>  Docs Center  >  IDL Reference  >  Advanced Math and Stats  >  IMSL_SP_LUSOL

IMSL_SP_LUSOL

IMSL_SP_LUSOL

The IMSL_SP_LUSOL function solves a sparse system of linear equations Ax = b. By using keywords, any of several related computations can be performed.

This routine requires an IDL Advanced Math and Stats license. For more information, contact your sales or technical support representative.

The IMSL_SP_LUSOL function solves a system of linear equations Ax = b, where A is sparse. In its default usage, it solves the so-called one off problem, by first performing an LU factorization of A using the improved generalized symmetric Markowitz pivoting scheme. The factor L is not stored explicitly because the saxpy operations performed during the elimination are extended to the right-hand side, along with any row interchanges. Thus, the system Ly = b is solved implicitly. The factor U is then passed to a triangular solver which computes the solution x from Ux = y.

If a sequence of systems Ax = b are to be solved where A is unchanged, it is usually more efficient to compute the factorization once, and perform multiple forward and back solves with the various right-hand sides. In this case the factor L is explicitly stored and a record of all row as well as column interchanges is made. The solve step then solves the two triangular systems Ly = b and Ux = y. In this case, you should first call IMSL_SP_LUFAC to compute the factorization, then use the keyword FACTOR_COORD with the IMSL_SP_LUSOL function.

If the solution to ATx = b is required, specify the keyword Transpose. This keyword only alters the forward elimination and back substitution so that the operations UTy = b and LTx = y are performed to obtain the solution. So, with one call to IMSL_SP_LUFAC to produce the factorization, solutions to both Ax = b and ATx = b can be obtained.

The keyword CONDITION is used to calculate and return an estimation of the L1 condition number of A. The algorithm used is due to Higham. Specifying CONDITION causes a complete L to be computed and stored, even if a one-off problem is being solved. This is due to the fact that Higham’s method requires a solution to problems of the form Az = r and ATz = b .

The default pivoting strategy is symmetric Markowitz (PIVOTING = 3). If a row or column oriented problem is encountered, there may be some reduction in fill-in by selecting either PIVOTING = 1 for Row Markowitz, or PIVOTING = 2 for column Markowitz. The Markowitz strategy will search a pre-elected number of rows or columns for pivot candidates. The default number is three, but this can be changed by using the keyword N_SEARCH_ROWS.

The keyword TOL_DROP can be used to set a tolerance which can reduce fill-in. This works by preventing any new fill element which has magnitude less than the specified drop tolerance from being added to the factorization. Since this can introduce substantial error into the factorization, it is recommended that the keyword ITER_REFINE be used to recover more accuracy in the final solution. The trade-off is between space savings from the drop tolerance and the extra time needed in repeated solve steps needed for refinement.

The IMSL_SP_LUSOL function provides the option of switching to a dense factorization method at some point during the decomposition. This option is enabled by specifying the keywords HYBRID_DENSITY and HYBRID_ORDER. HYBRID_DENSITY specifies a minimum density for the active submatrix before a format switch will occur. A density of 1.0 indicates complete fill-in. HYBRID_ORDER places an upper bound of the order of the active submatrix which will be converted to dense format. This is used to prevent a switch from occurring too early, possibly when the O(n3) nature of the dense factorization will cause performance degradation. Note that this option can significantly increase heap storage requirements.

Example


Consider the following matrix:

Let xT = (1, 2, 3, 4, 5, 6) so that Ax = (10, 7, 45, 33, –34, 31)T. The number of nonzeros in A is 15.

A = replicate(imsl_f_sp_elem, 15)
; Define the sparse matrix A using coordinate storage format.
a(*).row = [0, 1, 1, 1, 2, 3, 3, 3, 4, 4, 4, 4, 5, 5, 5]
a(*).col = [0, 1, 2, 3, 2, 0, 3, 4, 0, 3, 4, 5, 0, 1, 5]
a(*).val = [10, 10, -3, -1, 15, -2, 10, -1, -1, -5, $
1, -3, -1, -2, 6]
b = [10, 7, 45, 33, -34, 31]
; Define the right-hand side.
x = IMSL_SP_LUSOL(b, a)
; Call IMSL_SP_LUSOL, then print out result and the residual.
PM, x
  1.0000000
  2.0000000
  3.0000000
  4.0000000
  5.0000000
  6.0000000
PM, IMSL_SP_MVMUL(6, 6, a, x) - b
  0.0000000
  -8.8817842e-16
  0.0000000
  0.0000000
  0.0000000
  0.0000000

Syntax


Result = IMSL_SP_LUSOL(B, [, A] [, CONDITION=variable] [, /CSC_COL] [, /CSC_ROW] [, /CSC_VAL] [, FACTOR_COORD=variable] [, GWTH_FACTOR=value] [, GWTH_LIM=value] [, /HYBRID_DENSITY] [, /HYBRID_ORDER] [,/ITER_REFINE=value] [, MEMORY_BLOCK=value] [, N_NONZERO=variable] [, N_SEARCH_ROWS=value] [, PIVOTING=value] [, SMALLEST_PVT=variable] [, STABILITY=value] [, TOL_DROP=value] [, TRANSPOSE=value])

Return Value


Structure containing the LU factorization of A.

Arguments


B

One-dimensional matrix containing the right-hand side.

A (optional)

Sparse matrix stored as an array of structures containing the coefficient matrix A(i,j). See Direct Methods and its related sections for a description of structures used for sparse matrices.

Keywords


CONDITION (optional)

Named variable into which an estimate of the L1 condition number is stored.

CSC_COL (optional)

Accept the coefficient matrix in compressed sparse column (CSC) format. See Sparse Coordinate Storage Format for a discussion of this storage scheme. The keywords CSC_COL, CSC_ROW, and CSC_VAL must be used together.

CSC_ROW (optional)

Accept the coefficient matrix in compressed sparse column (CSC) format. See Sparse Coordinate Storage Format for a discussion of this storage scheme. The keywords CSC_COL, CSC_ROW, and CSC_VAL must be used together.

CSC_VAL (optional)

Accept the coefficient matrix in compressed sparse column (CSC) format. See Sparse Coordinate Storage Format for a discussion of this storage scheme. The keywords CSC_COL, CSC_ROW, and CSC_VAL must be used together.

FACTOR_COORD (optional)

The LU factorization of A as computed by IMSL_SP_LUFAC. If this keyword is used, then the argument A should not be used. This keyword is useful if solutions to systems involving the same coefficient matrix and multiple right-hand sides will be solved. The keywords FACTOR_COORD and CONDITION cannot be used together.

GWTH_FACTOR (optional)

Named variable into which the largest element in absolute value at any stage of the Gaussian elimination divided by the largest element in absolute value in A is stored.

GWTH_LIM (optional)

The computation stops if the growth factor exceeds GWTH_LIMIT. Default is 1.0e16.

HYBRID_DENSITY (optional)

Enable the function to switch to a dense factorization method when the density of the active submatrix reaches 0.0 ≤ HYBRID_DENSITY ≤ 1.0 and the order of the active submatrix is less than or equal to HYBRID_ORDER. The keywords HYBRID_DENSITY and HYBRID_ORDER must be used together.

HYBRID_ORDER (optional)

Enable the function to switch to a dense factorization method when the density of the active submatrix reaches 0.0 ≤ HYBRID_DENSITY ≤ 1.0 and the order of the active submatrix is less than or equal to HYBRID_ORDER. The keywords HYBRID_DENSITY and HYBRID_ORDER must be used together.

ITER_REFINE (optional)

If present and nonzero, iterative refinement will be applied.

MEMORY_BLOCK (optional)

Supply the number of non-zeros which will be added to the factor if current allocations are inadequate. Default:N_ELEMENTS(a)

N_NONZERO (optional)

Named variable into which the total number of non-zeros in the factor is stored.

N_SEARCH_ROWS (optional)

The number of rows which have the least number of non-zero elements that will be searched for a pivot element. Default: 3.

PIVOTING (optional)

Scalar value specifying the pivoting method to use. For Row Markowitz, set PIVOTING to 1; for Column Markowitz, set PIVOTING to 2; and for Symmetric Markowitz, set PIVOTING to 3. Default: 3.

SMALLEST_PVT (optional)

Named variable into which the value of the pivot element of smallest magnitude that occurred during the factorization is stored.

STABILITY (optional)

The absolute value of the pivot element must be bigger than the largest element in absolute value in its row divided by STABILITY. Default: 10.0.

TOL_DROP (optional)

Possible fill-in is checked against this tolerance. If the absolute value of the new element is less than TOL_DROP, it will be discarded. Default: 0.0.

TRANSPOSE (optional)

If present and nonzero, ATx = b is solved.

Version History


6.4

Introduced

See Also


IMSL_SP_BDPDFAC, IMSL_SP_BDPDSOL, IMSL_SP_BDSOL, IMSL_SP_CG, IMSL_SP_GMRES, IMSL_SP_LUFAC, IMSL_SP_MVMUL, IMSL_SP_PDFAC, IMSL_SP_PDSOL



© 2019 Harris Geospatial Solutions, Inc. |  Legal
My Account    |    Store    |    Contact Us