Aggregating what the matrix calls strong

The depth phase ended the anisotropy field on a failure and a diagnosis. Four methods on the 45°-rotated anisotropic operator — a point smoother, a line smoother, semi-coarsening, and classical algebraic multigrid — returning convergence factors of 0.784, 0.883 and 0.844, none of them usable. And the reason, measured before any solver ran: the standard nine-point discretisation of the rotated operator has axis couplings of −0.5005 and a coupling along the strong direction of −0.2498, so the matrix does not contain the anisotropy the equation has.

That diagnosis makes a prediction, and the prediction is the reason this essay exists. If the information was lost in the discretisation, then no better hierarchy can recover it, because every algebraic method reads the matrix. Smoothed aggregation is the standard answer to anisotropy — the method a practitioner reaches for after classical AMG — so it is the sharpest available test of a claim about where the failure lives.

It fails, at 0.789, on the operator where it returns 0.193 unrotated.

The method, in four steps, and only one of them matters #

Strength. The same graph classical AMG builds: j is strong for i when −aᵢⱼ is at least θ times the largest off-diagonal magnitude in row i. Reused rather than rewritten, so the two methods disagree about the hierarchy and not about the matrix.

Aggregation. Cover the points with disjoint aggregates. Vaněk’s two passes: seed an aggregate at any point whose whole strong neighbourhood is still free, then hand each leftover point to the neighbouring aggregate it is most strongly attached to. Where classical AMG chooses a subset of the points to keep, this chooses a partition of them — which is the difference the two names are about. That the result is a partition is checked rather than assumed: every point appears in exactly one aggregate on all three operators, and the tentative prolongator’s columns are orthonormal to 10⁻¹².

A tentative prolongator. One column per aggregate, holding a normalised constant. A constant is what these operators annihilate, and an interpolation that cannot represent it cannot represent the error a smoother leaves behind.

And one weighted Jacobi sweep applied to the columns of that prolongator, with ω = 4/(3ρ(D⁻¹A)). This is the method.

The sweep is the method #

Same aggregates, same cycle, same smoother, same number of levels. One weighted Jacobi sweep applied to each column of P, and a factor of 2.3 in the convergence rate.

What the sweep does is widen the columns. A piecewise constant that is one on an aggregate and zero on its neighbour becomes a function that decays across the boundary — so the coarse space contains something with a gradient in it, and the error a smoother leaves behind has gradients. Unsmoothed aggregation asks a piecewise-constant space to represent a smooth error, and 0.813 is what that costs.

The constant is worth a sentence because this site has met it from the other side. ω = 4/(3ρ) is Vaněk’s choice and it damps the top two thirds of the spectrum; on the one-dimensional Laplacian, ρ(D⁻¹A) = 2 and 4/(3ρ) is exactly ⅔ — which is the smoother weight the multigrid field measured as optimal, where its optimum is exactly ⅓ and a 2,000-point search misses it and reports 0.3335.

What the aggregates know #

With ε = 10⁻³ the x-couplings are 10⁻³ and the y-couplings are 1, so the strength test sees the y direction and nothing else. Measured on a 31×31 grid: the aggregates’ mean extent is 2.82 points in y and 1.00 in x, over 341 aggregates.

That is semi-coarsening — every grid column kept whole or dropped whole — derived from the entries of the matrix, by a method with no coordinate in it anywhere. The depth phase found the same thing for classical AMG and recorded it as the phase’s organising sentence: the algebraic method rediscovers what the geometric one was told. It survives the change of method, which is what makes it a property of the strength graph rather than of one coarsening algorithm.

Which direction counts as strong is read from coefficients(ε, 0) rather than assumed, and that is not pedantry: the two-dimensional operator’s coefficients are a = ε along x and c = 1 along y, so ε ≪ 1 makes y the strong direction, and the first version of this measurement asserted the extents the other way round and failed. A measurement of the right number attached to the wrong axis is the easiest mistake in this field to make and the hardest to see afterwards.

And the interpolation carries a constant, in the interior #

P₀ carries a constant exactly by construction — the coarse vector whose entries are √(aggregate size) maps to the all-ones vector, checked to 1.1·10⁻¹⁶.

The smoothed P does not, quite, and the measurement is worth reporting rather than glossing: (I − ωD⁻¹A) applied to a constant leaves it alone wherever the row sum of A is zero, which is every interior row, and changes it where the row sum is not — the boundary. Measured: the deviation from one is 1.1·10⁻¹⁶ over the interior points and 0.34 at the worst boundary point.

So the near-null-space property the method rests on is exact where the operator is homogeneous and broken at the edge, by an amount that is a boundary effect. That is the same shape as the finding two essays away in the structure field, where a clustering result thins because its exception set is the boundary of the grid — and it is worth noticing that two unrelated methods have their one-dimensional exactness spoiled by the same geometric fact.

The rotation, and the prediction it tests #

The rotated operator’s coefficients are a = c = (ε + 1)/2 = 0.5005 and b = (ε − 1)/2 = −0.4995, so the stencil has four axis entries of −0.5005 and four diagonal entries of ∓0.2498. Two consequences, and both of them break the strength test rather than the solver:

The largest off-diagonal in every interior row is an axis coupling. The strength test compares −aᵢⱼ against the largest in the row, so it marks the axis directions strong and the diagonal directions weak — and the diagonal directions are where the anisotropy runs.

And two of the four diagonal entries are positive. Every classical strength measure ignores positive off-diagonals entirely, for a good reason amg.js states: for a matrix with zero row sums, a positive off-diagonal is a coupling that pushes two values apart rather than together, and treating it as strong produces interpolation weights of the wrong sign. So half of the couplings along the strong direction are not merely under-weighted, they are invisible.

A method reading that matrix has no way to find the direction, because the direction is not in the matrix. Smoothed aggregation reads that matrix. It returns 0.789.

What a failure that was predicted is worth #

Two things, and it is worth separating them from the disappointment.

The diagnosis is confirmed rather than repeated. The depth phase’s claim was that the failure at 45° belongs to the discretisation. Its evidence was that four methods failed and that the stencil’s couplings pointed the wrong way. Adding a fifth failing method is weak evidence on its own; adding the one method whose entire purpose is anisotropy, which succeeds on the same anisotropy when it is axis-aligned, is much stronger. The variable that changed is the rotation and not the method.

And the fix is now locatable. If the information is missing from the matrix, the repair has to put it back before any solver sees it — which means the discretisation. A nine-point stencil aligned with the anisotropy rather than with the grid, or a mesh that follows the flow, or a stabilisation that knows the direction, all do that; a better hierarchy cannot. That is a conclusion about where to spend effort, and it is the kind of conclusion a failure is good for.

This is the same move the previous essay makes in the convection field, and the two are worth reading together: there, a scheme succeeds spectacularly because the answer’s shape was built into it, and fails on the problem next door. Here, a method fails because the operator’s shape was lost before it arrived. Both are the depth phase’s sentence — where a method’s competence was put into it by hand — measured from the two sides.

What it costs when it works #

The honest accounting, since a method that converges at 0.193 is not free:

Operator complexity is the total nonzeros across the hierarchy divided by the fine matrix’s, and it is the quantity that decides whether an optimal-complexity method is actually cheap. The smoothing widens P, which densifies every coarse operator below it, and 1.34 is a good number — classical AMG on the same grid runs higher, and the depth phase measured 17.7 on a graph Laplacian, with one level of forty-one unknowns entirely dense.

So this method is cheaper than classical AMG per cycle and converges more slowly per cycle, which is the standard trade between them and is measured here rather than repeated.

What the two passes each do, in counts #

The aggregation is two passes and the split between them is a measurement rather than an implementation detail — a method whose second pass is doing most of the work is a method whose first pass found little structure.

Read the last row against the first two. On the rotated operator the first pass opens thirty aggregates and the second pass has to place forty-five leftover points — half again as many as were seeded — where on the aligned operator it places fifteen. The strength graph the first pass walks is finding fewer clean neighbourhoods, because the couplings it can see are the axis ones and those do not agree about which points belong together.

That is a diagnostic available before any solve, and it is the second of the two the stencil supplies. A method whose aggregation is mostly second-pass adoption is a method operating on a graph that does not describe the problem.

The near-null space, and what “carries a constant” means precisely #

The tentative prolongator carries a constant exactly: the coarse vector whose entries are √(aggregate size) maps to the all-ones vector, to 1.1·10⁻¹⁶. That is by construction — the columns are normalised indicators of a partition — and it is why the partition property is checked rather than assumed.

The smoothed prolongator does not, quite, and the measurement is worth stating rather than glossing. (I − ωD⁻¹A) applied to a constant leaves it alone wherever the row sum of A is zero, which is every interior row of these operators, and changes it where the row sum is not — the boundary. Measured on the aligned operator: the deviation from one is 1.1·10⁻¹⁶ over the interior and 0.34 at the worst boundary point.

So the property the method’s theory rests on is exact in the interior and broken at the edge by a bounded amount. Both halves are asserted, and the reason to assert the second is that it is the difference between a claim that is true of an infinite grid and a claim that is true of the matrix in front of the method.

That is the same shape as the structure field’s loss two fields away, where a clustering result weakens in two dimensions because its exception set is the boundary of the grid rather than the ends of a vector. Two unrelated methods, the same geometric fact.

What the aggregation would have to see #

If the diagnosis is right — and the essay above is an experiment designed to test it rather than an argument for it — then the repair belongs in one of two places, and neither is the hierarchy.

A strength measure that sees positive couplings. Every classical measure compares −aᵢⱼ against the largest in the row and ignores positive off-diagonals, for a reason amg.js states: for a matrix with zero row sums a positive off-diagonal pushes two values apart. The rotated stencil has two positive diagonal entries per interior row, so half of the couplings along the anisotropy are not under-weighted but invisible. Affinity and algebraic-distance measures do see them, and whether they recover the 45° direction on this stencil is a measurement rather than an argument.

Or a discretisation that keeps the direction. The nine-point stencil used here is the standard one and it is where the information is lost. A stencil aligned with the anisotropy, or a mesh that follows it, would hand every method in this field a matrix that contains what the equation contains — and the measurement that matters is whether the four failing methods then return to the 0.19 they reach when the anisotropy happens to lie along an axis.

What is left #

Energy-minimising interpolation, the other half of the depth phase’s deferral, which chooses P by minimising the energy of its columns subject to reproducing the near-null space rather than by applying one Jacobi sweep. It is the natural next step from the smoothing sweep above — a sweep is one iteration of exactly that minimisation — and it is the method most likely to do better on the rotated operator, because it is the least dependent on the strength graph.

A strength measure that sees positive couplings. The affinity and algebraic-distance measures do, and whether they find the 45° direction in this stencil is a measurement rather than an argument. If they do, this essay’s conclusion needs narrowing to classical strength measures; if they do not, the diagnosis is stronger than stated.

And the discretisation itself, which is where the fix has to go and which this whole field has so far treated as given. The nine-point stencil used everywhere here is the standard one, and the measurement that matters is whether a stencil built to align with the anisotropy restores any of the four methods — including this one — to the 0.19 they reach when the anisotropy happens to lie along an axis.