Project 15 looks at the problem of solving a linear system, typically written A*x=b. You're probably familiar with a method like Gauss elimination for handling such problems. Gauss elimination works fine, so why should we look any further?
In your introductory classes on numerical methods, the matrix A was probably relatively small. If we let N stand for the number of columns of A, then the examples you looked at, and thought about, might have had values of N like 3, 4 or 5. If you used a Gauss elimination routine, you may have generated random linear systems with a value of N as high as 100, say.
The linear systems of interest to computational scientists, however, can easily have a value of N that is 10,000 or 1,000,000. With a matrix size this great, the methods appropriate for small systems break down.
In particular, the matrix may actually be too large to store. And even if we can store it, the time taken to apply Gauss elimination grows as N*N*N, so we can easily define problems whose solution would take hours, or days, or months!
Where do such huge matrices come from? Recall, for instance, that in Project 13, we discovered that we could approximate a partial differential equation using a grid, and that to get more accuracy, we needed to use a finer grid. Each time we double that particular grid, the order N of the associated matrix increases by a factor of 4.
What do such matrices look like? Unlike textbook examples, the huge matrices that occur in computational science are typically sparse, that is, a huge proportion of the entries of the matrix will be zero. This means that if we can think of a suitable data structure, we can save all the nonzero information in a huge matrix of order N, even if NxN, the size of the matrix, exceeds our available memory.
There are a number of methods for storing sparse matrices efficiently, and for solving the associated linear systems. Fortunately for us, MATLAB has these techniques built in, in such a way that we almost don't realize what is happening.
To start with, we will look at direct methods for solving the huge linear systems that are associated with the solution of the Poisson equation in 2D or 3D. Here, the word "direct" means that the method is not iterative. It carries out a series of steps with a definite termination at the end of which an (approximate) solution is produced.
In Project 20, we will look at the advantages of solving these kinds of systems using an iterative method instead of the direct approach employed here.
You can go up one level to the Computational Science Projects page.