← Back to branch summary

~ubuntu-branches/ubuntu/wily/julia/wily

« back to all changes in this revision

Viewing changes to deps/openlibm/slatec/sppfa.f

  • Committer: Package Import Robot
  • Author(s): Sébastien Villemot
  • Date: 2013-01-16 12:29:42 UTC
  • mfrom: (1.1.2)
  • Revision ID: package-import@ubuntu.com-20130116122942-x86e42akjq31repw
Tags: 0.0.0+20130107.gitd9656f41-1
* New upstream snashot
* No longer try to rebuild helpdb.jl.
   + debian/rules: remove helpdb.jl from build-arch rule
   + debian/control: move back python-sphinx to Build-Depends-Indep
* debian/copyright: reflect upstream changes
* Add Build-Conflicts on libatlas3-base (makes linalg tests fail)
* debian/rules: replace obsolete USE_DEBIAN makeflag by a list of
  USE_SYSTEM_* flags
* debian/rules: on non-x86 systems, use libm instead of openlibm
* dpkg-buildflags.patch: remove patch, applied upstream
* Refreshed other patches

Show diffs side-by-side

added added

removed removed

Lines of Context:
deps/openlibm/slatec/sppfa.f
 
1
*DECK SPPFA
 
2
      SUBROUTINE SPPFA (AP, N, INFO)
 
3
C***BEGIN PROLOGUE  SPPFA
 
4
C***PURPOSE  Factor a real symmetric positive definite matrix stored in
 
5
C            packed form.
 
6
C***LIBRARY   SLATEC (LINPACK)
 
7
C***CATEGORY  D2B1B
 
8
C***TYPE      SINGLE PRECISION (SPPFA-S, DPPFA-D, CPPFA-C)
 
9
C***KEYWORDS  LINEAR ALGEBRA, LINPACK, MATRIX FACTORIZATION, PACKED,
 
10
C             POSITIVE DEFINITE
 
11
C***AUTHOR  Moler, C. B., (U. of New Mexico)
 
12
C***DESCRIPTION
 
13
C
 
14
C     SPPFA factors a real symmetric positive definite matrix
 
15
C     stored in packed form.
 
16
C
 
17
C     SPPFA is usually called by SPPCO, but it can be called
 
18
C     directly with a saving in time if  RCOND  is not needed.
 
19
C     (Time for SPPCO) = (1 + 18/N)*(Time for SPPFA) .
 
20
C
 
21
C     On Entry
 
22
C
 
23
C        AP      REAL (N*(N+1)/2)
 
24
C                the packed form of a symmetric matrix  A .  The
 
25
C                columns of the upper triangle are stored sequentially
 
26
C                in a one-dimensional array of length  N*(N+1)/2 .
 
27
C                See comments below for details.
 
28
C
 
29
C        N       INTEGER
 
30
C                the order of the matrix  A .
 
31
C
 
32
C     On Return
 
33
C
 
34
C        AP      an upper triangular matrix  R , stored in packed
 
35
C                form, so that  A = TRANS(R)*R .
 
36
C
 
37
C        INFO    INTEGER
 
38
C                = 0  for normal return.
 
39
C                = K  if the leading minor of order  K  is not
 
40
C                     positive definite.
 
41
C
 
42
C
 
43
C     Packed Storage
 
44
C
 
45
C          The following program segment will pack the upper
 
46
C          triangle of a symmetric matrix.
 
47
C
 
48
C                K = 0
 
49
C                DO 20 J = 1, N
 
50
C                   DO 10 I = 1, J
 
51
C                      K = K + 1
 
52
C                      AP(K) = A(I,J)
 
53
C             10    CONTINUE
 
54
C             20 CONTINUE
 
55
C
 
56
C***REFERENCES  J. J. Dongarra, J. R. Bunch, C. B. Moler, and G. W.
 
57
C                 Stewart, LINPACK Users' Guide, SIAM, 1979.
 
58
C***ROUTINES CALLED  SDOT
 
59
C***REVISION HISTORY  (YYMMDD)
 
60
C   780814  DATE WRITTEN
 
61
C   890831  Modified array declarations.  (WRB)
 
62
C   890831  REVISION DATE from Version 3.2
 
63
C   891214  Prologue converted to Version 4.0 format.  (BAB)
 
64
C   900326  Removed duplicate information from DESCRIPTION section.
 
65
C           (WRB)
 
66
C   920501  Reformatted the REFERENCES section.  (WRB)
 
67
C***END PROLOGUE  SPPFA
 
68
      INTEGER N,INFO
 
69
      REAL AP(*)
 
70
C
 
71
      REAL SDOT,T
 
72
      REAL S
 
73
      INTEGER J,JJ,JM1,K,KJ,KK
 
74
C***FIRST EXECUTABLE STATEMENT  SPPFA
 
75
         JJ = 0
 
76
         DO 30 J = 1, N
 
77
            INFO = J
 
78
            S = 0.0E0
 
79
            JM1 = J - 1
 
80
            KJ = JJ
 
81
            KK = 0
 
82
            IF (JM1 .LT. 1) GO TO 20
 
83
            DO 10 K = 1, JM1
 
84
               KJ = KJ + 1
 
85
               T = AP(KJ) - SDOT(K-1,AP(KK+1),1,AP(JJ+1),1)
 
86
               KK = KK + K
 
87
               T = T/AP(KK)
 
88
               AP(KJ) = T
 
89
               S = S + T*T
 
90
   10       CONTINUE
 
91
   20       CONTINUE
 
92
            JJ = JJ + J
 
93
            S = AP(JJ) - S
 
94
            IF (S .LE. 0.0E0) GO TO 40
 
95
            AP(JJ) = SQRT(S)
 
96
   30    CONTINUE
 
97
         INFO = 0
 
98
   40 CONTINUE
 
99
      RETURN
 
100
      END