# LUdecomposition -- LU decomposition

## Synopsis

• Usage:
(P,L,U) = LUdecomposition A
• Inputs:
• A, , or ,of size m by n, over a a finite field ZZ/p, RR or CC
• Outputs:
• P, a list, of integers in [0, ..., m-1] forming a permutation
• L, , or , an m by min(m,n) lower triangular matrix, with 1's on the diagonal
• U, , or , a min(m,n) by n upper triangular matrix

## Description

The output matrices are mutable exactly when the input matrix is, but the matrix A is not modified

If Q is the m by m permutation matrix such that Q_(P_i,i) = 1, and all other entries are zero, then A = QLU.

There are several restrictions. The first is that there are only a limited number of rings for which this function is implemented. Second, if A is a mutable matrix defined over RR or CC, then A must be densely encoded. This restriction is not present if A is .

 i1 : kk = ZZ/101 o1 = kk o1 : QuotientRing i2 : A = matrix"1,2,3,4;1,3,6,10;19,7,11,13" ** kk o2 = | 1 2 3 4 | | 1 3 6 10 | | 19 7 11 13 | 3 4 o2 : Matrix kk <--- kk i3 : (P,L,U) = LUdecomposition A o3 = ({0, 1, 2}, | 1 0 0 |, | 1 2 3 4 |) | 1 1 0 | | 0 1 3 6 | | 19 -31 1 | | 0 0 47 22 | o3 : Sequence i4 : Q = id_(kk^3) _ P o4 = | 1 0 0 | | 0 1 0 | | 0 0 1 | 3 3 o4 : Matrix kk <--- kk i5 : Q * L * U == matrix A o5 = true
For matrices over RR or CC, this function calls the lapack routines, which are restricted to 53 bits of precision.
 i6 : A = matrix"1,2,3,4,5,6;1,3,6,12,13,16;19,7,11,47,48,21" ** RR o6 = | 1 2 3 4 5 6 | | 1 3 6 12 13 16 | | 19 7 11 47 48 21 | 3 6 o6 : Matrix RR <--- RR 53 53 i7 : (P,L,U) = LUdecomposition A o7 = ({2, 1, 0}, | 1 0 0 |, | 19 7 11 47 48 | .0526316 1 0 | | 0 2.63158 5.42105 9.52632 10.4737 | .0526316 .62 1 | | 0 0 -.94 -4.38 -4.02 ------------------------------------------------------------------------ 21 |) 14.8947 | -4.34 | o7 : Sequence i8 : Q = id_ (RR^3) _ P o8 = | 0 0 1 | | 0 1 0 | | 1 0 0 | 3 3 o8 : Matrix RR <--- RR 53 53 i9 : Q * L * U - A o9 = | 0 -2.22045e-16 0 0 0 0 | | 0 0 0 0 0 0 | | 0 0 0 0 0 0 | 3 6 o9 : Matrix RR <--- RR 53 53 i10 : clean(1e-15,oo) o10 = | 0 0 0 0 0 0 | | 0 0 0 0 0 0 | | 0 0 0 0 0 0 | 3 6 o10 : Matrix RR <--- RR 53 53
Mutable matrices can sometimes be useful for speed, and/or space. If A is a mutable matrix, it must be densely encoded (which is the default).
 i11 : A = mutableMatrix(CC,5,10, Dense=>true) o11 = | 0 0 0 0 0 0 0 0 0 0 | | 0 0 0 0 0 0 0 0 0 0 | | 0 0 0 0 0 0 0 0 0 0 | | 0 0 0 0 0 0 0 0 0 0 | | 0 0 0 0 0 0 0 0 0 0 | o11 : MutableMatrix i12 : printingPrecision = 2 o12 = 2 i13 : setRandomSeed 0 o13 = 0 i14 : fillMatrix A o14 = | .89+.67ii .91+.31ii .35+.56ii .19+.4ii .28+.61ii .034+.51ii .44+.64ii .47+.4ii .48+.82ii .5+.15ii | | .29+.63ii .074+.81ii .25+.15ii .62+.015ii .97+.68ii .15+.66ii .19+.52ii .16+.71ii .98+.06ii .47+.77ii | | .026+.71ii .36+.71ii .83+.54ii .22+.39ii .91+.89ii .17+.63ii .99+.57ii .91+.57ii .065+.28ii .31+.54ii | | .89+.23ii .13+.25ii .87+.42ii .56+.87ii .17+.97ii .35+.38ii .18+.37ii .31+.73ii .98+.21ii .21+.28ii | | .46+.78ii .74+.11ii .61+.85ii .7+.68ii .065+.88ii .24+.12ii .34+.062ii .56+.63ii .59+.83ii .096+.61ii | o14 : MutableMatrix i15 : (P,L,U) = LUdecomposition A; i16 : Q = id_(CC^5) _ P o16 = | 1 0 0 0 0 | | 0 0 0 0 1 | | 0 0 1 0 0 | | 0 1 0 0 0 | | 0 0 0 1 0 | 5 5 o16 : Matrix CC <--- CC 53 53 i17 : matrix Q * matrix L * matrix U - matrix A o17 = | 0 0 0 0 0 | 5.6e-17 -5.6e-17 0 -1.1e-16+1.7e-16ii -2.2e-16 | 5.6e-17 0 0 -2.8e-17 0 | -1.1e-16ii 0 0 0 0 | -1.1e-16ii 0 0 0 -2.8e-17+1.1e-16ii ----------------------------------------------------------------------- 0 0 0 0 0 | -2.8e-17 0 0 0 1.1e-16ii | 0 1.1e-16 1.1e-16 0 0 | 0 0 0 0 0 | -2.8e-17-5.6e-17ii 0 -1.1e-16ii 0 0 | 5 10 o17 : Matrix CC <--- CC 53 53 i18 : clean(1e-15,oo) o18 = | 0 0 0 0 0 0 0 0 0 0 | | 0 0 0 0 0 0 0 0 0 0 | | 0 0 0 0 0 0 0 0 0 0 | | 0 0 0 0 0 0 0 0 0 0 | | 0 0 0 0 0 0 0 0 0 0 | 5 10 o18 : Matrix CC <--- CC 53 53

## Caveat

This function is limited in scope, but is sometimes useful for very large matrices

• solve -- solve linear equation(s)
• SVD -- singular value decomposition of a matrix
• MutableMatrix -- the class of all mutable matrices
• fillMatrix -- fill a mutable matrix with random numbers
• clean -- Set to zero elements that are approximately zero
• norm

## Ways to use LUdecomposition :

• "LUdecomposition(Matrix)"
• "LUdecomposition(MutableMatrix)"

## For the programmer

The object LUdecomposition is .