National Center For Theoretical Sciences Mathematics Division, Taiwan High-Performance Numerical Solvers Dec. 2017 - Jan. 2018 Assignment 1

Download as pdf or txt
Download as pdf or txt
You are on page 1of 1

National Center for Theoretical Sciences

Mathematics Division, Taiwan

High-Performance Numerical Solvers


Dec. 2017 – Jan. 2018

Assignment 1

1. For solving Ax = b, write Matlab functions for:


(a) Classical conjugate gradient (CG) method
(b) CG method of Chronopoulos and Gear
(c) CG method of Ghysels and Vanroose (pipelined CG)
which were discussed in the lecture on Dec. 22. Implement each method such that the
known solution is one of the input parameters, and the method returns a vector containing
the A-norm of the error at each iteration.
2. Compare the three CG methods using some SPD matrices from the Suite Sparse Matrix
Collection (https://sparse.tamu.edu). For a matrix A, choose the known solution sol
as
sol = rand(length(A),1)-.5;
and then construct the right-hand side vector b as b=A*sol;. Be sure to use the same
right-hand side vector when comparing the three methods. At each step i of the CG
methods, compute and store the A-norm of the error:
error_norm(i) = (x-sol)’*A*(x-sol);
Also store the residual norm at each step. Plot the error norm and the residual norm as a
function of step number for several example matrices. Plot the error norms for the three
methods on the same graph; separately plot the residual norms for the three methods
on the same graph. Be sure to label the different methods on each graph. Try to find
matrices in which the behavior among the three methods is different. You may need to
choose ill-conditioned matrices, e.g., those with condition number around 109 . What can
you conclude about the accuracy of the methods compared to each other?
3. Experiment with the idea of residual replacement, i.e., to try to improve the stability of
the methods, occasionally (at some interval) replace the residual that was computed via
recurrences, by the residual computed explicitly by b − Axk . Show how this changes your
results above.
4. Prove the identity used to derive the Chronopoulos and Gear CG method:
βk−1
(Apk , pk ) = (Ark , rk ) − (rk , rk ).
αk−1
There may be many ways to prove this, but, as a hint, the following may or may not be
useful: (i) αk can also be written as (rk , rk )/(Apk , rk ); (ii) Start with the formula in CG
for updating pk .

You might also like