NAME
zspco - compute the UDU factorization and condition number
of a symmetric matrix A in packed storage. If the condition
number is not needed then xSPFA is slightly faster. It is
typical to follow a call to xSPCO with a call to xSPSL to
solve Ax = b or to xSPDI to compute the determinant,
inverse, and inertia of A.
SYNOPSIS
SUBROUTINE DSPCO (DA, N, IPIVOT, DRCOND, DWORK)
SUBROUTINE SSPCO (SA, N, IPIVOT, SRCOND, SWORK)
SUBROUTINE ZSPCO (ZA, N, IPIVOT, DRCOND, ZWORK)
SUBROUTINE CSPCO (CA, N, IPIVOT, SRCOND, CWORK)
#include <sunperf.h>
void dspco(double *dap, int n, int *kpvt, double *rcond) ;
void sspco(float *sap, int n, int *kpvt, float *rcond) ;
void zspco(doublecomplex *zap, int n, int *kpvt, double
*rcond) ;
void cspco(complex *cap, int n, int *kpvt, float *rcond) ;
ARGUMENTS
xA On entry, the upper triangle of the matrix A. On
exit, a UDU factorization of the matrix A.
N Order of the matrix A. N >= 0.
IPIVOT On exit, a vector of pivot indices.
xRCOND On exit, an estimate of the reciprocal condition
number of A. 0.0 <= RCOND <= 1.0. As the value
of RCOND gets smaller, operations with A such as
solving Ax = b may become less stable. If RCOND
satisfies RCOND + 1.0 = 1.0 then A may be singular
to working precision.
xWORK Scratch array with a dimension of N.
SAMPLE PROGRAM
PROGRAM TEST
IMPLICIT NONE
C
INTEGER LENGTA, N
PARAMETER (N = 3)
PARAMETER (LENGTA = (N * N + N) / 2)
C
DOUBLE PRECISION A(LENGTA), B(N), RCOND, WORK(N)
INTEGER IPIVOT(N)
C
EXTERNAL DSPCO, DSPSL
C
C Initialize the array A to store in packed symmetric format
C the matrix A shown below. Initialize the array B to store
C the vector b shown below.
C
C 1 0 4 30
C A = 0 2 0 b = 4
C 4 0 1 15
C
DATA A / 1.0D0, 0.0D0, 2.0D0, 4.0D0, 0.0D0, 1.0D0 /
DATA B / 3.0D1, 4.0D0, 1.5D1 /
C
PRINT 1000
PRINT 1010, A(1), A(2), A(4)
PRINT 1010, A(2), A(3), A(5)
PRINT 1010, A(4), A(5), A(6)
PRINT 1020
PRINT 1030, B
CALL DSPCO (A, N, IPIVOT, RCOND, WORK)
IF ((RCOND + 1.0D0) .EQ. RCOND) THEN
PRINT 1040
END IF
CALL DSPSL (A, N, IPIVOT, B)
PRINT 1050
PRINT 1030, B
C
1000 FORMAT (1X, 'A:')
1010 FORMAT (3(3X, F4.1))
1020 FORMAT (/1X, 'b:')
1030 FORMAT (3X, F4.1)
1040 FORMAT (1X, 'A may be singular to working precision.')
1050 FORMAT (/1X, 'A**(-1) * b:')
C
END
SAMPLE OUTPUT
A:
1.0 0.0 4.0
0.0 2.0 0.0
4.0 0.0 1.0
b:
30.0
4.0
15.0
A**(-1) * b:
2.0
2.0
7.0
|
Закладки на сайте Проследить за страницей |
Created 1996-2025 by Maxim Chirkov Добавить, Поддержать, Вебмастеру |