home
***
CD-ROM
|
disk
|
FTP
|
other
***
search
/
Simtel MSDOS 1992 September
/
Simtel20_Sept92.cdr
/
msdos
/
fortran
/
linpklib.arc
/
SPBFA.FOR
< prev
next >
Wrap
Text File
|
1984-01-06
|
3KB
|
95 lines
SUBROUTINE SPBFA(ABD,LDA,N,M,INFO)
INTEGER LDA,N,M,INFO
REAL ABD(LDA,1)
C
C SPBFA FACTORS A DOUBLE PRECISION SYMMETRIC POSITIVE DEFINITE
C MATRIX STORED IN BAND FORM.
C
C SPBFA IS USUALLY CALLED BY SPBCO, BUT IT CAN BE CALLED
C DIRECTLY WITH A SAVING IN TIME IF RCOND IS NOT NEEDED.
C
C ON ENTRY
C
C ABD REAL(LDA, N)
C THE MATRIX TO BE FACTORED. THE COLUMNS OF THE UPPER
C TRIANGLE ARE STORED IN THE COLUMNS OF ABD AND THE
C DIAGONALS OF THE UPPER TRIANGLE ARE STORED IN THE
C ROWS OF ABD . SEE THE COMMENTS BELOW FOR DETAILS.
C
C LDA INTEGER
C THE LEADING DIMENSION OF THE ARRAY ABD .
C LDA MUST BE .GE. M + 1 .
C
C N INTEGER
C THE ORDER OF THE MATRIX A .
C
C M INTEGER
C THE NUMBER OF DIAGONALS ABOVE THE MAIN DIAGONAL.
C 0 .LE. M .LT. N .
C
C ON RETURN
C
C ABD AN UPPER TRIANGULAR MATRIX R , STORED IN BAND
C FORM, SO THAT A = TRANS(R)*R .
C
C INFO INTEGER
C = 0 FOR NORMAL RETURN.
C = K IF THE LEADING MINOR OF ORDER K IS NOT
C POSITIVE DEFINITE.
C
C BAND STORAGE
C
C IF A IS A SYMMETRIC POSITIVE DEFINITE BAND MATRIX,
C THE FOLLOWING PROGRAM SEGMENT WILL SET UP THE INPUT.
C
C M = (BAND WIDTH ABOVE DIAGONAL)
C DO 20 J = 1, N
C I1 = MAX0(1, J-M)
C DO 10 I = I1, J
C K = I-J+M+1
C ABD(K,J) = A(I,J)
C 10 CONTINUE
C 20 CONTINUE
C
C LINPACK. THIS VERSION DATED 08/14/78 .
C CLEVE MOLER, UNIVERSITY OF NEW MEXICO, ARGONNE NATIONAL LAB.
C
C SUBROUTINES AND FUNCTIONS
C
C BLAS SDOT
C FORTRAN MAX0,SQRT
C
C INTERNAL VARIABLES
C
REAL SDOT,T
REAL S
INTEGER IK,J,JK,K,MU
C BEGIN BLOCK WITH ...EXITS TO 40
C
C
DO 30 J = 1, N
INFO = J
S = 0.0E0
IK = M + 1
JK = MAX0(J-M,1)
MU = MAX0(M+2-J,1)
IF (M .LT. MU) GO TO 20
DO 10 K = MU, M
T = ABD(K,J) - SDOT(K-MU,ABD(IK,JK),1,ABD(MU,J),1)
T = T/ABD(M+1,JK)
ABD(K,J) = T
S = S + T*T
IK = IK - 1
JK = JK + 1
10 CONTINUE
20 CONTINUE
S = ABD(M+1,J) - S
C ......EXIT
IF (S .LE. 0.0E0) GO TO 40
ABD(M+1,J) = SQRT(S)
30 CONTINUE
INFO = 0
40 CONTINUE
RETURN
END