The LA_LUDC procedure computes the LU decomposition of an n-column by m-row array as:

A = P L U

where P is a permutation matrix, L is lower trapezoidal with unit diagonal elements (lower triangular if n = m), and U is upper trapezoidal (upper triangular if n = m).

LA_LUDC is based on the following LAPACK routines:

 Output Type LAPACK Routine Float sgetrf Double dgetrf Complex cgetrf Double complex zgetrf

## Examples

The following example uses the LU decomposition on a given array, then determines the residual error of using the resulting lower and upper arrays to recompute the original array:

`PRO ExLA_LUDC`
`; Create a random array:`
`n = 20`
`seed = 12321`
`array = RANDOMN(seed, n, n, /RAN1)`
` `
`; Compute LU decomposition.`
`aludc = array           ; make a copy`
`LA_LUDC, aludc, index`
` `
`; Extract the lower and upper triangular arrays.`
`l = IDENTITY(n)`
`u = FLTARR(n, n)`
`FOR j = 1,n - 1 DO l[0:j-1,j] = aludc[0:j-1,j]`
`FOR j=0,n - 1 DO u[j:*,j] = aludc[j:*,j]`
` `
`; Reconstruct array, but with rows permuted.`
`arecon = l ## u`
`; Adjust from LAPACK back to IDL indexing.`
`Index = Index - 1`
`; Permute the array rows back into correct order.`
`; Note that we need to loop in reverse order.`
`FOR i = n - 1,0,-1 DO BEGIN & \$`
`   temp = arecon[*,i]`
`   arecon[*, i] = arecon[*,index[i]]`
`   arecon[*, index[i]] = temp`
`ENDFOR`
`PRINT, 'LA_LUDC Error:', MAX(ABS(arecon - array))`
`END`
` `
`ExLA_LUDC`

When this program is compiled and run, IDL prints:

`LA_LUDC error: 4.76837e-007`

## Syntax

LA_LUDC, Array, Index [, /DOUBLE] [, STATUS=variable]

## Arguments

### Array

A named variable containing the real or complex array to decompose. This procedure returns Array as its LU decomposition.

### Index

An output vector with MIN(m, n) elements that records the row permutations which occurred as a result of partial pivoting. For 1 < j < MIN(m,n), row j of the matrix was interchanged with row Index[j].

Note: Row numbers within Index start at one rather than zero.

## Keywords

### DOUBLE

Set this keyword to use double-precision for computations and to return a double-precision (real or complex) result. Set DOUBLE = 0 to use single-precision for computations and to return a single-precision (real or complex) result. The default is /DOUBLE if Array is double precision, otherwise the default is DOUBLE = 0.

### STATUS

Set this keyword to a named variable that will contain the status of the computation. Possible values are:

• STATUS = 0: The computation was successful.
• STATUS > 0: One of the diagonal elements of U is zero. The STATUS value specifies which value along the diagonal (starting at one) is zero.

Note: If STATUS is not specified, any error messages will output to the screen.

## Version History

 5.6 Introduced

## Resources and References

For details see Anderson et al., LAPACK Users' Guide, 3rd ed., SIAM, 1999.