230 likes | 345 Views
This document outlines a teaching session on the Conjugate Gradient Method (CGM) aimed at solving large sparse linear algebraic systems represented as Ax = b. It covers key topics such as the method of Gaussian elimination, basic iterative methods, and the derivation and theory behind the Conjugate Gradient Method. Additionally, it explores convergence criteria, optimization methods like steepest descent, and practical implementation strategies. Homework and projects are included for hands-on learning, emphasizing the application of the CGM using MATLAB.
E N D
2010年教学实践周7.12-7.16 Cojugate Gradient Method Zhengru Zhang (张争茹) zrzhang@bnu.edu.cn Office: Math. Building 413(West)
Outline • Aim • Method of Gauss Elimination • Basic Iterative Methods • Conjugate Gradient Method • Derivation • Theory • Algorithm • References • Homework & Project
Aim Solve linear algebraic system like a11 x1 + a12 x2 + ... + a1n xn = b1 a21 x1 + a22 x2 + ... + a2n xn = b2 ... an1 x1 + an2 x2 + ... + ann xn = bn Using matrix, the above system can be written as Ax=b A is a N x N matrix, b is a N x 1 vector Consider the case: A is large and sparse
Algorithm of Gaussian Elimination without Pivoting U = A, L = I for k = 1 to N-1 for j = k +1 to N ljk= ujk/ukk uj,k:m = uj,k:m – ljkuk,k:m • LU Factorization, let A=LU • Solve Ly=b • Solve Ux=y
Operation Count of Gauss Elimination • Gauss Elimination and Back Substitution • There are 3 loops • There are 2 flops per entry • For eachk, the inner loop is repeated for rows k +1, …, N • Cost: about About N 3flops
Instability of Gaussian Elimination without Pivoting Examples A2= A1= • Pivoting • Partial Pivoting • Complete Pivoting Remedy
Algorithm of Gaussian Elimination with Partial Pivoting U = A, L = I For k = 1 to N-1 for j = k +1 to N ljk= ujk/ukk uj,k:m = uj,k:m – ljkuk,k:m
Basic Iterative Methods • How to construct iterative sequence? • Convergence? Conditions? • Convergence rate?
Jacobi iteration Gauss Seidel iteration X[k+1] = D-1(L+U) X[k] + D-1 b B = D-1(L+U) X[k+1] = (D-L)-1 U X[k] + (D-L)-1 b B = (D-L)-1 U • Iterative method X[k+1] = BX[k]+g converges if and only if • (B) < 1 • Convergence rate • ||X[k]-X*|| ||X[1]-X[0]||, whereq =||B||<1
Steepest Decent Method • Consider the case: A is symmetric positive definite • Quadratic functional • (x)= xTAx - 2bTx • The solution of Ax=b is equivalent to find the minimizer • of the functional(x) • Method of optimization: find a direction pk and a step k
Steepest Decent Method Determine pkand k • Suppose that pk is determined. Let’s start from xk • Let f() =(xk + pk) • = (xk + pk)TA(xk + pk)-2bT(xk + pk) • = 2pkTApk - 2 rkTpk +(xk) • whererk = b - Axk(Residual) • By calculas f’() = 2pkTApk- 2rkTpk =0 • Then let xk+1 = xk + k pk
Algorithm for Steepest Decent Method • Verify(xk+1) - (xk) =(xk +k pk) - (xk) • = k2pkTApk - 2k rkTpk • How to determine the directionpk ? • take as the negative gradientpk = rk
Convergence Theorem Algorithm Suppose the eigenvalues of A then there holds where
Conjugate Gradient Method Derivation • Negative gradient direction rk is the locally steepest • decent direction, but it may not be the global one • Consider a new direction: combination of rk and pk-1 • Initially, take p0 = r0 , x1 = x0 + 0p0 • For step k +1, choose and to minimize • By calculas
The corresponding minimizer is • 0 and 0satisfy • take • Let • In summary, 0
Algorithm for CG method Where and are obtained in a simple form • Operations involved: • Transpose, • Scalar Multiply, • Matrix Add, • Matrix Multiply
Properties for CG method Orthogonal properties • Theoretically, CG method is an exact method. Actually, • works as an iterative method. • Convergence rate: • where
References • 徐树方,高立,张平文, • 数值线性代数,北京大学出版社,北京,2007 • 袁亚湘,孙文瑜, • 最优化理论与方法,科学出版社,北京,2000 • Yousef Saad, • Iterative Methods for Sparse Linear Systems, 2000
Home & Project Due at the end of this week Problem: Minimize the functional E(u)=∫(|u|2+u2-2fu )dx The corresponding Euler-Lagrange equation is E/u=-2u+2u-2f=0 or -u+u=f • Solve the following linear systems using CG method -u xx + u =f 0<x<1 f=(1+42)sin2x u(0)=u(1)=0 where • Set n = 100, 200, 300, 400, 500 • Use Matlab to graph the solution (j, uj)
Home & Project Due at the end of this week • Solve the following linear systems using CG method • The unknowns can be ordered as below
The coefficient matrix is of the block tridiagonal form Where S Tridiagonal matrix with diagonal entry: other entry: • Set n=20,40,80,100. Find the solution • Use Matlab to graph the solution (i, j, uij) -u+u=f (x,y)(0,1)(0,1) u(x,y)=100(x2-x)(y2-y) f=200(y-y2) + 200(x-x2) + 100(x2-x)(y2-y)