Digital

Name

dsskyf − Symmetric sparse matrix factorization using skyline storage scheme (Serial and Parallel Versions)

FORMAT

DSSKYF (n, au, iaudiag, nau, iparam, rparam, iwrk, rwrk, ierror)

Arguments

ninteger∗4
On entry, the order of the matrix A.
On exit, n is unchanged. 

aureal∗8
On entry, an array of length at least nau, containing the matrix A stored in the skyline storage scheme, using either the profile-in or the diagonal-out storage mode.
On exit,  au contains the transp(U)∗D∗U factorization of the matrix A.  au must remain unchanged between the call to the routine DSSKYF and any routine that uses the factors such as DSSKYS, DSSKYC, DSSKYR and DSSKYX. 

iaudiaginteger∗4
On entry, an array of length at least n for the profile-in storage mode and n+1 for the diagonal-out storage mode, containing the pointers to the locations of the diagonal elements in array AU.
On exit,  iaudiag is unchanged. 

nauinteger∗4
On entry, the number of elements in array AU. nau is also the envelope size of the symmetric part of the matrix A. For the profile-in storage mode, nau = IAUDIAG(n).  For the diagonal-out storage mode, nau = IAUDIAG(n+1) - 1. 
On exit,  nau is unchanged. 

iparaminteger∗4
An array of length at least 100, containing the integer parameters for the transp(U)∗D∗U factorization.

iparam(1): niparam
On entry, defines the length of the array IPARAM. niparam >= 100.
On exit, iparam(1) is unchanged. 

iparam(2): nrparam
On entry, defines the length of the array RPARAM. nrparam >= 100.
On exit,  iparam(2) is unchanged. 

iparam(3): niwrk
On entry, defines the size of the integer work array, IWRK.  niwrk >=2n.
On exit,  iparam(3) is unchanged. 

iparam(4): nrwrk
On entry, defines the size of the real work array, RWRK.  As the real work array is not used at present, nrwrk can be unspecified.
On exit,  iparam(4) is unchanged. 

iparam(5): iounit
On entry, defines the I/O unit number for printing error messages and information from the routine DSSKYF. The I/O unit must be opened in the calling subprogram. If iounit <= 0, no output is generated. 
On exit,  iparam(5) is unchanged. 

iparam(6): iolevel
On entry, defines the message level that determines the amount of information printed out to iounit, when iounit > 0. 

iolevel = 0 : fatal error messages only

iolevel = 1 : error messages and minimal information

iolevel = 2 : error messages and detailed information

On exit,  iparam(6) is unchanged. 

iparam(7): idefault
On entry, defines if the default values should be used in arrays IPARAM and RPARAM. If idefault = 0, then the following default values are assigned:

IPARAM(1) = niparam = 100

IPARAM(2) = nrparam = 100

IPARAM(6) = iolevel = 0

IPARAM(8) = istore = 1

IPARAM(9) = ibeg = 0

IPARAM(10) = idet = 0

IPARAM(11) = ipvt = 0

IPARAM(13) = inertia = 0

RPARAM(1) = pvt_sml = 10∗∗(-12)

If idefault = 1, then you must assign values to the above variables before the call to the DSSKYF routine. 
On exit,  iparam(7) is unchanged. 

iparam(8): istore
On entry, defines the type of storage scheme used for the skyline matrix. If istore = 1, the matrix A is stored using the profile-in storage mode; if istore = 2, the matrix A is stored using the diagonal-out storage mode.  Default: istore = 1. 
On exit,  iparam(8) is unchanged. 

iparam(9):  ibeg
On entry, defines if full or partial factorization is to be performed. If ibeg = 0, then a full factorization is performed for rows and columns 1 through n. If ibeg > 0, then a partial factorization is performed starting from rows and columns ibeg+1 through n, that is, rows and columns from 1 through ibeg have already been factorized.  Default: ibeg = 0. 
On exit,  iparam(9) is unchanged. 

iparam(10):  idet
On entry, defines if the determinant of the matrix A is to be calculated. If idet = 0, then the determinant is not calculated; if idet = 1, the determinant is calculated as det_base ∗ 10∗∗det_pwr See RPARAM(4) and RPARAM(5).  Default: idet = 0. 
On exit,  iparam(10) is unchanged. 

iparam(11): ipvt
On entry, defines if the factorization should continue when a small pivot (defined by RPARAM(1)) is encountered. If ipvt = 0 and the absolute value of the pivot element is smaller than pvt_sml = RPARAM(1), then the factorization process is stopped and control returned to the calling sub-program.  If ipvt = 1 and a pivot smaller than RPARAM(1) in absolute value is encountered in the factorization, the process continues.  If ipvt = 2 and a pivot smaller than RPARAM(1) in absolute value is encountered in the factorization, it is replaced by a predetermined value pvt_new = RPARAM(2), and the factorization is continued.  Default: ipvt = 0. 
On exit,  iparam(11) is unchanged. 

iparam(12): ipvt_loc
On entry, an unspecified variable.
On exit,  iparam(12) contains the location of the first pivot element smaller in absolute value than pvt_sml. The pivot element is returned in pvt_val = RPARAM(3). If iparam(12) = 0, then no such pivot element exists. 

iparam(13): inertia
On entry, defines if the inertia of the matrix A should be calculated. The inertia of the symmetric matrix A is the triplet of integers (ipeigen, ineigen, izeigen), consisting of the number of positive, negative and zero eigenvalues, respectively. If inertia = 0, then the inertia is not calculated; if inertia = 1, then the number of positive and negative eigenvalues are returned in ipeigen = IPARAM(14) and ineigen = IPARAM(15), respectively. An indication of the existence of zero eigenvalues is returned in izeigen = IPARAM(16).  Default: inertia  = 0. 
On exit,  iparam(13) is unchanged. 

iparam(14): ipeigen
On entry, an unspecified variable.
On exit, if inertia = 1,  iparam(14) contains the number of positive eigenvalues of the matrix A. 

iparam(15): ineigen
On entry, an unspecified variable.
On exit, if inertia = 1,  iparam(15) contains the number of negative eigenvalues of the matrix A. 

iparam(16): izeigen
On entry, an unspecified variable.
On exit, if inertia = 1,  iparam(16) indicates if the matrix A has any zero eigenvalues. If izeigen = 0, then the matrix A does not have a zero eigenvalue; if izeigen = 1, then the matrix A has at least one zero eigenvalue. 

rparamreal∗8
An array of length at least 100, containing the real parameters for the transp(U)∗D∗U factorization.

rparam(1): pvt_sml
On entry, defines the value of the pivot  element which is considered to be small. If a pivot element smaller than pvt_sml, in absolute value, is encountered in the factorization process, then, depending on the value of ipvt = IPARAM(11), the process either stops, continues or continues after the pivot is set equal to pvt_new = RPARAM(2). pvt_sml  > 0.  Recommended value: 10∗∗(-15) <= pvt_sml <= 1.  Default: pvt_sml = 10∗∗(-12). 
On exit,  rparam(1) is unchanged. 

rparam(2): pvt_new
On entry, defines the value to which the pivot element must be set if ipvt = 2 and the pivot element is less than pvt_sml in absolute value. pvt_new should be large enough to avoid overflow when calculating the reciprocal of the pivot element. 
On exit,  rparam(2) is unchanged. 

rparam(3): pvt_val
On entry, an unspecified variable.
On exit,  rparam(3) contains the value of the first pivot element smaller than pvt_sml in absolute value. This element occurs at the location returned in IPARAM(12). If no such pivot element is found, the value of pvt_val is unspecified. 

rparam(4): det_base
On entry, an unspecified variable.
On exit, defines the base for the determinant of the matrix A. If idet = 1, the determinant is calculated as det_base ∗ 10∗∗(det_pwr). 1.0 <=  det_base < 10.0. 

rparam(5): det_pwr
On entry, an unspecified variable.
On exit, defines the power for the determinant of the matrix A. If idet = 1, the determinant is calculated as det_base ∗ 10∗∗det_pwr. 

iwrkinteger∗4
On entry, an array of length at least 2n used for integer workspace.
On exit,  iwrk contains information for use by routines that use the factorization such as DSSKYS, DSSKYC, DSSKYR and DSSKYX. The first 2n elements of  iwrk should therefore be passed unchanged to these routines. 

rwrk
 real∗8
On entry, an array used for real workspace.
On exit,  rwrk is unchanged.  Presently, rwrk is not used by the routine DSSKYF.  It can be a dummy variable. 

ierrorinteger∗4
On entry, an unspecified variable.
On exit,  ierror contains the error flag.  A value of zero indicates a normal exit from the routine DSSKYF. 

Description

DSSKYF obtains the factorization of the symmetric matrix A as:

A = transp(U)∗D∗U

where D is a diagonal matrix, and U is a unit upper triangular matrix.  The matrix A is stored in a skyline form, using either the profile-in storage mode or the diagonal-out storage mode. 

The routine DSSKYF does not perform any pivoting to preserve the numerical stability of the transp(U)∗D∗U factorization. It is therefore primarily intended for the solution of symmetric positive (or negative) definite systems as they do not require pivoting for numerical stability. Caution is urged when using this routine for symmetric indefinite systems. 

If a small pivot, in absolute value, pvt_sml, is encountered in the process of factorization, you have the option of either stopping the factorization process and returning to the calling sub-program, continuing the factorization process with the small value of the pivot, or continuing after setting the pivot equal to some predetermined value, pvt_new. The location of the first occurrence of a small pivot is returned in ipvt_loc and its value in pvt_val. 

In addition to the transp(U)∗D∗U factorization, you can also obtain the determinant of A, the number of positive and negative eigenvalues of the matrix A and an indication of the existence of zero eigenvalues. A partial factorization can also be obtained by appropriately setting the value of ibeg. If ibeg > 0, then factorization begins at row and column ibeg + 1; the rows and columns from 1 to ibeg are assumed to have been already factorized. When ibeg > 0, the determinant of A, the inertia of A, and the statistics on the matrix are calculated from rows and columns ibeg + 1 through n.  If the factorization process is stopped at row i due to a small pivot, then the inertia, determinant and statistics on the matrix are evaluated for rows ibeg + 1 through (i-1). 

The data in the first 2n elements of the integer workspace, IWRK, are used in routines that use the transp(U)∗D∗U factorization, such as DSSKYS, DSSKYC and DSSKYR. This data must therefore remain unchanged between the call to DSSKYF and any one of these routines. The real workspace array, RWRK, is not used at present. 

This routine is available in both serial and parallel versions. The routine names and parameter list are identical for both versions. For information about linking to the serial or to the parallel library, refer to the DXML Reference Manual.