MA/CS/EGR 537: Assignment 3

This is due Thursday, October 2 in class. Please e-mail your Matlab script(s) to

douglas@ccs.uky.edu

Theory

Assume that A is real, symmetric, positive definite, and N by N. Use the definition of the conjugate gradient method for solving Ax+b=0 given in class:

1.	x0 arbitrary	(approximate solution)
	r0 = Ax0 + b	(approximate residual)
	w0 = r0	(search direction)
2.	For k=0,1,...
		alphak = -(rk,wk)/(wk,Awk)
		xk+1 = xk + alphak wk
		rk+1 = rk + alphak Awk
		betak = -(rk+1,Awk)/(wk,Awk)
		wk+1 = rk+1 + betak wk
	End

Rigorously prove the results below, where indicated.

Lemma 2: The parameter beta_k is chosen so that (w_k+1,Aw_k)=0.

Proof: Problem 1.

Lemma 3: For 0<=q<=k,

(a)	(rk+1,wq)=0
(b)	(rk+1,rq)=0
(c)	(wk+1,Awq)=0

Proof: Problem 2.

Lemma 4: alpha_k = (r_k+1,r_k+1) / (r_k+1,Ar_k).

Proof: Problem 3.

Lemma 5: beta_k = (r_k+1,r_k+1) / (r_k,r_k).

Proof: Problem 4.

Computation

Assume that M=n², where n>1. For your calculations, use n=4. Use a matrix A that is block tridiagonal. Both block matrices are n by n. Let

A = [ -I, T, -I ], where T = [ -1, 4, -1 ] and I = identity matrix.

I = eye(n);
T = diag(ones(n-1,1),1);
T = 2 * I - ( T + T' );
A = kron(T,I) + kron(I,T);

M = n * n;
rand('seed', number );
b = rand(M,1);

You should run the CG algorithm until the residual norm has been reduced by a factor of 10^-6, i.e., run until

|| r_k || / || r₀ || <= 10^-6

Cheers,
Craig C. Douglas