101
0/318
Loading content...
The Conjugate Gradient (CG) method is one of the most elegant and powerful iterative algorithms in numerical linear algebra. It is specifically designed to solve systems of linear equations of the form Ax = b, where A is a symmetric positive-definite (SPD) matrix. Unlike direct methods such as Gaussian elimination or LU decomposition, CG operates iteratively, making it highly efficient for large, sparse systems where direct methods become computationally prohibitive.
The CG method is rooted in optimization theory. Solving Ax = b when A is SPD is mathematically equivalent to minimizing the quadratic functional:
$$f(x) = \frac{1}{2}x^\top A x - b^\top x$$
The minimum of this function occurs precisely at the solution x* = A⁻¹b. The CG algorithm navigates through the solution space by generating a sequence of A-conjugate search directions, ensuring that each iteration makes progress toward the solution without undoing previous work.
The algorithm maintains several key vectors at each iteration:
At each iteration, the algorithm:
For an n×n SPD matrix, the CG method is guaranteed to converge to the exact solution (in exact arithmetic) in at most n iterations. In practice, the method often converges much faster, especially when the matrix has clustered eigenvalues or when a preconditioner is applied. The convergence rate depends on the condition number κ(A) = λ_max / λ_min, where λ_max and λ_min are the largest and smallest eigenvalues of A.
Implement a function that solves the linear system Ax = b using the Conjugate Gradient method. Your implementation should:
A = [[4.0, 1.0], [1.0, 3.0]]
b = [1.0, 2.0]
max_iterations = 5
convergence_threshold = 1e-8[0.09090909, 0.63636364]The 2×2 matrix A is symmetric positive-definite (both eigenvalues are positive). Starting from the zero vector, the CG method computes:
Iteration 1: • Initial residual r₀ = b - A·x₀ = [1.0, 2.0] (since x₀ = [0, 0]) • Initial search direction p₀ = r₀ = [1.0, 2.0] • Compute A·p₀ = [4·1 + 1·2, 1·1 + 3·2] = [6, 7] • Step length α₀ = (r₀·r₀) / (p₀·A·p₀) = 5 / 20 = 0.25 • Update x₁ = x₀ + α₀·p₀ = [0.25, 0.5] • Update r₁ = r₀ - α₀·A·p₀ = [1, 2] - 0.25·[6, 7] = [-0.5, 0.25]
Iteration 2: • The algorithm continues, refining the solution further • After 2 iterations (n = 2), the exact solution is reached
The final solution x = [0.09090909, 0.63636364] satisfies Ax = b with the given matrix.
A = [[4.0, 1.0, 0.0], [1.0, 5.0, 1.0], [0.0, 1.0, 4.0]]
b = [1.0, 2.0, 3.0]
max_iterations = 10
convergence_threshold = 1e-8[0.19444444, 0.22222222, 0.69444444]This is a tridiagonal SPD matrix, which is common in finite difference discretizations. The matrix structure:
A = | 4 1 0 | | 1 5 1 | | 0 1 4 |
The off-diagonal elements represent coupling between adjacent variables, while the dominant diagonal ensures positive-definiteness. The CG method converges in at most 3 iterations (the matrix dimension), with the solution vector balancing the right-hand side values according to the system's structure.
A = [[1.0, 0.0, 0.0], [0.0, 1.0, 0.0], [0.0, 0.0, 1.0]]
b = [5.0, 10.0, 15.0]
max_iterations = 5
convergence_threshold = 1e-8[5.0, 10.0, 15.0]When A is the identity matrix, the system Ax = b simplifies to Ix = b, meaning x = b. The CG method converges in exactly 1 iteration:
• Starting from x₀ = [0, 0, 0], the residual r₀ = b = [5, 10, 15] • Since A = I, the product A·r₀ = r₀, and the step length α₀ = 1 • Update x₁ = x₀ + 1·r₀ = [5, 10, 15] — the exact solution
This demonstrates that for well-conditioned systems (condition number = 1 for identity), CG converges immediately.
Constraints