IDL

LA_LUDC

LA_LUDC

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:

LAPACK Routine Basis for LA_LUDC

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.

See Also


LA_LUMPROVE, LA_LUSOL, LUDC



Notes


This page has no user notes yet. Be the first one!


This information is not subject to the controls of the International Traffic in Arms Regulations (ITAR) or the Export Administration Regulations (EAR). However, it may be restricted from transfer to various embargoed countries under U.S. laws and regulations.
© 2014 Exelis Visual Information Solutions