## Lower and Upper Triangular Systems.

A lower triangular matrix is a square matrix in which the elements to the right of the diagonal are all zero. An upper triangular matrix is a square matrix in which the elements to the left of the diagonal are all zero. A unit lower triangular matrix is a lower triangular matrix in which the diagonal elements are all ones. A unit upper triangular matrix is an upper triangular matrix in which the diagonal elements are all ones. The determinant of an upper or lower triangular matrix is simply the product of its diagonal elements. In particular, the determinant of a unit upper or lower triangular matrix is 1.A lower triangular system of equations,

**L x**=

**B**where

**L**= (

*l*

_{i,j}) is a nonsingular

*n × n*lower triangular matrix and

**B**= (

*b*

_{i}) is an

*n*-dimensional vector can be solved for the

*n*-dimensional vector

**x**= (

*x*

_{i}) by forward substitution:

*x*

_{1}=

*b*

_{1}/

*l*

_{1,1}

*x*

_{i}= [

*b*

_{i}- (

*l*

_{i,1}

*x*

_{1}+ ··· +

_{i-1})] /

*l*

_{i,i},

*i*= 2, . . . ,

*n*.

The algorithm for the inverse of the nonsingular

*n × n*lower triangular matrix

**L**can be obtained by expressing the product

**LL**=

^{-1}**I**in terms of the components of

**L**= (

*l*

_{ij}) and

**L**= (

^{-1}*l*

^{ ij}):

For

*i*= 1, . . . ,

*n*

*l*

^{ i,i}= 1 /

*l*

_{ i,i}

*j*=

*i*+1, . . . ,

*n*,

*l*

^{ i,j}= - (

*l*

_{ i,j}

*l*

^{ j,j}+ · · · +

*l*

_{i,i-1}

*l*

^{ i-1,j}) /

*l*

_{ i,i}

Conversely, an upper triangular system of equations,

**U x**=

**B**where

**U**= (

*u*

_{i,j}) is a nonsingular

*n × n*upper triangular matrix and

**B**= (

*b*

_{i}) is an

*n*-dimensional vector can be solved for the

*n*-dimensional vector

**x**= (

*x*

_{i}) by back substitution:

*x*

_{n}=

*b*

_{n}/

*u*

_{n,n}

*x*

_{i}= [

*b*

_{i}- (

*u*

_{i,i+1}

*x*

_{i+1}+ ... +

*u*

_{i,n}

*x*

_{n})] /

*u*

_{i,i},

*i*=

*n*-1, . . . , 1.

The algorithm for the inverse of the nonsingular

*n × n*upper triangular matrix

**U**can be obtained by expressing the product

**UU**=

^{-1}**I**in terms of the components of

**U**= (

*u*

_{ij}) and

**U**= (

^{-1}*u*

^{ ij}):

For

*i*=

*n*, . . . ,1

*u*

^{ i,i}= 1 /

*u*

_{ i,i}

*j*=

*n*, . . . ,

*i*+1,

*u*

^{ i,j}= - (

*u*

_{ i,i+1}

*u*

^{ i+1,j}+ · · · +

*u*

_{i,j}

*u*

^{ j,j}) /

*u*

_{ i,i}

Two versions of source code are provided for solving lower, upper, unit lower, and unit upper triangular linear systems of equations and two versions of source code are provided for calculating the inverse of lower, upper, unit lower, and unit upper triangular matrices. The different versions reflect alternative storage methods for triangular matrices.

### Function List

- int Lower_Triangular_Solve( double *L, double *B, double x[ ], int n )

Given an*n*×*n*nonsingular lower triangular matrix*L*, the function Lower_Triangular_Solve solves the linear equation*L x*=*B*given the*n*-dimensional vector*B*for the*n*-dimensional vector*x*. The function returns 0 if successful and -1 if the matrix*L*is singular.

- int Lower_Triangular_Inverse( double *L, int n )

Given an*n*×*n*nonsingular lower triangular matrix*L*, the function Lower_Triangular_Inverse calculates the inverse of*L*, the inverse overwrites the lower triangular part of*L*, the superdiagonal part of the matrix*L*is untouched. The function returns 0 if successful and -1 if the matrix*L*is singular.

- int Lower_Triangular_Solve_lt( double *L, double *B, double x[ ], int n )

Given an*n×n*lower triangular matrix*L*, one need only store the lower triangular part of*L*where the upper triangular part of*L*is assumed to be zero. For this routine it is assumed that only the lower triangular part of*L*is stored, i.e.*L*points to a memory location which contains*L[0][0]*immediately followed in memory by*L[1][0], L[1][1], … ,L[n-1][0], L[n-1][1], … ,L[n-1][n-1]*. Given an*n*×*n*nonsingular lower triangular matrix*L*, the function Lower_Triangular_Solve_lt solves the linear equation*L x*=*B*given the*n*-dimensional vector*B*for the*n*-dimensional vector*x*. The function returns 0 if successful and -1 if the matrix*L*is singular.

- int Lower_Triangular_Inverse_lt( double *L, int n )

Given an*n×n*lower triangular matrix*L*, one need only store the lower triangular part of*L*where the upper triangular part of*L*is assumed to be zero. For this routine it is assumed that only the lower triangular part of*L*is stored, i.e.*L*points to a memory location which contains*L[0][0]*immediately followed in memory by*L[1][0], L[1][1], … ,L[n-1][0], L[n-1][1], … ,L[n-1][n-1]*. Given an*n*×*n*nonsingular lower triangular matrix*L*, the function Lower_Triangular_Inverse_lt calculates the inverse of*L*, the inverse overwrites the lower triangular part of*L*. The function returns 0 if successful and -1 if the matrix*L*is singular.

- int Upper_Triangular_Solve( double *U, double *B, double x[ ], int n )

Given an*n*×*n*nonsingular upper triangular matrix*U*, the function Upper_Triangular_Solve solves the linear equation*U x*=*B*given the*n*-dimensional vector*B*for the*n*-dimensional vector*x*. The function returns 0 if successful and -1 if the matrix*U*is singular.

- int Upper_Triangular_Inverse( double *U, int n )

Given an*n*×*n*nonsingular upper triangular matrix*U*, the function Upper_Triangular_Inverse calculates the inverse of*U*, the inverse overwrites the upper triangular part of*U*, the subdiagonal part of the matrix*U*is untouched. The function returns 0 if successful and -1 if the matrix*U*is singular.

- int Upper_Triangular_Solve_ut( double *U, double *B, double x[ ], int n )

Given an*n×n*upper triangular matrix*U*, one need only store the upper triangular part of*U*where the lower triangular part of*U*is assumed to be zero. For this routine it is assumed that only the upper triangular part of*U*is stored, i.e.*U*points to a memory location which contains*U[0][0]*immediately followed in memory by*U[0][1], . . . , U[0][n-1], U[1][1], . . . . ,U[1][n-1], . . . ,U[n-1][n-1]*. Given an*n*×*n*nonsingular upper triangular matrix*U*, the function Upper_Triangular_Solve_ut solves the linear equation*U x*=*B*given the*n*-dimensional vector*B*for the*n*-dimensional vector*x*. he function returns 0 if successful and -1 if the matrix*U*is singular.

- int Upper_Triangular_Inverse_ut( double *U, int n )

Given an*n×n*upper triangular matrix*U*, one need only store the upper triangular part of*U*where the lower triangular part of*U*is assumed to be zero. For this routine it is assumed that only the upper triangular part of*U*is stored, i.e.*U*points to a memory location which contains*U[0][0]*immediately followed in memory by*U[0][1], . . . , U[0][n-1], U[1][1], . . . . ,U[1][n-1], . . . ,U[n-1][n-1]*. Given an*n*×*n*nonsingular upper triangular matrix*U*, the function Upper_Triangular_Inverse_ut calculates the inverse of*U*, the inverse overwrites the upper triangular part of*U*. The function returns 0 if successful and -1 if the matrix*U*is singular.

- void Unit_Lower_Triangular_Solve( double *L, double *B, double x[ ], int n )

Given an*n*×*n*unit lower triangular matrix*L*, the function Unit_Lower_Triangular_Solve solves the linear equation*L x*=*B*given the*n*-dimensional vector*B*for the*n*-dimensional vector*x*.

- void Unit_Lower_Triangular_Inverse( double *L, int n )

Given an*n*×*n*unit lower triangular matrix*L*, the function Unit_Lower_Triangular_Inverse calculates the inverse of*L*, the inverse overwrites the lower triangular part of*L*, the diagonal and superdiagonal part of the matrix*L*is untouched.

- void Unit_Lower_Triangular_Solve_lt( double *L, double *B, double x[ ], int n )

Given an*n × n*unit lower triangular matrix*L*, one need only store the lower triangular part of*L*where the upper triangular part of*L*is assumed to be zero. For this set of routines it is assumed that only the lower triangular part of*L*is stored, i.e.*L*points to a memory location which contains*L[0][0]*immediately followed in memory by*L[1][0], L[1][1], . . . . ,L[n-1][0], L[n-1][1], . . . ,L[n-1][n-1]*. Given an*n*×*n*unit lower triangular matrix*L*, the function Unit_Lower_Triangular_Solve_lt solves the linear equation*L x*=*B*given the*n*-dimensional vector*B*for the*n*-dimensional vector*x*.

- void Unit_Lower_Triangular_Inverse_lt(double *L, int n)

Given an*n × n*unit lower triangular matrix*L*, one need only store the lower triangular part of*L*where the upper triangular part of*L*is assumed to be zero. For this set of routines it is assumed that only the lower triangular part of*L*is stored, i.e.*L*points to a memory location which contains*L[0][0]*immediately followed in memory by*L[1][0], L[1][1], . . . . ,L[n-1][0], L[n-1][1], . . . ,L[n-1][n-1]*. Given an*n*×*n*unit lower triangular matrix*L*, the function Unit_Lower_Triangular_Inverse_lt calculates the inverse of*L*, the inverse overwrites the lower triangular part of*L*.

- void Unit_Upper_Triangular_Solve( double *U, double *B, double x[ ], int n )

Given an*n*×*n*unit upper triangular matrix*U*, the function Unit_Upper_Triangular_Solve solves the linear equation*U x*=*B*given the*n*-dimensional vector*B*for the*n*-dimensional vector*x*.

- void Unit_Upper_Triangular_Inverse( double *U, int n )

Given an*n*×*n*unit upper triangular matrix*U*, the function Unit_Upper_Triangular_Inverse calculates the inverse of*U*, the inverse overwrites the upper triangular part of*U*.

- void Unit_Upper_Triangular_Solve_ut( double *U, double *B, double x[ ], int n )

Given an*n × n*unit upper triangular matrix*U*, one need only store the upper triangular part of*U*where the lower triangular part of*U*is assumed to be zero. For this set of routines it is assumed that only the upper triangular part of*U*is stored, i.e.*A*points to a memory location which contains*U[0][0]*immediately followed in memory by*U[0][1], . . . , U[0][n-1], U[1][1], . . . . ,U[1][n-1], . . . ,U[n-1][n-1]*. Given an*n*×*n*unit upper triangular matrix*U*, the function Unit_Upper_Triangular_Solve_ut solves the linear equation*U x*=*B*given the*n*-dimensional vector*B*for the*n*-dimensional vector*x*.

- void Unit_Upper_Triangular_Inverse_ut( double *U, int n )

Given an*n × n*unit upper triangular matrix*U*, one need only store the upper triangular part of*U*where the lower triangular part of*U*is assumed to be zero. For this set of routines it is assumed that only the upper triangular part of*U*is stored, i.e.*A*points to a memory location which contains*U[0][0]*immediately followed in memory by*U[0][1], . . . , U[0][n-1], U[1][1], . . . . ,U[1][n-1], . . . ,U[n-1][n-1]*. Given an*n*×*n*unit upper triangular matrix*U*, the function Unit_Upper_Triangular_Inverse_ut calculates the inverse of*U*, the inverse overwrites the upper triangular part of*U*.

*C* Source Code

- The file, lower_triangular.c, contains the routines Lower_Triangular_Solve( ) and Lower_Triangular_Inverse( ).

- The file, lower_triangular_lt.c, contains the routines Lower_Triangular_Solve_lt( ) and Lower_Triangular_Inverse_lt( ).

- The file, upper_triangular.c, contains the routines Upper_Triangular_Solve( ) and Upper_Triangular_Inverse( ).

- The file, upper_triangular_ut.c, contains the routines Upper_Triangular_Solve_ut( ) and Upper_Triangular_Solve_ut( ).

- The file, unit_lower_triangular.c, contains the routines Unit_Lower_Triangular_Solve( ) and Unit_Lower_Triangular_Inverse( ).

- The file, unit_lower_triangular_lt.c, contains the routines Unit_Lower_Triangular_Solve_lt( ) and Unit_Lower_Triangular_Inverse_lt( ).

- The file, unit_upper_triangular.c, contains the routines Unit_Upper_Triangular_Solve( ) and Unit_Upper_Triangular_Inverse( ).

- The file, unit_upper_triangular_ut.c, contains the routines Unit_Upper_Triangular_Solve_ut( ) and Unit_Upper_Triangular_Inverse_ut( ).