How much faster are analytical derivatives with Ceres?
Lets try a simple optimization problem. Fitting a circle to 2D data
$$ r^2 = x^2 + y^2 $$
The measurement residual is defined as the distance to the circle boundary for each point i
$$ f(x, y, r, x_i, y_i) = \sqrt{(x - x_i)^2 + (y - y_i)^2} - r $$
Taking the partial derivatives for each parameter we get
$$ \frac{\partial f}{\partial x} = \frac{ x - x_i }{\sqrt{(x-x_i)^2 + (y-y_i)^2}} $$
$$ \frac{\partial f}{\partial y} = \frac{ y - y_i }{\sqrt{(x-x_i)^2 + (y-y_i)^2}} $$
$$ \frac{\partial f}{\partial r} = -1 $$
These derivatives are a bit complex. It would be best if we avoided the square root in the residual
We can redefine the residual as
$$ f(x, y, r, x_i, y_i) = r^2 - (x - x_i)^2 + (y - y_i)^2 $$
Taking the partial derivatives for each parameter we get
$$ \frac{\partial f}{\partial x} = 4r^2 $$
$$ \frac{\partial f}{\partial y} = -2(x-x_i) $$
$$ \frac{\partial f}{\partial r} = -2(y-y_i) $$
We employ a little trick. We define m
as the square root of r
.
By optimizing with parameter m
it is impossible to get a negative radius.
With this information we can implement two cost functions.
CircleFitAuto
which relies on Ceres automatic Jet derivatives, and CircleFitAnalytical
which analytically computes the derivatives each step.
class CircleFitAuto
{
public:
CircleFitAuto(const double x, const double y) : x_(x), y_(y) {}
template <typename T>
bool operator()(const T* const x, const T* const y,
const T* const m,
T* residual) const
{
const T r = *m * *m;
const T xp = x_ - *x;
const T yp = y_ - *y;
residual[0] = r * r - xp * xp - yp * yp;
return true;
}
private:
const double x_;
const double y_;
};
struct CircleFitAnalytical : public ceres::SizedCostFunction<1, 1, 1, 1>
{
CircleFitAnalytical(const double x, const double y) : x_(x), y_(y) {}
virtual ~CircleFitAnalytical() = default;
virtual bool Evaluate(double const* const* parameters, double* residuals,
double** jacobians) const
{
const double x = parameters[0][0];
const double y = parameters[1][0];My conclusion
const double m = parameters[2][0];
const double r = m * m;
const double xp = x_ - x;
const double yp = y_ - y;
residuals[0] = r * r - xp * xp - yp * yp;
if (jacobians != nullptr && jacobians[0] != nullptr) {
jacobians[0][0] = 2 * (x_ - x);
jacobians[1][0] = 2 * (y_ - y);
jacobians[2][0] = 4 * m * m * m;
}
return true;
}
private:
const double x_;
const double y_;
};
We can generate points around the unit circle with normal noise (stddev = 0.1)
std::random_device rd{};
std::mt19937 gen{ rd() };
std::normal_distribution<double> dist{ 0, 0.1 };
std::vector<std::array<double, 2>> data;
const int samples = 1e2;
const double angle_step = 2 * M_PI / samples;
for (int i = 0; i < samples; ++i) {
const double theta = angle_step * i;
const double radius = 1.0 + dist(gen);
data.push_back({ std::cos(theta) * radius, std::sin(theta) * radius });
}
Visualising the data
Running the optimization with different sample sizes
The results are compelling. For this problem the analytical derivative is 25-30x faster than Ceres auto diff Jets.
However…
I mentioned these results to a colleague and they suggested compiling with optimizations…
Sure enough, compiling with optimizations ('-O2') and the speedup completely disappears!
So then, how much faster are analytical derivatives with Ceres? In the end, not noticeably faster, at least not for the example above!