NAME
chpsl - solve the linear system Ax = b for a Hermitian
matrix A in packed storage, which has been UDU-factored by
xHPCO or xHPFA, and vectors b and x.
SYNOPSIS
SUBROUTINE ZHPSL (ZA, N, IPIVOT, ZB)
SUBROUTINE CHPSL (CA, N, IPIVOT, CB)
#include <sunperf.h>
void zhpsl(doublecomplex *za, int n, int *ipivot , doub-
lecomplex *b) ;
void chpsl(complex *ca, int n, int *ipivot , complex *b) ;
ARGUMENTS
xA On entry, the UDU factorization of the matrix A,
as computed xHPCO or xHPFA.
N Order of the matrix A. N >= 0.
IPIVOT Pivot vector as computed by xHPCO or xHPFA.
xB On entry, the right-hand side vector b. On exit,
the solution vector x.
SAMPLE PROGRAM
PROGRAM TEST
IMPLICIT NONE
C
INTEGER LENGTA, N
PARAMETER (N = 3)
PARAMETER (LENGTA = (N * N + N) / 2)
C
REAL RCOND
COMPLEX A(LENGTA), B(N), WORK(N)
INTEGER IPIVOT(N)
C
EXTERNAL CHPCO, CHPSL
INTRINSIC CONJG
C
C Initialize the array A to store the matrix A shown below.
C Initialize the array B to store the vector b shown below.
C
C 1 1+2i 1+2i 95-180i
C A = 1+2i 6 -2+6i b = 545-118i
C 1+2i -2+6i 11 865+ 62i
C
C
DATA A / (1.0,0.0), (1.0,-2.0), (6.0,0.0),
$ (1.0,-2.0), (6.0,-2.0), (11.0,0.0) /
DATA B / (95.0,-180.0), (545.0,-118.0), (865.0,62.0) /
C
PRINT 1000
PRINT 1010, A(1), A(2), A(4)
PRINT 1010, CONJG(A(2)), A(3), A(5)
PRINT 1010, CONJG(A(4)), CONJG(A(5)), A(6)
PRINT 1020
PRINT 1030, B
CALL CHPCO (A, N, IPIVOT, RCOND, WORK)
PRINT 1040, RCOND
IF ((RCOND + 1.0) .EQ. 1.0) THEN
PRINT 1050
END IF
CALL CHPSL (A, N, IPIVOT, B)
PRINT 1060
PRINT 1030, B
C
1000 FORMAT (1X, 'A in full form:')
1010 FORMAT (3(3X, '(', F4.1, ',', F4.1, ')'))
1020 FORMAT (/1X, 'b:')
1030 FORMAT (3X, '(', F6.1, ',', F6.1, ')')
1040 FORMAT (/1X, 'Reciprocal condition number of A:', F6.3)
1050 FORMAT (1X, 'A may be singular to working precision.')
1060 FORMAT (/1X, 'A**(-1) * b:')
C
END
SAMPLE OUTPUT
A in full form:
( 1.0, 0.0) ( 1.0,-2.0) ( 1.0,-2.0)
( 1.0, 2.0) ( 6.0, 0.0) ( 6.0,-2.0)
( 1.0, 2.0) ( 6.0, 2.0) (11.0, 0.0)
b:
( 95.0,-180.0)
( 545.0,-118.0)
( 865.0, 62.0)
Reciprocal condition number of A: 0.001
A**(-1) * b:
( 5.0, 0.0)
( 26.0, 0.0)
( 64.0, 0.0)
|
Закладки на сайте Проследить за страницей |
Created 1996-2025 by Maxim Chirkov Добавить, Поддержать, Вебмастеру |