How to put the finite difference in a a matrix to solve the linear system

Question

After doing the finite difference approximation of a pde equation or ode, we have a linear equation of the fnite difference . how do we put this equation into a matrix to solve it? I want the actual steps of putting the numbers in the rows of the matrix because i don't find anything explains this!! How do we put the coefficients in the matrix? Is it like the usual linear equations? I how does the linear system become sparse?

Answer 1

The answer of your questions of course depends on the problem you are discretising and the particular finite difference stencil you use, however I'll try to illustrate with an example:

Say that we want to solve the ODE \begin{align} \frac{du}{dt} + \lambda u &= 0, \quad 0 < t \leq 1 \\ u &= f, \quad t=0. \end{align} We divide time into $N$ intervals of length $\Delta t$ and denote the numerical solution at time step $j$ as $v^j$. In other words, time is now discretised into points $t_j$, $j=0, \dots, N$ and the numerical solution is $v^j \approx u(t_j)$.

Say that we want to approximate the time derivative using a second order accurate central stencil: $$ \frac{du}{dt} \Big|_{t=t_j} \approx \frac{v_{j+1} - v_{j-1}}{2 \Delta t}. $$ Thus, we replace the ODE above with the $N-1$ equations $$ \frac{v_{j+1} - v_{j-1}}{2 \Delta t} + \lambda v_j = 0, \quad j=1, \dots N-1. $$ I have excluded $j=0$ and $j=N$ here because the central stencil cannot be applied at these points. Nonetheless, we now have a set of linear equations involving the unknowns $v_j$ and $v_{j \pm 1}$. On matrix form this becomes $$ \begin{pmatrix} \ddots & & & & & & \\ & \ddots & & & & & \\ & -\frac{1}{2\Delta t} & \lambda & \frac{1}{2\Delta t} & & & \\ & & -\frac{1}{2\Delta t} & \lambda & \frac{1}{2\Delta t} & & \\ & & & -\frac{1}{2\Delta t} & \lambda & \frac{1}{2\Delta t} & \\ & & & & & \ddots & \\ & & & & & & \ddots \end{pmatrix} % \begin{pmatrix} \vdots \\ v_{j-2} \\ v_{j-1} \\ v_j \\ v_{j+1} \\ v_{j+2} \\ \vdots \end{pmatrix} = \begin{pmatrix} \vdots \\ \vdots \\ 0 \\ 0 \\ 0 \\ \vdots \\ \vdots \end{pmatrix}. $$ The general procedure here is that everything that multiplies $v_j$ (in this case $\lambda$) ends up on the diagonal, everything that multiplies $v_{j-1}$ ends up on the first sub-diagonal and everything that multiplies $v_{j+1}$ ends up on the first super-diagonal. Everything else is zero, hence the matrix is sparse. The number of rows of the matrix corresponds to the number of points in the discretisation, in this case $N+1$. Of course, the stencil would have to be modified on the first and last rows since the central stencil simple does not fit.