Schur complement / Sparsity

Bundle adjustment looks intractable at first glance: a modest map with $m$ keyframes and $n$ points has $6m + 3n$ unknowns, and $n$ is easily in the hundreds of thousands. What makes it practical is that the problem is sparse, and the Schur complement is the trick that exploits that sparsity.

Where the sparsity comes from. Each reprojection-error term involves exactly one pose and one point. So in the Gauss-Newton normal equations $H \Delta\mathbf{x} = -\mathbf{b}$ with $H = J^T J$, the Hessian has the arrow-like block structure

$$H = \begin{bmatrix} B & E \\ E^T & C \end{bmatrix}$$

where $B$ ($6m \times 6m$) couples only poses, $C$ ($3n \times 3n$) couples only points, and $E$ holds the pose-point coupling. Crucially, two points never appear in the same residual, so $C$ is block-diagonal — one independent $3 \times 3$ block per point.

The Schur complement step (also called marginalizing out the points) eliminates the point variables from the linear system, leaving the reduced camera system:

$$\left(B - E C^{-1} E^T\right) \Delta\mathbf{x}_{\text{cam}} = -\mathbf{b}_{\text{cam}} + E C^{-1} \mathbf{b}_{\text{pts}}$$

Because $C$ is block-diagonal, $C^{-1}$ costs almost nothing (invert each $3 \times 3$ block). One solves the $6m \times 6m$ reduced system for the pose update, then back-substitutes to recover each point update independently. A $(6m + 3n)$-dimensional solve becomes a $6m$-dimensional one — orders of magnitude faster when $n \gg m$, which is always.

Two further layers of sparsity matter in practice. First, the reduced camera matrix $B - EC^{-1}E^T$ is itself sparse: entry $(i, j)$ is non-zero only if keyframes $i$ and $j$ observe a common point (the covisibility structure), so sparse Cholesky factorization with good variable ordering (COLAMD) applies. Second, the same elimination viewpoint generalises: in factor-graph terms, the Schur complement is just variable elimination, the operation underlying marginalization in sliding-window VIO and incremental smoothers like iSAM2.

Where the formula comes from

The Schur complement is nothing more exotic than block Gaussian elimination. Write out the two block rows of the normal equations:

$$B\,\Delta\mathbf{x}_{\text{cam}} + E\,\Delta\mathbf{x}_{\text{pts}} = -\mathbf{b}_{\text{cam}}$$ $$E^T \Delta\mathbf{x}_{\text{cam}} + C\,\Delta\mathbf{x}_{\text{pts}} = -\mathbf{b}_{\text{pts}}$$

Solve the second row for the point update, $\Delta\mathbf{x}_{\text{pts}} = -C^{-1}\left(\mathbf{b}_{\text{pts}} + E^T \Delta\mathbf{x}_{\text{cam}}\right)$, and substitute into the first — the reduced camera system above appears immediately, and the same expression is the back-substitution rule once $\Delta\mathbf{x}_{\text{cam}}$ is known. Statistically, the reduced system is the information-matrix form of the marginal over the cameras: eliminating the points does not approximate anything, it re-expresses the same Gaussian exactly.

Counting the win. A dense solve of the full system costs $O((6m + 3n)^3)$. With the Schur trick: inverting $C$ is $n$ independent $3\times 3$ inversions ($O(n)$), forming the reduced system is bounded by the number of observations, and the remaining solve is $O((6m)^3)$ worst-case — usually far less thanks to covisibility sparsity. With $m$ in the hundreds and $n$ in the hundreds of thousands, the point elimination is the difference between milliseconds and minutes. When $m$ itself grows large (city-scale SfM), even the reduced system stops being cheap to factor, and solvers switch to iterative methods (preconditioned conjugate gradient on the Schur complement — Ceres’ ITERATIVE_SCHUR).

In Ceres, all of this is one configuration choice:

ceres::Solver::Options options;
options.linear_solver_type = ceres::SPARSE_SCHUR;  // or DENSE_SCHUR, ITERATIVE_SCHUR

with the solver automatically detecting the pose/point elimination ordering (or taking an explicit ParameterBlockOrdering). Every serious solver — Ceres, g2o, GTSAM — implements the same trick; knowing it explains why BA scales, why marginalizing old states creates fill-in (dense blocks) in the remaining Hessian, and why solver choice and ordering can change runtime by orders of magnitude.

Common pitfalls

• Eliminating in the wrong order — the trick works because points are many, cheap, and mutually independent; eliminating poses first couples all their points together and densifies the system. Elimination order is a correctness-of-performance decision.
• Confusing marginalization with deletion — dropping old states and marginalizing them are different: marginalization keeps their information as a dense prior over their neighbours (fill-in). That density is why sliding-window VIO systems carefully choose what to marginalize versus discard.
• Gauge freedom — full BA has 7 unobservable DoF for monocular (6 + scale); without fixing a pose (or adding priors), the reduced camera matrix is singular and the solver limps or fails.
• Ill-conditioning from far points — points at near-infinite depth make their $3\times3$ blocks nearly singular; inverse-depth parameterisation or depth priors keep $C^{-1}$ well-behaved.

Why it matters for SLAM

Real-time SLAM exists because of structure-exploiting linear algebra: without the Schur complement, bundle adjustment over even a few hundred keyframes would be hopeless. The same idea — eliminate variables, keep the problem sparse — reappears in marginalization, sliding-window estimators, and incremental smoothing, so this is one of the highest-leverage pieces of math on the whole roadmap.

Related

• MAP inference as sparse nonlinear least squares
• Pose graph optimization
• Marginalization
• Math libraries