NAG Library Routine Document
F01BSF factorizes a real sparse matrix using the pivotal sequence previously obtained by F01BRF
when a matrix of the same sparsity pattern was factorized.
|SUBROUTINE F01BSF (
||N, NZ, A, LICN, IVECT, JVECT, ICN, IKEEP, IW, W, GROW, ETA, RPMIN, ABORT, IDISP, IFAIL)
||N, NZ, LICN, IVECT(NZ), JVECT(NZ), ICN(LICN), IKEEP(5*N), IW(5*N), IDISP(2), IFAIL
||A(LICN), W(N), ETA, RPMIN
F01BSF accepts as input a real sparse matrix of the same sparsity pattern as a matrix previously factorized by a call of F01BRF
. It first applies to the matrix the same permutations as were used by F01BRF
, both for permutation to block triangular form and for pivoting, and then performs Gaussian elimination to obtain the
factorization of the diagonal blocks.
Extensive data checks are made; duplicated nonzeros can be accumulated.
The factorization is intended to be used by F04AXF
to solve sparse systems of linear equations
F01BSF is much faster than F01BRF
and in some applications it is expected that there will be many calls of F01BSF for each call of F01BRF
The method is fully described in Duff (1977)
A more recent algorithm for the same calculation is provided by F11MEF
Duff I S (1977) MA28 – a set of Fortran subroutines for sparse unsymmetric linear equations AERE Report R8730 HMSO
- 1: N – INTEGERInput
On entry: , the order of the matrix .
- 2: NZ – INTEGERInput
On entry: the number of nonzero elements in the matrix .
- 3: A(LICN) – REAL (KIND=nag_wp) arrayInput/Output
On entry: , for , must contain the nonzero elements of the sparse matrix . They can be in any order since F01BSF will reorder them.
: the nonzero elements in the
factorization. The array must not
be changed by you between a call of F01BSF and a call of F04AXF
- 4: LICN – INTEGERInput
: the dimension of the arrays A
as declared in the (sub)program from which F01BSF is called. It should have the same value as it had for F01BRF
- 5: IVECT(NZ) – INTEGER arrayInput
- 6: JVECT(NZ) – INTEGER arrayInput
On entry: and , for , must contain the row index and the column index respectively of the nonzero element stored in .
- 7: ICN(LICN) – INTEGER arrayInput
contains, on entry, the same information as output by F01BRF
. It must not be changed by you between a call of F01BSF and a call of F04AXF
is used as internal workspace prior to being restored on exit and hence is unchanged.
- 8: IKEEP() – INTEGER arrayCommunication Array
: the same indexing information about the factorization as output in IKEEP
You must not
between a call of F01BSF and subsequent calls to F04AXF
- 9: IW() – INTEGER arrayWorkspace
- 10: W(N) – REAL (KIND=nag_wp) arrayOutput
contains an estimate (an upper bound) of the increase in size of elements encountered during the factorization (see GROW
); the rest of the array is used as workspace.
If , the array is not used.
- 11: GROW – LOGICALInput
, then on exit
contains an estimate (an upper bound) of the increase in size of elements encountered during the factorization. If the matrix is well-scaled (see Section 8
), then a high value for
indicates that the
factorization may be inaccurate and you should be wary of the results and perhaps increase the parameter PIVOT
for subsequent runs (see Section 7
- 12: ETA – REAL (KIND=nag_wp)Input
: the relative pivot threshold below which an error diagnostic is provoked and IFAIL
is set to
. If ETA
is greater than
, then no check on pivot size is made.
- 13: RPMIN – REAL (KIND=nag_wp)Output
: if ETA
is less than
, then RPMIN
gives the smallest ratio of the pivot to the largest element in the row of the corresponding upper triangular factor thus monitoring the stability of the factorization. If RPMIN
is very small it may be advisable to perform a new factorization using F01BRF
- 14: ABORT – LOGICALInput
, F01BSF exits immediately (with
) if it finds duplicate elements in the input matrix.
If , F01BSF proceeds using a value equal to the sum of the duplicate elements.
In either case details of each duplicate element are output on the current advisory message unit (see X04ABF
), unless suppressed by the value of IFAIL
- 15: IDISP() – INTEGER arrayCommunication Array
must be as output in IDISP
by the previous call of F01BRF
- 16: IFAIL – INTEGERInput/Output
For this routine, the normal use of IFAIL
is extended to control the printing of error and warning messages as well as specifying hard or soft failure (see Section 3.3
in the Essential Introduction).
must be set to a value with the decimal expansion
, where each of the decimal digits
must have a value of
||specifies hard failure, otherwise soft failure;
||suppresses error messages, otherwise error messages will be printed (see Section 6);
||suppresses warning messages, otherwise warning messages will be printed (see Section 6).
The recommended value for inexperienced users is (i.e., hard failure with all messages printed).
unless the routine detects an error or a warning has been flagged (see Section 6
6 Error Indicators and Warnings
If on entry
, explanatory error messages are output on the current error message unit (as defined by X04AAF
Errors or warnings detected by the routine:
On entry, an element of the input matrix has a row or column index (i.e., an element of IVECT
) outside the range
The input matrix is incompatible with the matrix factorized by the previous call of F01BRF
(see Section 8
The input matrix is numerically singular.
A very small pivot has been detected (see Section 5
). The factorization has been completed but is potentially unstable.
Duplicate elements have been found in the input matrix and the factorization has been abandoned ( on entry).
The factorization obtained is exact for a perturbed matrix whose th element differs from by less than where is the machine precision, is the growth value returned in if , and the number of Gaussian elimination operations applied to element .
is very large or RPMIN
is very small, then a fresh call of F01BRF
If you have a sequence of problems with the same sparsity pattern then F01BSF is recommended after F01BRF
has been called for one such problem. It is typically
times faster but is potentially unstable since the previous pivotal sequence is used. Further details on timing are given in the document for F01BRF
If growth estimation is performed (), then the time increases by between and . Pivot size monitoring () involves a similar overhead.
We normally expect this routine to be entered with a matrix having the same pattern of nonzeros as was earlier presented to F01BRF
. However there is no record of this pattern, but rather a record of the pattern including all fill-ins. Therefore we permit additional nonzeros in positions corresponding to fill-ins.
If singular matrices are being treated then it is also required that the present matrix be sufficiently like the previous one for the same permutations to be suitable for factorization with the same set of zero pivots.
This example factorizes the real sparse matrices
This example program simply prints the values of
returned by F01BSF. Normally the calls of F01BRF
and F01BSF would be followed by calls of F04AXF
9.1 Program Text
Program Text (f01bsfe.f90)
9.2 Program Data
Program Data (f01bsfe.d)
9.3 Program Results
Program Results (f01bsfe.r)