Performs the symbolic Cholesky factorization for a sparse symmetric matrix using a minimum degree ordering or a user-specified ordering, and set up the data structure for the numerical Cholesky factorization

Required Arguments

IROW — Vector of length NZ containing the row subscripts of the nonzeros in the lower triangular part of the matrix including the nonzeros on the diagonal.   (Input)

JCOL — Vector of length NZ containing the column subscripts of the nonzeros in the lower triangular part of the matrix including the nonzeros on the diagonal.   (Input)
(IROW (K), JCOL(K)) gives the row and column indices of the k-th nonzero element of the matrix stored in coordinate form. Note, IROW(K) ≥ JCOL(K).

NZSUB — Vector of length MAXSUB containing the row subscripts for the off-diagonal nonzeros in the Cholesky factor in compressed format.   (Output)

INZSUB — Vector of length N + 1 containing pointers for NZSUB. The row subscripts for the off-diagonal nonzeros in column J are stored in NZSUB from location INZSUB (J) to INZSUB(J + (ILNZ (J +1) −ILNZ(J) − 1).   (Output)

MAXNZ — Total number of off-diagonal nonzeros in the Cholesky factor.   (Output)

ILNZ — Vector of length N + 1 containing pointers to the Cholesky factor. The off-diagonal nonzeros in column J of the factor are stored from location ILNZ (J) to
ILNZ(J + 1) − 1.   (Output)
(ILNZ, NZSUB, INZSUB) sets up the data structure for the off-diagonal nonzeros of the Cholesky factor in column ordered form using compressed subscript format.

INVPER — Vector of length N containing the inverse permutation.   (Output)
INVPER (K) = I indicates that the original row K is the new row I.

Optional Arguments

N — Number of equations.   (Input)
Default: N = size (INVPER,1).

NZ — Total number of the nonzeros in the lower triangular part of the symmetric matrix, including the nonzeros on the diagonal.   (Input)
Default: NZ = size (IROW,1).

IJOB — Integer parameter selecting an ordering to permute the matrix symmetrically.   (Input)
IJOB = 0 selects the user ordering specified in IPER and reorders it so that the multifrontal method can be used in the numerical factorization.
IJOB = 1 selects the user ordering specified in IPER.
IJOB = 2 selects a minimum degree ordering.
IJOB = 3 selects a minimum degree ordering suitable for the multifrontal method in the numerical factorization.
Default: IJOB = 3.

ITWKSP — The total workspace needed.   (Input)
If the default is desired, set ITWKSP to zero.
Default: ITWKSP = 0.

MAXSUB — Number of subscripts contained in array NZSUB.   (Input/Output)
On input, MAXSUB gives the size of the array NZSUB.
Note that when default workspace (ITWKSP = 0) is used, set MAXSUB = 3 * NZ. Otherwise (ITWKSP > 0), set MAXSUB = (ITWKSP − 10 * N − 7) ∕ 4. On output, MAXSUB gives the number of subscripts used by the compressed subscript format.
Default: MAXSUB = 3*NZ.

IPER — Vector of length N containing the ordering specified by IJOB.   (Input/Output)
IPER (I) = K indicates that the original row K is the new row I.

ISPACE — The storage space needed for stack of frontal matrices.   (Output)

FORTRAN 90 Interface

Generic: Because the Fortran compiler cannot determine the precision desired from the required arguments, there is no generic Fortran 90 Interface for this routine. The specific Fortran 90 Interfaces are:

Single:     CALL LSCXD (IROW, JCOL, NZSUB, INZSUB, MAXNZ, ILNZ, INVPER [,…])

Or

                       CALL S_LSCXD (IROW, JCOL, NZSUB, INZSUB, MAXNZ, ILNZ, INVPER [,…])

Double:             CALL DLSCXD (IROW, JCOL, NZSUB, INZSUB, MAXNZ, ILNZ, INVPER [,…])

Or

          CALL D_LSCXD (IROW, JCOL, NZSUB, INZSUB, MAXNZ, ILNZ, INVPER [,…])

FORTRAN 77 Interface

Single:            CALL LSCXD (N, NZ, IROW, JCOL, IJOB, ITWKSP, MAXSUB, NZSUB, INZSUB, MAXNZ, ILNZ, IPER, INVPER, ISPACE)

Double:                              The double precision name is DLSCXD.

Description

Consider the linear equation

Ax = b

where A is sparse, positive definite and symmetric. The sparse coordinate format for the matrix A requires one real and two integer vectors. The real array a contains all the nonzeros in the lower triangle of A including the diagonal. Let the number of nonzeros be nz. The two integer arrays irow and jcol, each of length nz, contain the row and column indices for these entries in A. That is

^Airow(i),icol(i) ^{= a(i),           i = 1,  …, nz}

irow(i) ≥ jcol(i)          i = 1,  ^…, nz

with all other entries in the lower triangle of A zero.

The routine LSCXD computes a minimum degree ordering or uses a user-supplied ordering to set up the sparse data structure for the Cholesky factor, L. Then the routine LNFXD produces the numerical entries in L so that we have

P APT= LLT

Here, P is the permutation matrix determined by the ordering.

The numerical computations can be carried out in one of two ways. The first method performs the factorization using a multifrontal technique. This option requires more storage but in certain cases will be faster. The multifrontal method is based on the routines in Liu (1987). For detailed description of this method, see Liu (1990), also Duff and Reid (1983, 1984), Ashcraft (1987), Ashcraft et al. (1987), and Liu (1986, 1989). The second method is fully described in George and Liu (1981). This is just the standard factorization method based on the sparse compressed storage scheme.

Comments

1.         Workspace may be explicitly provided, if desired, by use of L2CXD. The reference is:

CALL L2CXD (N, NZ, IROW, JCOL, IJOB, MAXSUB, NZSUB, INZSUB, MAXNZ, ILNZ, IPER, INVPER, ISPACE, LIWK, IWK)

The additional arguments are as follows:

LIWK — The length of IWK, LIWK should be at least 10N + 12NZ + 7. Note that the argument MAXSUB should be set to (LIWK − 10N − 7)/4.

IWK — Integer work vector of length LIWK.

Note that the parameter ITWKSP is not an argument to this routine.

2.         Informational errors

Type   Code

4           1                  The matrix is structurally singular.

Example

As an example, the following matrix is symbolically factorized, and the result is printed:

The number of nonzeros in the lower triangle of A is nz= 10. The sparse coordinate form for the lower triangle of A is given by:

 

irow      1      2      3      3      4      4      5      5      5      5

jcol      1      2      1      3      3      4      1      2      4      5

or equivalently by

irow      4      5      5      5      1      2      3      3      4      5

jcol      4      1      2      4      1      2      1      3      3      5

 

      USE LSCXD_INT
      USE WRIRN_INT
      INTEGER    N, NZ

      PARAMETER  (N=5, NZ=10)

!

      INTEGER    ILNZ(N+1), INVPER(N), INZSUB(N+1), IPER(N),&

                 IROW(NZ), ISPACE, JCOL(NZ), MAXNZ, MAXSUB,&

                 NZSUB(3*NZ)

!

      DATA IROW/1, 2, 3, 3, 4, 4, 5, 5, 5, 5/

      DATA JCOL/1, 2, 1, 3, 3, 4, 1, 2, 4, 5/

      MAXSUB = 3 * NZ

      CALL LSCXD (IROW, JCOL, NZSUB, INZSUB, MAXNZ, ILNZ, INVPER,&
                 MAXSUB=MAXSUB, IPER=IPER)

!                                 Print results

      CALL WRIRN (' iper ', IPER, 1, N, 1)

      CALL WRIRN (' invper ',INVPER, 1, N, 1)

      CALL WRIRN (' nzsub ', NZSUB, 1, MAXSUB, 1)

      CALL WRIRN (' inzsub ', INZSUB, 1, N+1, 1)

      CALL WRIRN (' ilnz ', ILNZ, 1, N+1, 1)

      END

Output

 

          iper
1   2   3   4   5
2   1   5   4   3

        invper

1   2   3   4   5

2   1   5   4   3

      nzsub

1   2   3   4

3   5   4   5

          inzsub

1   2   3   4   5   6

1   1   3   4   4   4

           ilnz

1   2   3   4   5   6

1   2   4   6   7   7

Visual Numerics, Inc. Visual Numerics - Developers of IMSL and PV-WAVE http://www.vni.com/ PHONE: 713.784.3131 FAX:713.781.9260