Make an iteration matrix explain convergence

Try

Solve 5x+y=11 and x+4y=9 by rewriting x=(11−y)/5 and y=(9−x)/4. Jacobi uses both old coordinates to compute the next pair. Gauss–Seidel uses the freshly updated x when computing y. Starting at zero, Jacobi gives (11/5,9/4), while Gauss–Seidel gives (11/5,17/10). The order of updates is part of the algorithm, not merely a programming detail.

Predict

Predict whether repeated substitution always solves a nonsingular system. It does not. A unique exact solution can coexist with a divergent stationary iteration. The splitting of A determines an iteration matrix, and that matrix—not invertibility alone—controls convergence. Also predict whether a small update guarantees a small error: only with a suitable contraction estimate or residual-and-conditioning argument does it become a certificate.

Experiment

Compare the first sweep with an explicit dependency table.

Method	x update uses	y update uses	First pair
Jacobi	old y	old x	(11/5,9/4)
Gauss–Seidel	old y	new x	(11/5,17/10)

Perform a second sweep for each method. Keep separate old and new arrays for Jacobi; accidentally overwriting the old x turns the implementation into a different iteration.

Understand

Write A=M−N with M invertible. The stationary method is xₖ₊₁=M⁻¹Nxₖ+M⁻¹b, so B=M⁻¹N. Subtracting the exact fixed point gives eₖ=B^ke₀. Convergence for every start is equivalent to B^k→0, which in finite dimension is equivalent to spectral radius ρ(B)<1. Necessity follows from eigenvectors. Sufficiency follows from Jordan form: powers of a block involve polynomial factors in k multiplied by λ^k, tending to zero for |λ|<1, while nilpotent blocks eventually vanish. A simpler sufficient test is ‖B‖<1 in an induced norm; it need not be necessary for that particular norm.

Apply

For A=D+L+U, Jacobi uses B=−D⁻¹(L+U); Gauss–Seidel uses B=−(D+L)⁻¹U. These expressions describe triangular solves, not a recommendation to form explicit inverses. Strict row diagonal dominance gives Jacobi infinity norm below one because each row sum is Σⱼ≠ᵢ|aᵢⱼ|/|aᵢᵢ|. Other conditions, including symmetric positive definiteness for Gauss–Seidel, can establish convergence by different arguments. Convergence speed depends on the spectrum and possible transient behavior, not just whether the eventual limit exists.

A diagonal-dominance test is sufficient, not necessary. An iteration matrix can have spectral radius below one while its infinity norm exceeds one, so failure of that easy bound does not prove divergence. Conversely, observing ten decreasing errors does not prove eventual convergence for an arbitrary nonnormal matrix. Use a theorem appropriate to the splitting or report the observed behavior without promoting it to a guarantee. The residual can also be checked cheaply during each sweep, with its extra work included in the cost model.

Continue

Implement bounded sweeps, residual monitoring and explicit nonconvergence status. Repeated iterates are development evidence, not a theorem about every right-hand side. For a known contraction factor q, the latest-step bound q/(1−q) times the step norm estimates remaining error; otherwise use residual and conditioning carefully. Large sparse systems often need Krylov methods and preconditioning beyond this stationary foundation. Next, eigenvalue iteration uses similar repeated matrix actions but normalizes vectors and targets invariant directions rather than a fixed solution vector.

Review calculus and linear algebra when a prerequisite needs practice. Continue to Extract spectral information from repeated actions.