Conjugate gradient method without agonizing pain

64 Pages·1994·0.492 MB·English

by  Schewchuk.

An Introduction to
the Conjugate Gradient Method
Without the Agonizing Pain
Edition11
4
JonathanRichardShewchuk
August4,1994
SchoolofComputerScience
CarnegieMellonUniversity
Pittsburgh,PA15213
Abstract
TheConjugateGradientMethodisthemostprominentiterativemethodforsolvingsparsesystemsoflinearequations.
Unfortunately, many textbook treatments of the topic are written with neither illustrations nor intuition, and their
victims can be found to this day babbling senselessly in the corners of dusty libraries. For this reason, a deep,
geometricunderstandingofthemethodhasbeenreservedfortheelitebrilliantfewwhohavepainstakinglydecoded
themumblingsoftheirforebears. Nevertheless,theConjugateGradientMethodisacompositeofsimple,elegantideas
thatalmostanyonecanunderstand. Ofcourse,areaderasintelligentasyourselfwilllearnthemalmosteffortlessly.
TheideaofquadraticformsisintroducedandusedtoderivethemethodsofSteepestDescent,ConjugateDirections,
and Conjugate Gradients. Eigenvectorsare explained and used to examine the convergenceof the Jacobi Method,
SteepestDescent,andConjugateGradients. OthertopicsincludepreconditioningandthenonlinearConjugateGradient
Method. I have taken pains to make this article easy to read. Sixty-six illustrations are provided. Dense prose is
avoided. Conceptsareexplainedinseveraldifferentways. Mostequationsarecoupledwithanintuitiveinterpretation.
SupportedinpartbytheNaturalSciencesandEngineeringResearchCouncilofCanadaundera1967ScienceandEngineering
Scholarshipand by the National ScienceFoundationunder Grant ASC-9318163. The views and conclusionscontainedin this
documentarethoseoftheauthorandshouldnotbeinterpretedasrepresentingtheofficialpolicies,eitherexpressorimplied,of
NSERC,NSF,ortheU.S.Government.
Keywords: conjugategradientmethod,preconditioning,convergenceanalysis,agonizingpain
Contents
1. Introduction 1
2. Notation 1
3. TheQuadraticForm 2
4. TheMethodofSteepestDescent 6
5. ThinkingwithEigenvectorsandEigenvalues 9
(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:3)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)
5.1. EigendoitifItry 9
(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:3)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)
5.2. Jacobiiterations 11
(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:3)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)
5.3. AConcreteExample 12
6. ConvergenceAnalysisofSteepestDescent 13
(cid:0)(cid:3)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:3)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)
6.1. InstantResults 13
(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:3)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)
6.2. GeneralConvergence 17
7. TheMethodofConjugateDirections 21
(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:3)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:3)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)
7.1. Conjugacy 21
(cid:0)(cid:3)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)
7.2. Gram-SchmidtConjugation 25
(cid:0)(cid:3)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)
7.3. OptimalityoftheErrorTerm 26
8. TheMethodofConjugateGradients 30
9. ConvergenceAnalysisofConjugateGradients 32
(cid:0)(cid:3)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)
9.1. PickingPerfectPolynomials 33
(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:3)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)
9.2. ChebyshevPolynomials 35
10.Complexity 37
11.StartingandStopping 38
(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:3)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:3)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)
11.1. Starting 38
(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:3)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:3)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)
11.2. Stopping 38
12.Preconditioning 39
13.ConjugateGradientsontheNormalEquations 41
14.TheNonlinearConjugateGradientMethod 42
(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)
14.1. OutlineoftheNonlinearConjugateGradientMethod 42
(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:3)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)
14.2. GeneralLineSearch 43
(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:3)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)
14.3. Preconditioning 47
A Notes 48
B CannedAlgorithms 49
(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:3)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)
B1. SteepestDescent 49
(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:3)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)
B2. ConjugateGradients 50
(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)
B3. PreconditionedConjugateGradients 51
i
(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)
B4. NonlinearConjugateGradientswithNewton-RaphsonandFletcher-Reeves 52
(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)
B5. PreconditionedNonlinearConjugateGradientswithSecantandPolak-Ribie`re 53
(cid:4)(cid:6)(cid:5)(cid:8)(cid:7)(cid:10)(cid:9)
C UglyProofs 54
(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)
C1. TheSolutionto MinimizestheQuadraticForm 54
(cid:11) (cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)
C2. ASymmetricMatrixHas OrthogonalEigenvectors. 54
(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:1)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)(cid:3)(cid:0)(cid:1)(cid:0)(cid:2)(cid:0)
C3. OptimalityofChebyshevPolynomials 55
D Homework 56
ii
AboutthisArticle
An electronic copy of this article is available by anonymous FTP to WARP.CS.CMU.EDU (IP address
128.2.209.103) under the filename quake-papers/painless-conjugate-gradient.ps. A
PostScript file containing full-page copies of the figures herein, suitable for transparencies, is available
as quake-papers/painless-conjugate-gradient-pics.ps. Most of the illustrations were
createdusingMathematica.
(cid:12)
c1994byJonathanRichardShewchuk. Thisarticlemaybefreelyduplicatedanddistributedsolongas
noconsiderationisreceivedinreturn,andthiscopyrightnoticeremainsintact.
ThisguidewascreatedtohelpstudentslearnConjugateGradientMethodsaseasilyaspossible. Please
mail me ([email protected]) comments, corrections, and any intuitions I might have missed; some of
these will be incorporated into a second edition. I am particularly interested in hearing about use of this
guideforclassroomteaching.
Forthosewhowishtolearnmoreaboutiterativemethods,IrecommendWilliamL.Briggs’“AMultigrid
Tutorial”[2],oneofthebest-writtenmathematicalbooksIhaveread.
Special thanks to Omar Ghattas, who taught me much of what I know about numerical methods, and
provided me with extensive comments on the first draft of this article. Thanks also to James Epperson,
DavidO’Hallaron,JamesStichnoth,NickTrefethen,andDanielTunkelangfortheircomments.
Tohelpyouskipchapters,here’sadependencegraphofthesections:
1 Introduction 10 Complexity
2 Notation 13 Normal Equations
3 Quadratic Form
4 Steepest Descent 5 Eigenvectors
7 Conjugate Directions 6 SD Convergence
8 Conjugate Gradients 9 CG Convergence
11 Stop & Start 12 Preconditioning 14 Nonlinear CG
ThisarticleisdedicatedtoeverymathematicianwhousesfiguresasabundantlyasIhaveherein.
iii
iv
1. Introduction
WhenIdecidedtolearntheConjugateGradientMethod(henceforth,CG),Ireadfourdifferentdescriptions,
whichIshallpolitelynotidentify. Iunderstoodnoneofthem. Mostofthemsimplywrotedownthemethod,
thenproveditspropertieswithoutanyintuitiveexplanationorhintofhowanybodymighthaveinventedCG
inthefirstplace. Thisarticlewasbornofmyfrustration,withthewishthatfuturestudentsofCGwilllearn
arichandelegantalgorithm,ratherthanaconfusingmassofequations.
CG is the most popular iterative methodfor solving large systemsof linear equations. CG is effective
(cid:4)(cid:13)(cid:5)(cid:14)(cid:7)(cid:15)(cid:9)
forsystemsoftheform
(1)
(cid:5) (cid:9) (cid:4)
where isanunknownvector, isaknownvector,and isaknown,square,symmetric,positive-definite
(or positive-indefinite)matrix. (Don’t worry if you’ve forgotten what “positive-definite” means; we shall
review it.) These systems arise in many important settings, such as finite difference and finite element
methodsforsolvingpartialdifferentialequations,structuralanalysis,circuitanalysis,andmathhomework.
(cid:4)
(cid:4)
Iterative methods like CG are suited for use with sparse matrices. If is dense, your best course of
(cid:4) (cid:4)
actionisprobablytofactor andsolvetheequationbybacksubstitution. Thetimespentfactoringadense
(cid:9)
isroughlyequivalenttothetimespentsolvingthesystemiteratively;andonce isfactored,thesystem
(cid:4)
can be backsolved quickly for multiple values of . Compare this dense matrix with a sparse matrix of
(cid:4)
larger size that fills the same amount of memory. The triangular factors of a sparse usually have many
more nonzero elements than itself. Factoring may be impossible due to limited memory, and will be
time-consuming as well; even the backsolving step may be slower than iterative solution. On the other
hand,mostiterativemethodsarememory-efficientandrunquicklywithsparsematrices.
I assume that you have taken a first course in linear algebra, and that you have a solid understanding
of matrix multiplication and linear independence, although you probably don’t remember what those
eigenthingieswereallabout. Fromthisfoundation,IshallbuildtheedificeofCGasclearlyasIcan.
2. Notation
Letusbeginwithafewdefinitionsandnotesonnotation.
(cid:4) (cid:5) (cid:9)
Withafewexceptions,Ishallusecapitalletterstodenotematrices,lowercaseletterstodenotevectors,
(cid:11)(cid:17)(cid:16)(cid:18)(cid:11) (cid:11)(cid:17)(cid:16)
andGreekletterstodenotescalars. isan matrix,and and arevectors—thatis, 1matrices.
Writtenoutfully,Equation1is
(cid:19)(cid:20) (cid:4) (cid:4) (cid:4) (cid:25)(cid:27)(cid:26) (cid:19)(cid:20) (cid:5) (cid:25)(cid:27)(cid:26) (cid:19)(cid:20) (cid:9) (cid:25)(cid:27)(cid:26)
(cid:20)(cid:20) (cid:0)(cid:22)(cid:0)(cid:22)(cid:0) (cid:23) (cid:26)(cid:26) (cid:20)(cid:20) (cid:26)(cid:26) (cid:20)(cid:20) (cid:26)(cid:26)
(cid:20) (cid:4) (cid:4) (cid:4) (cid:26) (cid:20) (cid:5) (cid:26) (cid:20) (cid:9) (cid:26)
11 12 1(cid:23) 1 (cid:7) 1
(cid:0)
21 22 2 2 2
(cid:21) (cid:4) ... (cid:4) ... (cid:4) ... (cid:28) (cid:21) (cid:5) ... (cid:28) (cid:21) (cid:9) ... (cid:28)
(cid:23) (cid:23) (cid:0)(cid:22)(cid:0)(cid:22)(cid:0) (cid:23)(cid:24)(cid:23) (cid:23) (cid:23)
1 2
(cid:5)(cid:30)(cid:29) (cid:31) (cid:23)"$# (cid:5) " (cid:31) "
!
(cid:5)(cid:30)(cid:29) (cid:31)(cid:18)(cid:7)%(cid:31)&(cid:29) (cid:5) (cid:5) (cid:31) (cid:5)(cid:30)(cid:29) (cid:31)(cid:18)(cid:7)
The inner product of two vectors is written , and represents the scalar sum . Note that
1 (cid:16)
(cid:5)(cid:30)(cid:29) (cid:31) (cid:5)’(cid:29)((cid:4)(cid:6)(cid:5)
. If and areorthogonal,then 0. Ingeneral,expressionsthatreduceto1 1matrices,
suchas and ,aretreatedasscalarvalues.
(cid:5)
2 JonathanRichardShewchuk
2
4
2
(cid:5) (cid:5) (cid:7)
)
3 2 2 (cid:5)
1 2
1
-4 -2 2 4 6
(cid:5) (cid:5) (cid:7)+*
)
-2 2 1 6 2 8
-4
-6
Figure1: Sampletwo-dimensionallinearsystem. Thesolutionliesattheintersectionofthelines.
(cid:4) (cid:5)
Amatrix ispositive-definiteif,foreverynonzerovector ,
(cid:5) (cid:29) (cid:4)(cid:6)(cid:5)-,
(cid:0)
0 (2)
This may mean littleto you, butdon’t feel bad; it’s not a veryintuitive idea, and it’s hard to imaginehow
a matrixthat is positive-definitemight look differently fromone that isn’t. We will get a feeling for what
positive-definitenessisaboutwhenweseehowitaffectstheshapeofquadraticforms.
(cid:4)0/(cid:18)12(cid:29)3(cid:7)4/5(cid:29) (cid:4)6(cid:29) (cid:4)0/7198 (cid:7)4/:8 (cid:4)(cid:3)8
. .
Finally,don’tforgettheimportantbasicidentities and 1 1 1.
3. TheQuadraticForm
Aquadraticformissimplyascalar,quadraticfunctionofavectorwiththeform
; (cid:5)<1=(cid:7) (cid:5) (cid:29) (cid:4)(cid:6)(cid:5)>*?(cid:9) (cid:29) (cid:5)
. 1 )A@
(3)
(cid:4) (cid:5) (cid:9) 2 (cid:4)
@
; (cid:5)<1 (cid:4)(cid:13)(cid:5)(cid:14)(cid:7)(cid:15)(cid:9)
where isamatrix, and arevectors,and isascalarconstant. Ishallshowshortlythatif issymmetric
.
andpositive-definite, isminimizedbythesolutionto .
Throughoutthispaper,Iwilldemonstrateideaswiththesimplesampleproblem
(cid:4)4(cid:7)CB (cid:9)F(cid:7)CB (cid:7)
3 2 * 2 @ (cid:0)
D(cid:2)E DGE 0 (4)
2 6 8
(cid:4)(cid:13)(cid:5)H(cid:7)I(cid:9) (cid:5)
* (cid:5)%(cid:7)KJ * (cid:29)
The system is illustrated in Figure 1. In general, the solution lies at the intersection point
(cid:11) (cid:11) L
of hyperplanes, each having; d(cid:5)Mim1 ension 1. For this problem, the so; lu(cid:5)<ti1on is 2E 2 . The
. .
correspondingquadratic form appears in Figure 2. A contour plot of is illustrated in Figure 3.
TheQuadraticForm 3
150
; (cid:5)<1
100 .
4
50
2
0
(cid:5) 0
66
2 -2 44
(cid:5)(cid:5) 22
-4 00
11
--22
-6
--44
N<OQPSR TUPWVYX
Figure2: Graphofaquadraticform . Theminimumpointofthissurfaceisthesolutionto .
(cid:5)
2
4
2
(cid:5)
1
-4 -2 2 4 6
-2
-4
-6
N<O$PZR
Figure3: Contoursofthequadraticform. Eachellipsoidalcurvehasconstant .
(cid:5)
4 JonathanRichardShewchuk
2
6
4
2
(cid:5)
-4 -2 2 4 6 1
-2
-4
-6
-8
NZ[\O$PZR P
Figure4: Gradient ofthequadraticform. Forevery ,thegradientpointsinthedirectionofsteepest
N<O$PZR
increaseof ,andisorthogonaltothecontourlines.
(cid:4) ; (cid:5)M1
.
Because ispositive-definite,thesurfacedefinedby isshapedlikeaparaboloidbowl. (I’llhavemore
tosayaboutthisinamoment.)
Thegradientofaquadraticformisdefinedtobe
(cid:19)(cid:20) _ ; (cid:5)<1 (cid:25)(cid:26)
(cid:20) _(cid:22)‘ . (cid:26)
;^] (cid:5)M1=(cid:7) (cid:20)(cid:20) _(cid:22)_‘ 1; . (cid:5)<1 (cid:26)(cid:26)
. (cid:0)
(cid:21) _a‘a_ b2;... . (cid:5)M1 (cid:28) (5)
(cid:5) ; (cid:5)<1
.
The gradient is a vectorfield that, for a given point , pointsin the directionof greatestincrease of .
Figure4illustratesthegradientvectorsforEquation3withtheco; n(cid:5)<st1antsgivenin; E] q(cid:5)Mu1 ation4. Atthebottom
. .
oftheparaboloidbowl,thegradientiszero. Onecanminimize bysetting equaltozero.
Withalittlebitoftediousmath,onecanapplyEquation5toEquation3,andderive
; ] (cid:5)M1c(cid:7) (cid:4) (cid:29) (cid:5) (cid:4)(cid:6)(cid:5):*?(cid:9)
. 1 ) 1 (cid:0)
(6)
(cid:4) 2 2
If issymmetric,thisequationreducesto
; ] (cid:5)M1c(cid:7)4(cid:4)(cid:13)(cid:5)d*?(cid:9)
. (cid:0)
(7)
(cid:4)(cid:13)(cid:5)%(cid:7)e(cid:9) ; (cid:5)M1 (cid:4)
Setting the gradient to zero, we obtain Equation 1, the linear system we wish to solve. Therefore, the
.
solution to is a critical point of . If is positive-definite as well as symmetric, then this