Chapter 3 Numerical linear algebra
Up to this point, our investigations have been mostly theoretical. We will now turn to some important questions in numerical linear algebra.
Throughout our study of linear algebra, two computational issues have repeatedly presented themselves:
-
Solve the equation \(A\xvec=\bvec\text{.}\)
-
Find the eigenvalues and eigenvectors of a matrix.
In fact, these issues are important throughout a wide range of scientific applications. For instance, studying a lightning strike on an airplane involves solving \(A\xvec=\bvec\) where \(A\) is a very large matrix. This problem is actually modeled by a partial differential equation that proves to be too difficult to solve directly. Instead, we βdiscretizeβ the equation, replacing the derivatives with difference quotients on a discrete mesh of points. In this way, the partial differential equation is approximated by a linear system \(A\xvec=\bvec\) where there is one equation for each point in our mesh.
Just as in your first calculus course, the finer the mesh, the closer the points are to one another and the better the approximation. That is, the more points we have, the bigger the matrix and the more confidence we have in our results. Throughout this chapter, we will usually think of an \(m\times
m\) matrix \(A\text{,}\) but we actually want to think of \(m\) as being infinite. In practice, this could mean that \(m=100\) billion.
As we saw in our earlier course, the effort involved in solving \(A\xvec = \bvec\) using Gaussian elimination is proportional to \(m^3\text{.}\) This means that solving a system where \(m=100\) is roughly 1000 times more work than when \(m=10\text{.}\) So if \(m=100\) billion and want to know a solution in our lifetime, Gaussian elimination is completely unrealistic.
Fortunately, our matrices \(A\) typically have some kind of structure that we can exploit. For instance, \(A\) could be sparse, meaning that most of its entries are zero. This often happens when approximating the solution to a partial differential equation as mesh points only interact with nearby mesh points. For instance, if \(A\) is approximating a two-dimensional differential operator on a discrete rectangular mesh, each row may only have four entries since each point interacts with four adjacent points, two horizontally and two vertically.
In fact, in this case, we typically donβt want to actually store the matrix \(A\) in a computer because the matrix is huge and most of the entries are zero. Instead, we just want to use the fact that we can compute \(A\xvec\) for any vector \(\xvec\text{,}\) which is certainly the case when \(A\) is approximating a differential operator.