## Factoring block tridiagonal symmetric positive definite matrices

### Goal

Perform Cholesky factorization of a symmetric positive definite block tridiagonal matrix.

### Solution

To perform Cholesky factorization of a symmetric positive definite block tridiagonal matrix, with `N` square blocks of size `NB` by `NB`:

Perform Cholesky factorization of the first diagonal block.

Repeat

`N`- 1 times moving down along the diagonal:Compute the off-diagonal block of the triangular factor.

Update the diagonal block with newly computed off-diagonal block.

Perform Cholesky factorization of a diagonal block.

Source code: see the BlockTDS_SPD/source/dpbltrf.f file in the samples archive available at https://software.intel.com/content/dam/develop/external/us/en/documents/mkl-cookbook-samples-120115.zip.

### Cholesky factorization of a symmetric positive definite block tridiagonal matrix

```
…
CALL DPOTRF('L', NB, D, LDD, INFO)
…
DO K=1,N-1
CALL DTRSM('R', 'L', 'T', 'N', NB, NB, 1D0, D(1,(K-1)*NB+1), LDD, B(1,(K-1)*NB+1), LDB)
CALL DSYRK('L', 'N', NB, NB, -1D0, B(1,(K-1)*NB+1), LDB, 1D0, D(1,K*NB+1), LDD)
CALL DPOTRF('L', NB, D(1,K*NB+1), LDD, INFO)
…
END DO
```

### Routines Used

Task |
Routine |
Description |
---|---|---|

Perform Cholesky factorization of diagonal blocks |
DPOTRF |
Computes the Cholesky factorization of a symmetric (Hermitian) positive-definite matrix. |

Compute off-diagonal blocks of the triangular factor |
DTRSM |
Solves a triangular matrix equation. |

Update the diagonal blocks |
DSYRK |
Performs a symmetric rank-k update. |

### Discussion

A symmetric positive definite block tridiagonal matrix, with `N` diagonal blocks `D`_{i} and `N` - 1 sub-diagonal blocks `B`_{i} of size `NB` by `NB` is factored as:

Multiplying the blocks of the matrices on the right gives:

Equating the elements of the original block tridiagonal matrix to the elements of the multiplied factors yields:

Solving for `C`_{i} and `L`_{i}`L`_{i}^{T}:

Note that the right-hand sides of the equations for `L`_{i}`L`_{i}^{T} is a Cholesky factorization. Therefore a routine chol() for performing Cholesky factorization can be applied to this problem using code such as:

L_{1}=chol(D_{1}) do i=1,N-1 C_{i}=B_{i}∙L_{i}^{-T}//trsm() D_{i + 1}:=D_{i + 1}- C_{i}∙C_{i}^{T}//syrk() L_{i + 1}=chol(D_{i + 1}) end do

**Parent topic:**Intel® oneAPI Math Kernel Library Cookbook