Computing Tridiagonal Forms

Introduction and description of problem

We begin with a few basic definitions of terms, one of them being informal.
Given a square matrix A, denote its (i,j) element as aij, for all i,j. We say that A is a lower unit triangular matrix, if aii=1 for all i and aij=0 for all i,j with i>j. Furthermore, we say that A is a tridiagonal matrix, if aij=0 for all i,j with |i-j|>1.
Informally, we say that A is sparse, if it has just a few non-zero elements. The term "just a few" is what we don't define here formally. One can treat it as some constant part of the elements (e.g. 1/10). The other can treat it as 1/n of them, where n is the size of the matrix A.
There are reasons for sparse matrices being handful. One advantage they offer is relatively simple mathematical calculations. For instance, it is common that the sum of two sparse matrices is sparse. Additionally, in the scope of computing, there are ways of storing sparse matrices in memory efficiently (see section 2.1 in [1]).

As for the main research topic, we were seeking for an efficient algorithm that, given a symmetric sparse matrix A, partitions it into a form A=LTLt, where L is a lower unit triangular sparse matrix and T is a tridiagonal matrix. That is, if such a partitioning exists. Ignoring the requirement that L is sparse, there is already an algorithm which solves this named after Parlett-Raid, along with a recursive version[2].

Progress and results

During the quest of solving, we were experimenting with implementations of a few mathematical functions, the main candidate being a recursive Parlett-Reid algorithm[2]. Writing A=LTLt as before, during the execution of a recursive step in the algorithm we are basically filling the unknown elements of a sub-matrix of L column-by-column. It is common that there is one and unique column to fill in at a certain stage, while there are also times we may not continue and have to halt (e.g. if there's no partitioning). Similarly, at a certain stage we may have a freedom of choice for a specific column. This is where it can be wise to try and set it to a zero column. After all, our goal is to have the matrix L in the partitioning as sparse as we can.
What was not found is an efficient algorithm proven to always return a partitioning A = LTLt with the mostly sparse L (assuming one exists). What can be shown to be found as seen above is a way of playing with a known algorithm in order to try and get a sparse L.

Attached file

Project document

1. Timothy A. Davis. Direct Methods for Sparse Linear Systems. Published: Philadelphia : SIAM, Society for Industrial and Applied Mathematics, c2006
2. Gil Shklarski and Sivan Toledo. Blocked and Recursive Algorithms for Triangular Tridiagonalization. Accepted to ACM Transactions on Mathematical Software, Nov 2009.