Solution Methods in Numerical Linear Algebra

Abstract

This thesis investigated the influence of matrix structure, arising from various boundary conditions, on the performance of both direct and iterative solvers for linear systems resulting from discretized partial differential equations (PDEs). Using the negative one-dimensional Poisson equation as a model problem, system matrices were constructed for Dirichlet, Neumann, and mixed boundary conditions, including symmetrized variants.

The study classified the resulting matrices in terms of singularity and positive definiteness. Direct solvers such as LU decomposition and Cholesky factorization were tested for accuracy and applicability. Iterative solvers (Jacobi, Gauss-Seidel, and Successive Over-Relaxation or SOR) were analysed using the spectral radius of their iteration matrices to predict convergence.

Key findings include:
- Matrices from Dirichlet and mixed boundary conditions were non-singular; Neumann matrices were singular but still allowed convergence in iterative solvers under specific conditions.
- Cholesky factorization only applied to positive definite systems, while LU decomposition was more generally applicable to non-singular systems.
- Gauss-Seidel converged faster than Jacobi, with twice the convergence speed in tridiagonal cases.
- SOR achieved the fastest convergence when using the optimal relaxation parameter, though this optimal ω is generally not known a priori.
- Even when the spectral radius ρ(B) approached 1, as in Neumann systems, convergence was observed due to the structure of the eigenvalue spectrum.
\end{itemize}

The thesis concludes that the structural properties of the matrix, such as singularity and diagonal dominance, are critical in determining the suitability and efficiency of numerical solvers. It also highlights the fragility of these methods when applied to nearly singular systems, motivating future research into more robust solution techniques and convergence criteria.