~ubuntu-branches/ubuntu/vivid/atlas/vivid

« back to all changes in this revision

Viewing changes to src/blas/f77reference/stpsv.f

  • Committer: Package Import Robot
  • Author(s): Sébastien Villemot, Sylvestre Ledru, Sébastien Villemot
  • Date: 2013-06-11 15:58:16 UTC
  • mfrom: (1.1.4) (25 sid)
  • mto: This revision was merged to the branch mainline in revision 26.
  • Revision ID: package-import@ubuntu.com-20130611155816-8xeeiziu1iml040c
Tags: 3.10.1-1
[ Sylvestre Ledru ]
* New upstream release (Closes: #609287)

[ Sébastien Villemot ]
* Provide architectural defaults (i.e. precomputed timings) for all
  release archs (except armel and mips for the time being, due to slow
  porterboxes). This will make the package build much faster and should
  eliminate transient build failures due to excessive variance in the
  timings.
* Move symlinks for lib{cblas,f77blas,atlas,lapack_atlas} out of the
  libblas.so.3 alternative and make them always present, so that
  software relying on these libs do not break when another alternative
  is selected for BLAS
* ATLAS now has improved ARM support with native asm constructs. This required
  the following tunes:
  + armel-is-v4t.diff: new patch, prevents FTBFS on armel; otherwise,
    ATLAS uses asm constructs too recent for the platform (armel is only v4t)
  + debian/rules: on armhf, define the ATL_ARM_HARDFP flag; otherwise the asm
    constructs use the soft-float ABI for passing floating points
  + on armhf, ensure that -mfloat-abi=softfp and -mcpu=vfpv3 flags are never
    used; this is implemented via a patch (armhf.diff) and by the use of fixed
    archdefs
* The generic package is now built without multi-threading, because otherwise
  the package fails to build on some single-processor machines (this required
  the introduction of a patch: fix-non-threaded-build.diff). As a side effect,
  the build of the custom package gracefully handles non-threaded
  builds. (Closes: #602524)
* Add libblas.a as slave in the libblas.so alternative (Closes: #701921)
* Add symlinks for lib{f77blas,atlas}.a in /usr/lib (Closes: #666203)
* Modify shlibs file of libatlas3-base, such that packages using
  libblas/liblapack depend on any BLAS/LAPACK alternative, while packages
  depending on ATLAS-specific libraries (e.g. libatlas.so) depend specifically
  on libatlas3-base.
* corei1.diff: remove patch, applied upstream
* Use my @debian.org email address
* Remove obsolete DM-Upload-Allowed flag
* Switch VCS to git
* Remove Conflicts/Replaces against pre-squeeze packages
* libatlas-base-dev now provides libblas.so, as libblas-dev
* No longer use -Wa,--noexecstack in CFLAGS, it makes the package FTBFS
* Do not use POWER3 arch for powerpcspe port (Closes: #701068)
* Bump to debhelper compat level 9
* README.Debian: mention that devscripts is needed to compile the custom
  package (Closes: #697431)
* Bump Standards-Version to 3.9.4. As a consequence, add Built-Using
  fields because the package embeds stuff from liblapack-pic

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
      SUBROUTINE STPSV(UPLO,TRANS,DIAG,N,AP,X,INCX)
 
2
*     .. Scalar Arguments ..
 
3
      INTEGER INCX,N
 
4
      CHARACTER DIAG,TRANS,UPLO
 
5
*     ..
 
6
*     .. Array Arguments ..
 
7
      REAL AP(*),X(*)
 
8
*     ..
 
9
*
 
10
*  Purpose
 
11
*  =======
 
12
*
 
13
*  STPSV  solves one of the systems of equations
 
14
*
 
15
*     A*x = b,   or   A'*x = b,
 
16
*
 
17
*  where b and x are n element vectors and A is an n by n unit, or
 
18
*  non-unit, upper or lower triangular matrix, supplied in packed form.
 
19
*
 
20
*  No test for singularity or near-singularity is included in this
 
21
*  routine. Such tests must be performed before calling this routine.
 
22
*
 
23
*  Arguments
 
24
*  ==========
 
25
*
 
26
*  UPLO   - CHARACTER*1.
 
27
*           On entry, UPLO specifies whether the matrix is an upper or
 
28
*           lower triangular matrix as follows:
 
29
*
 
30
*              UPLO = 'U' or 'u'   A is an upper triangular matrix.
 
31
*
 
32
*              UPLO = 'L' or 'l'   A is a lower triangular matrix.
 
33
*
 
34
*           Unchanged on exit.
 
35
*
 
36
*  TRANS  - CHARACTER*1.
 
37
*           On entry, TRANS specifies the equations to be solved as
 
38
*           follows:
 
39
*
 
40
*              TRANS = 'N' or 'n'   A*x = b.
 
41
*
 
42
*              TRANS = 'T' or 't'   A'*x = b.
 
43
*
 
44
*              TRANS = 'C' or 'c'   A'*x = b.
 
45
*
 
46
*           Unchanged on exit.
 
47
*
 
48
*  DIAG   - CHARACTER*1.
 
49
*           On entry, DIAG specifies whether or not A is unit
 
50
*           triangular as follows:
 
51
*
 
52
*              DIAG = 'U' or 'u'   A is assumed to be unit triangular.
 
53
*
 
54
*              DIAG = 'N' or 'n'   A is not assumed to be unit
 
55
*                                  triangular.
 
56
*
 
57
*           Unchanged on exit.
 
58
*
 
59
*  N      - INTEGER.
 
60
*           On entry, N specifies the order of the matrix A.
 
61
*           N must be at least zero.
 
62
*           Unchanged on exit.
 
63
*
 
64
*  AP     - REAL             array of DIMENSION at least
 
65
*           ( ( n*( n + 1 ) )/2 ).
 
66
*           Before entry with  UPLO = 'U' or 'u', the array AP must
 
67
*           contain the upper triangular matrix packed sequentially,
 
68
*           column by column, so that AP( 1 ) contains a( 1, 1 ),
 
69
*           AP( 2 ) and AP( 3 ) contain a( 1, 2 ) and a( 2, 2 )
 
70
*           respectively, and so on.
 
71
*           Before entry with UPLO = 'L' or 'l', the array AP must
 
72
*           contain the lower triangular matrix packed sequentially,
 
73
*           column by column, so that AP( 1 ) contains a( 1, 1 ),
 
74
*           AP( 2 ) and AP( 3 ) contain a( 2, 1 ) and a( 3, 1 )
 
75
*           respectively, and so on.
 
76
*           Note that when  DIAG = 'U' or 'u', the diagonal elements of
 
77
*           A are not referenced, but are assumed to be unity.
 
78
*           Unchanged on exit.
 
79
*
 
80
*  X      - REAL             array of dimension at least
 
81
*           ( 1 + ( n - 1 )*abs( INCX ) ).
 
82
*           Before entry, the incremented array X must contain the n
 
83
*           element right-hand side vector b. On exit, X is overwritten
 
84
*           with the solution vector x.
 
85
*
 
86
*  INCX   - INTEGER.
 
87
*           On entry, INCX specifies the increment for the elements of
 
88
*           X. INCX must not be zero.
 
89
*           Unchanged on exit.
 
90
*
 
91
*
 
92
*  Level 2 Blas routine.
 
93
*
 
94
*  -- Written on 22-October-1986.
 
95
*     Jack Dongarra, Argonne National Lab.
 
96
*     Jeremy Du Croz, Nag Central Office.
 
97
*     Sven Hammarling, Nag Central Office.
 
98
*     Richard Hanson, Sandia National Labs.
 
99
*
 
100
*
 
101
*     .. Parameters ..
 
102
      REAL ZERO
 
103
      PARAMETER (ZERO=0.0E+0)
 
104
*     ..
 
105
*     .. Local Scalars ..
 
106
      REAL TEMP
 
107
      INTEGER I,INFO,IX,J,JX,K,KK,KX
 
108
      LOGICAL NOUNIT
 
109
*     ..
 
110
*     .. External Functions ..
 
111
      LOGICAL LSAME
 
112
      EXTERNAL LSAME
 
113
*     ..
 
114
*     .. External Subroutines ..
 
115
      EXTERNAL XERBLA
 
116
*     ..
 
117
*
 
118
*     Test the input parameters.
 
119
*
 
120
      INFO = 0
 
121
      IF (.NOT.LSAME(UPLO,'U') .AND. .NOT.LSAME(UPLO,'L')) THEN
 
122
          INFO = 1
 
123
      ELSE IF (.NOT.LSAME(TRANS,'N') .AND. .NOT.LSAME(TRANS,'T') .AND.
 
124
     +         .NOT.LSAME(TRANS,'C')) THEN
 
125
          INFO = 2
 
126
      ELSE IF (.NOT.LSAME(DIAG,'U') .AND. .NOT.LSAME(DIAG,'N')) THEN
 
127
          INFO = 3
 
128
      ELSE IF (N.LT.0) THEN
 
129
          INFO = 4
 
130
      ELSE IF (INCX.EQ.0) THEN
 
131
          INFO = 7
 
132
      END IF
 
133
      IF (INFO.NE.0) THEN
 
134
          CALL XERBLA('STPSV ',INFO)
 
135
          RETURN
 
136
      END IF
 
137
*
 
138
*     Quick return if possible.
 
139
*
 
140
      IF (N.EQ.0) RETURN
 
141
*
 
142
      NOUNIT = LSAME(DIAG,'N')
 
143
*
 
144
*     Set up the start point in X if the increment is not unity. This
 
145
*     will be  ( N - 1 )*INCX  too small for descending loops.
 
146
*
 
147
      IF (INCX.LE.0) THEN
 
148
          KX = 1 - (N-1)*INCX
 
149
      ELSE IF (INCX.NE.1) THEN
 
150
          KX = 1
 
151
      END IF
 
152
*
 
153
*     Start the operations. In this version the elements of AP are
 
154
*     accessed sequentially with one pass through AP.
 
155
*
 
156
      IF (LSAME(TRANS,'N')) THEN
 
157
*
 
158
*        Form  x := inv( A )*x.
 
159
*
 
160
          IF (LSAME(UPLO,'U')) THEN
 
161
              KK = (N* (N+1))/2
 
162
              IF (INCX.EQ.1) THEN
 
163
                  DO 20 J = N,1,-1
 
164
                      IF (X(J).NE.ZERO) THEN
 
165
                          IF (NOUNIT) X(J) = X(J)/AP(KK)
 
166
                          TEMP = X(J)
 
167
                          K = KK - 1
 
168
                          DO 10 I = J - 1,1,-1
 
169
                              X(I) = X(I) - TEMP*AP(K)
 
170
                              K = K - 1
 
171
   10                     CONTINUE
 
172
                      END IF
 
173
                      KK = KK - J
 
174
   20             CONTINUE
 
175
              ELSE
 
176
                  JX = KX + (N-1)*INCX
 
177
                  DO 40 J = N,1,-1
 
178
                      IF (X(JX).NE.ZERO) THEN
 
179
                          IF (NOUNIT) X(JX) = X(JX)/AP(KK)
 
180
                          TEMP = X(JX)
 
181
                          IX = JX
 
182
                          DO 30 K = KK - 1,KK - J + 1,-1
 
183
                              IX = IX - INCX
 
184
                              X(IX) = X(IX) - TEMP*AP(K)
 
185
   30                     CONTINUE
 
186
                      END IF
 
187
                      JX = JX - INCX
 
188
                      KK = KK - J
 
189
   40             CONTINUE
 
190
              END IF
 
191
          ELSE
 
192
              KK = 1
 
193
              IF (INCX.EQ.1) THEN
 
194
                  DO 60 J = 1,N
 
195
                      IF (X(J).NE.ZERO) THEN
 
196
                          IF (NOUNIT) X(J) = X(J)/AP(KK)
 
197
                          TEMP = X(J)
 
198
                          K = KK + 1
 
199
                          DO 50 I = J + 1,N
 
200
                              X(I) = X(I) - TEMP*AP(K)
 
201
                              K = K + 1
 
202
   50                     CONTINUE
 
203
                      END IF
 
204
                      KK = KK + (N-J+1)
 
205
   60             CONTINUE
 
206
              ELSE
 
207
                  JX = KX
 
208
                  DO 80 J = 1,N
 
209
                      IF (X(JX).NE.ZERO) THEN
 
210
                          IF (NOUNIT) X(JX) = X(JX)/AP(KK)
 
211
                          TEMP = X(JX)
 
212
                          IX = JX
 
213
                          DO 70 K = KK + 1,KK + N - J
 
214
                              IX = IX + INCX
 
215
                              X(IX) = X(IX) - TEMP*AP(K)
 
216
   70                     CONTINUE
 
217
                      END IF
 
218
                      JX = JX + INCX
 
219
                      KK = KK + (N-J+1)
 
220
   80             CONTINUE
 
221
              END IF
 
222
          END IF
 
223
      ELSE
 
224
*
 
225
*        Form  x := inv( A' )*x.
 
226
*
 
227
          IF (LSAME(UPLO,'U')) THEN
 
228
              KK = 1
 
229
              IF (INCX.EQ.1) THEN
 
230
                  DO 100 J = 1,N
 
231
                      TEMP = X(J)
 
232
                      K = KK
 
233
                      DO 90 I = 1,J - 1
 
234
                          TEMP = TEMP - AP(K)*X(I)
 
235
                          K = K + 1
 
236
   90                 CONTINUE
 
237
                      IF (NOUNIT) TEMP = TEMP/AP(KK+J-1)
 
238
                      X(J) = TEMP
 
239
                      KK = KK + J
 
240
  100             CONTINUE
 
241
              ELSE
 
242
                  JX = KX
 
243
                  DO 120 J = 1,N
 
244
                      TEMP = X(JX)
 
245
                      IX = KX
 
246
                      DO 110 K = KK,KK + J - 2
 
247
                          TEMP = TEMP - AP(K)*X(IX)
 
248
                          IX = IX + INCX
 
249
  110                 CONTINUE
 
250
                      IF (NOUNIT) TEMP = TEMP/AP(KK+J-1)
 
251
                      X(JX) = TEMP
 
252
                      JX = JX + INCX
 
253
                      KK = KK + J
 
254
  120             CONTINUE
 
255
              END IF
 
256
          ELSE
 
257
              KK = (N* (N+1))/2
 
258
              IF (INCX.EQ.1) THEN
 
259
                  DO 140 J = N,1,-1
 
260
                      TEMP = X(J)
 
261
                      K = KK
 
262
                      DO 130 I = N,J + 1,-1
 
263
                          TEMP = TEMP - AP(K)*X(I)
 
264
                          K = K - 1
 
265
  130                 CONTINUE
 
266
                      IF (NOUNIT) TEMP = TEMP/AP(KK-N+J)
 
267
                      X(J) = TEMP
 
268
                      KK = KK - (N-J+1)
 
269
  140             CONTINUE
 
270
              ELSE
 
271
                  KX = KX + (N-1)*INCX
 
272
                  JX = KX
 
273
                  DO 160 J = N,1,-1
 
274
                      TEMP = X(JX)
 
275
                      IX = KX
 
276
                      DO 150 K = KK,KK - (N- (J+1)),-1
 
277
                          TEMP = TEMP - AP(K)*X(IX)
 
278
                          IX = IX - INCX
 
279
  150                 CONTINUE
 
280
                      IF (NOUNIT) TEMP = TEMP/AP(KK-N+J)
 
281
                      X(JX) = TEMP
 
282
                      JX = JX - INCX
 
283
                      KK = KK - (N-J+1)
 
284
  160             CONTINUE
 
285
              END IF
 
286
          END IF
 
287
      END IF
 
288
*
 
289
      RETURN
 
290
*
 
291
*     End of STPSV .
 
292
*
 
293
      END