240 likes | 407 Views
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.
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)