~ubuntu-branches/ubuntu/lucid/pdl/lucid

« back to all changes in this revision

Viewing changes to Lib/Slatec/slatec/dpofa.f

  • Committer: Bazaar Package Importer
  • Author(s): Ben Gertzfield
  • Date: 2002-04-08 18:47:16 UTC
  • Revision ID: james.westby@ubuntu.com-20020408184716-0hf64dc96kin3htp
Tags: upstream-2.3.2
ImportĀ upstreamĀ versionĀ 2.3.2

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
*DECK DPOFA
 
2
      SUBROUTINE DPOFA (A, LDA, N, INFO)
 
3
C***BEGIN PROLOGUE  DPOFA
 
4
C***PURPOSE  Factor a real symmetric positive definite matrix.
 
5
C***LIBRARY   SLATEC (LINPACK)
 
6
C***CATEGORY  D2B1B
 
7
C***TYPE      DOUBLE PRECISION (SPOFA-S, DPOFA-D, CPOFA-C)
 
8
C***KEYWORDS  LINEAR ALGEBRA, LINPACK, MATRIX FACTORIZATION,
 
9
C             POSITIVE DEFINITE
 
10
C***AUTHOR  Moler, C. B., (U. of New Mexico)
 
11
C***DESCRIPTION
 
12
C
 
13
C     DPOFA factors a double precision symmetric positive definite
 
14
C     matrix.
 
15
C
 
16
C     DPOFA is usually called by DPOCO, but it can be called
 
17
C     directly with a saving in time if  RCOND  is not needed.
 
18
C     (time for DPOCO) = (1 + 18/N)*(time for DPOFA) .
 
19
C
 
20
C     On Entry
 
21
C
 
22
C        A       DOUBLE PRECISION(LDA, N)
 
23
C                the symmetric matrix to be factored.  Only the
 
24
C                diagonal and upper triangle are used.
 
25
C
 
26
C        LDA     INTEGER
 
27
C                the leading dimension of the array  A .
 
28
C
 
29
C        N       INTEGER
 
30
C                the order of the matrix  A .
 
31
C
 
32
C     On Return
 
33
C
 
34
C        A       an upper triangular matrix  R  so that  A = TRANS(R)*R
 
35
C                where  TRANS(R)  is the transpose.
 
36
C                The strict lower triangle is unaltered.
 
37
C                If  INFO .NE. 0 , the factorization is not complete.
 
38
C
 
39
C        INFO    INTEGER
 
40
C                = 0  for normal return.
 
41
C                = K  signals an error condition.  The leading minor
 
42
C                     of order  K  is not positive definite.
 
43
C
 
44
C***REFERENCES  J. J. Dongarra, J. R. Bunch, C. B. Moler, and G. W.
 
45
C                 Stewart, LINPACK Users' Guide, SIAM, 1979.
 
46
C***ROUTINES CALLED  DDOT
 
47
C***REVISION HISTORY  (YYMMDD)
 
48
C   780814  DATE WRITTEN
 
49
C   890531  Changed all specific intrinsics to generic.  (WRB)
 
50
C   890831  Modified array declarations.  (WRB)
 
51
C   890831  REVISION DATE from Version 3.2
 
52
C   891214  Prologue converted to Version 4.0 format.  (BAB)
 
53
C   900326  Removed duplicate information from DESCRIPTION section.
 
54
C           (WRB)
 
55
C   920501  Reformatted the REFERENCES section.  (WRB)
 
56
C***END PROLOGUE  DPOFA
 
57
      INTEGER LDA,N,INFO
 
58
      DOUBLE PRECISION A(LDA,*)
 
59
C
 
60
      DOUBLE PRECISION DDOT,T
 
61
      DOUBLE PRECISION S
 
62
      INTEGER J,JM1,K
 
63
C***FIRST EXECUTABLE STATEMENT  DPOFA
 
64
         DO 30 J = 1, N
 
65
            INFO = J
 
66
            S = 0.0D0
 
67
            JM1 = J - 1
 
68
            IF (JM1 .LT. 1) GO TO 20
 
69
            DO 10 K = 1, JM1
 
70
               T = A(K,J) - DDOT(K-1,A(1,K),1,A(1,J),1)
 
71
               T = T/A(K,K)
 
72
               A(K,J) = T
 
73
               S = S + T*T
 
74
   10       CONTINUE
 
75
   20       CONTINUE
 
76
            S = A(J,J) - S
 
77
            IF (S .LE. 0.0D0) GO TO 40
 
78
            A(J,J) = SQRT(S)
 
79
   30    CONTINUE
 
80
         INFO = 0
 
81
   40 CONTINUE
 
82
      RETURN
 
83
      END