CME 211: Project Part 1
This project was designed by Patrick LeGresley and modified for the purposes of this course.
For the final project you will be developing a program to solve the 2D heat equation on a simple geometry
using a sparse matrix solver written in C++. As a first step you will be implementing a sparse matrix solver
in C++. In Part 2, you will be forming the system of equations for a specified geometry, solving the system,
and writing Python code to visualize your results.
When solving partial dierential equations, a commonly used method to solve the equation numerically is the
finite-dierence method. In this technique, the solution domain of interest is discretized into a grid of points
and derivative terms are replaced by appropriate finite dierence approximations, using stencils that involve
information from neighboring points. Consider the steady-state heat equation in two dimensions:
ˆy2 = 0.
To discretize this equation on a Cartesian grid, we would use the following finite-dierence approximations to
obtain the 2nd derivatives at any grid location (indexed with integers i and j):
x2 (ui≠1,j ≠ 2ui,j + ui+1,j )
y2 (ui,j≠1 ≠ 2ui,j + ui,j+1)
Substituting these approximations into the PDE and assuming constant grid spacing (x = y = h) gives
the following equation for each point (i, j) in the solution domain:
h2 (ui≠1,j + ui+1,j + ui,j≠1 + ui,j+1 ≠ 4ui,j )=0
Assembling the equations for the grid points in your domain will lead to a linear system of equations Au = b.
Due to the compact nature of the finite dierence stencils used in the discretization, systems of equations
derived using FD are generally quite sparse. Additionally, for certain geometries (like the one you will be
dealing with in Part 2), the resulting system will be symmetric positive/negative definite. To solve systems
with this structure, an ecient iterative algorithm to use is the Conjugate Gradient (CG) method. The
pseudo-code for the CG algorithm is given here:
For this first part of the project your task is to implement the CG algorithm in C++ as a function within
the file CGSolver.cpp. A prototype of this function can be found in the provided CGSolver.hpp header
file. You should use CSR as the sparse matrix format in your solver, as this format is better suited to the
matrix-vector operations that arise in CG.
However, building a matrix in COO format is much more straightforward and that is the format you will use
in the second part of the project. To assist you with the required conversion, a header file COO2CSR.hpp has
ft’T 1 1 start _o
colo j in 0 I
r0 = b – A u0
L2normr0 = L2norm(r0)
p0 = r0
niter = 0
while (niter < nitermax)
niter = niter + 1
alpha = (rn
T A pn)
un+1= un + alphan pn
rn+1= rn – alphan A pn
L2normr = L2norm(rn+1)
if L2normr/L2normr0 < threshold
betan = (rn+1
pn+1 = rn+1 + betan pn
For$ this$ first part$ of$ the$ project your$ task$ is$ to$ implement$ the$ CG$ algorithm$ in$ C++ as$ a$
function$within$ the file$CGSolver.cpp.$A$ prototype$ of$ this$ function$ can$ be$ found$in$ the$
provided$CGSolver.hpp header$ file.$$You$should use$CSR$as$ the$sparse$matrix$ format$in$
your$solver,$as$this$format$is$better$suited$to$the matrixSvector$operations that$arise$in$CG.
However,$ building$ a$ matrix$ in$ COO$ format$ is$ much$ more$ straightforward$ and$ that$ is$ the$
format$ you$ will$ use$ in$ the$ second$ part$ of$ the$ project.$ $ To$ assist$ you$ with$ the$ required$
conversion,$ a$ header$ file$ COO2CSR.hpp has$ been$ provided$ for$ you$ which$ contains$ a$
function$COO2CSR() that$can$be$used$to$convert$an$existing$COO$matrix to$CSR$format$inS
place$ (i.e. the$ provided$ function$ will$ overwrite$ your$ input vectors$ and$ fill$ them$ in$ with$
When$ implementing$ your$ solver$ you$ should develop additional$ functions$ to$ perform$
common$ vector$ and$ matrix$ operations$ that$ occur in$ the$ CG$ algorithm.$ $ These$ functions$
should$ be$ implemented$ in$ a$ file$ matvecops.cpp with$ a$ corresponding$ include$ file$
To$ test$ and$ debug$ your$ solver implementation$ two test$ matrices in$ COO$ format$ are
provided in$the$files$matrix1.txt and$matrix2.txt.$$In$each$file$the$first$line contains
To$ demonstrate$ the$ use$ of$ your$ solver$write$a$main() in$main.cpp that$loads a$matrix$
Figure 1: Conjugate Gradient pseudo-code
y yf fx wit l t I A value 24 I I 33 sure
b double zero setarat.IE 5AA E Y403
E o o y g
Mxn matrix I IX
matrix1 txt l e 5 to
I I I
a d iL II
zf f 4
nonzero in row 1
C a S
been provided for you which contains a function COO2CSR() that can be used to convert an existing COO
matrix to CSR format in-place (i.e. the provided function will overwrite your input vectors and fill them in
with appropriate data for the CSR format).
When implementing your solver you should develop additional functions to perform common vector and matrix
operations that occur in the CG algorithm. These functions should be implemented in a file matvecops.cpp
with a corresponding include file matvecops.hpp for the prototypes.
To test and debug your solver implementation two test matrices in COO format are provided in the files
matrix1.txt and matrix2.txt. In each file the first line contains the number of rows and columns in the
matrix. For subsequent lines, the values are the row index, column index, and the value of the matrix entry.
To demonstrate the use of your solver write a main() in main.cpp that loads a matrix from a file specified at
the command line, converts the matrix to CSR format, runs your CG solver function with a starting guess of
ones for the solution and zeros for the right hand side, and writes the solution to the specified file name. Use
a tolerance of 1.e-5 and if your solver is running properly, the solution should be a vector of zeros:
$ ./main matrix2.txt solution2.txt
SUCCESS: CG solver converged in 25 iterations.
Adherence to Prototypes
We provide you with header files, and we ask that you use them as a starting place. You are encouraged to
consider adding use of keyword const to arguments when appropriate.
Preliminary write up
In LATEX, create a writeup.tex that uses the algorithm environment to provide the pseudo-code for the CG
algorithm. Also write a short (no more than one page) discussion of how you implemented your CG solver in
terms of functions to eliminate redundant code.
Summary of requirements
1. Create a directory called project in your CME 211 homework repository. All source files and
documentation will go in this directory.
2. Write a function in CGSolver.cpp that solves a linear system in CSR format using the CG method.
The implementation of this function must be consistent with the provided prototype in CGSolver.hpp.
3. For common operations in your solver algorithm, develop vector and matrix functions in matvecops.cpp
and create the associated header file matvecops.hpp with the prototypes.
4. Write a main() in main.cpp that loads a matrix from a file in COO format, converts the matrix to
CSR format using the provided COO2CSR() function, and solves the system using your CGSolver()
function. Write the solution vector to the specified solution file: one value per line, scientific notation, 4
5. Write a makefile to compile all source code and produce the main executable.
• $ make: must compile all source code and produce main
• $ make clean: must remove all object (*.o) files and remove main
Command line interface
Your build system must produce an executable named main that operates according to the following examples:
A µ seavectorcanine
• provide a usage message if there are no command line arguments:
./main <input matrix file name> <output solution file name>
• with files specified:
$ ./main matrix2.txt solution2.txt
SUCCESS: CG solver converged in 25 iterations.
Please be very careful with directory and file names. To be clear you should have a minimum of these files in
the project directory in your GitHub repository:
Do not commit:
• temporary files produced by your text editor
• object files produced by compiler (*.o files)
• the binary executable (main)
Whatever files you have in your GitHub repository at the deadline will be considered your final submission.
I I S
cosowercdata.it if if I e 5
A fly b X
I b I x fi
r f I Po III
L 2 norm to 25 Jf
alpha ET 4T
Iff 120 168
re III alpha t
I I Est III I L
co 2 CSR app Te coozcSR hpp
CG solver app maturecops cpp
main.o mainapp CGSolver.hppcoozcsr.tn pco02CSR.o c002CSR.cppCGSolver o i
CGSolver.cppmatvecops.hppmatuecops.co i matrecops cpp