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

« back to all changes in this revision

Viewing changes to src/blas/f77reference/zhpr2.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 ZHPR2(UPLO,N,ALPHA,X,INCX,Y,INCY,AP)
 
2
*     .. Scalar Arguments ..
 
3
      DOUBLE COMPLEX ALPHA
 
4
      INTEGER INCX,INCY,N
 
5
      CHARACTER UPLO
 
6
*     ..
 
7
*     .. Array Arguments ..
 
8
      DOUBLE COMPLEX AP(*),X(*),Y(*)
 
9
*     ..
 
10
*
 
11
*  Purpose
 
12
*  =======
 
13
*
 
14
*  ZHPR2  performs the hermitian rank 2 operation
 
15
*
 
16
*     A := alpha*x*conjg( y' ) + conjg( alpha )*y*conjg( x' ) + A,
 
17
*
 
18
*  where alpha is a scalar, x and y are n element vectors and A is an
 
19
*  n by n hermitian matrix, supplied in packed form.
 
20
*
 
21
*  Arguments
 
22
*  ==========
 
23
*
 
24
*  UPLO   - CHARACTER*1.
 
25
*           On entry, UPLO specifies whether the upper or lower
 
26
*           triangular part of the matrix A is supplied in the packed
 
27
*           array AP as follows:
 
28
*
 
29
*              UPLO = 'U' or 'u'   The upper triangular part of A is
 
30
*                                  supplied in AP.
 
31
*
 
32
*              UPLO = 'L' or 'l'   The lower triangular part of A is
 
33
*                                  supplied in AP.
 
34
*
 
35
*           Unchanged on exit.
 
36
*
 
37
*  N      - INTEGER.
 
38
*           On entry, N specifies the order of the matrix A.
 
39
*           N must be at least zero.
 
40
*           Unchanged on exit.
 
41
*
 
42
*  ALPHA  - COMPLEX*16      .
 
43
*           On entry, ALPHA specifies the scalar alpha.
 
44
*           Unchanged on exit.
 
45
*
 
46
*  X      - COMPLEX*16       array of dimension at least
 
47
*           ( 1 + ( n - 1 )*abs( INCX ) ).
 
48
*           Before entry, the incremented array X must contain the n
 
49
*           element vector x.
 
50
*           Unchanged on exit.
 
51
*
 
52
*  INCX   - INTEGER.
 
53
*           On entry, INCX specifies the increment for the elements of
 
54
*           X. INCX must not be zero.
 
55
*           Unchanged on exit.
 
56
*
 
57
*  Y      - COMPLEX*16       array of dimension at least
 
58
*           ( 1 + ( n - 1 )*abs( INCY ) ).
 
59
*           Before entry, the incremented array Y must contain the n
 
60
*           element vector y.
 
61
*           Unchanged on exit.
 
62
*
 
63
*  INCY   - INTEGER.
 
64
*           On entry, INCY specifies the increment for the elements of
 
65
*           Y. INCY must not be zero.
 
66
*           Unchanged on exit.
 
67
*
 
68
*  AP     - COMPLEX*16       array of DIMENSION at least
 
69
*           ( ( n*( n + 1 ) )/2 ).
 
70
*           Before entry with  UPLO = 'U' or 'u', the array AP must
 
71
*           contain the upper triangular part of the hermitian matrix
 
72
*           packed sequentially, column by column, so that AP( 1 )
 
73
*           contains a( 1, 1 ), AP( 2 ) and AP( 3 ) contain a( 1, 2 )
 
74
*           and a( 2, 2 ) respectively, and so on. On exit, the array
 
75
*           AP is overwritten by the upper triangular part of the
 
76
*           updated matrix.
 
77
*           Before entry with UPLO = 'L' or 'l', the array AP must
 
78
*           contain the lower triangular part of the hermitian matrix
 
79
*           packed sequentially, column by column, so that AP( 1 )
 
80
*           contains a( 1, 1 ), AP( 2 ) and AP( 3 ) contain a( 2, 1 )
 
81
*           and a( 3, 1 ) respectively, and so on. On exit, the array
 
82
*           AP is overwritten by the lower triangular part of the
 
83
*           updated matrix.
 
84
*           Note that the imaginary parts of the diagonal elements need
 
85
*           not be set, they are assumed to be zero, and on exit they
 
86
*           are set to zero.
 
87
*
 
88
*
 
89
*  Level 2 Blas routine.
 
90
*
 
91
*  -- Written on 22-October-1986.
 
92
*     Jack Dongarra, Argonne National Lab.
 
93
*     Jeremy Du Croz, Nag Central Office.
 
94
*     Sven Hammarling, Nag Central Office.
 
95
*     Richard Hanson, Sandia National Labs.
 
96
*
 
97
*
 
98
*     .. Parameters ..
 
99
      DOUBLE COMPLEX ZERO
 
100
      PARAMETER (ZERO= (0.0D+0,0.0D+0))
 
101
*     ..
 
102
*     .. Local Scalars ..
 
103
      DOUBLE COMPLEX TEMP1,TEMP2
 
104
      INTEGER I,INFO,IX,IY,J,JX,JY,K,KK,KX,KY
 
105
*     ..
 
106
*     .. External Functions ..
 
107
      LOGICAL LSAME
 
108
      EXTERNAL LSAME
 
109
*     ..
 
110
*     .. External Subroutines ..
 
111
      EXTERNAL XERBLA
 
112
*     ..
 
113
*     .. Intrinsic Functions ..
 
114
      INTRINSIC DBLE,DCONJG
 
115
*     ..
 
116
*
 
117
*     Test the input parameters.
 
118
*
 
119
      INFO = 0
 
120
      IF (.NOT.LSAME(UPLO,'U') .AND. .NOT.LSAME(UPLO,'L')) THEN
 
121
          INFO = 1
 
122
      ELSE IF (N.LT.0) THEN
 
123
          INFO = 2
 
124
      ELSE IF (INCX.EQ.0) THEN
 
125
          INFO = 5
 
126
      ELSE IF (INCY.EQ.0) THEN
 
127
          INFO = 7
 
128
      END IF
 
129
      IF (INFO.NE.0) THEN
 
130
          CALL XERBLA('ZHPR2 ',INFO)
 
131
          RETURN
 
132
      END IF
 
133
*
 
134
*     Quick return if possible.
 
135
*
 
136
      IF ((N.EQ.0) .OR. (ALPHA.EQ.ZERO)) RETURN
 
137
*
 
138
*     Set up the start points in X and Y if the increments are not both
 
139
*     unity.
 
140
*
 
141
      IF ((INCX.NE.1) .OR. (INCY.NE.1)) THEN
 
142
          IF (INCX.GT.0) THEN
 
143
              KX = 1
 
144
          ELSE
 
145
              KX = 1 - (N-1)*INCX
 
146
          END IF
 
147
          IF (INCY.GT.0) THEN
 
148
              KY = 1
 
149
          ELSE
 
150
              KY = 1 - (N-1)*INCY
 
151
          END IF
 
152
          JX = KX
 
153
          JY = KY
 
154
      END IF
 
155
*
 
156
*     Start the operations. In this version the elements of the array AP
 
157
*     are accessed sequentially with one pass through AP.
 
158
*
 
159
      KK = 1
 
160
      IF (LSAME(UPLO,'U')) THEN
 
161
*
 
162
*        Form  A  when upper triangle is stored in AP.
 
163
*
 
164
          IF ((INCX.EQ.1) .AND. (INCY.EQ.1)) THEN
 
165
              DO 20 J = 1,N
 
166
                  IF ((X(J).NE.ZERO) .OR. (Y(J).NE.ZERO)) THEN
 
167
                      TEMP1 = ALPHA*DCONJG(Y(J))
 
168
                      TEMP2 = DCONJG(ALPHA*X(J))
 
169
                      K = KK
 
170
                      DO 10 I = 1,J - 1
 
171
                          AP(K) = AP(K) + X(I)*TEMP1 + Y(I)*TEMP2
 
172
                          K = K + 1
 
173
   10                 CONTINUE
 
174
                      AP(KK+J-1) = DBLE(AP(KK+J-1)) +
 
175
     +                             DBLE(X(J)*TEMP1+Y(J)*TEMP2)
 
176
                  ELSE
 
177
                      AP(KK+J-1) = DBLE(AP(KK+J-1))
 
178
                  END IF
 
179
                  KK = KK + J
 
180
   20         CONTINUE
 
181
          ELSE
 
182
              DO 40 J = 1,N
 
183
                  IF ((X(JX).NE.ZERO) .OR. (Y(JY).NE.ZERO)) THEN
 
184
                      TEMP1 = ALPHA*DCONJG(Y(JY))
 
185
                      TEMP2 = DCONJG(ALPHA*X(JX))
 
186
                      IX = KX
 
187
                      IY = KY
 
188
                      DO 30 K = KK,KK + J - 2
 
189
                          AP(K) = AP(K) + X(IX)*TEMP1 + Y(IY)*TEMP2
 
190
                          IX = IX + INCX
 
191
                          IY = IY + INCY
 
192
   30                 CONTINUE
 
193
                      AP(KK+J-1) = DBLE(AP(KK+J-1)) +
 
194
     +                             DBLE(X(JX)*TEMP1+Y(JY)*TEMP2)
 
195
                  ELSE
 
196
                      AP(KK+J-1) = DBLE(AP(KK+J-1))
 
197
                  END IF
 
198
                  JX = JX + INCX
 
199
                  JY = JY + INCY
 
200
                  KK = KK + J
 
201
   40         CONTINUE
 
202
          END IF
 
203
      ELSE
 
204
*
 
205
*        Form  A  when lower triangle is stored in AP.
 
206
*
 
207
          IF ((INCX.EQ.1) .AND. (INCY.EQ.1)) THEN
 
208
              DO 60 J = 1,N
 
209
                  IF ((X(J).NE.ZERO) .OR. (Y(J).NE.ZERO)) THEN
 
210
                      TEMP1 = ALPHA*DCONJG(Y(J))
 
211
                      TEMP2 = DCONJG(ALPHA*X(J))
 
212
                      AP(KK) = DBLE(AP(KK)) +
 
213
     +                         DBLE(X(J)*TEMP1+Y(J)*TEMP2)
 
214
                      K = KK + 1
 
215
                      DO 50 I = J + 1,N
 
216
                          AP(K) = AP(K) + X(I)*TEMP1 + Y(I)*TEMP2
 
217
                          K = K + 1
 
218
   50                 CONTINUE
 
219
                  ELSE
 
220
                      AP(KK) = DBLE(AP(KK))
 
221
                  END IF
 
222
                  KK = KK + N - J + 1
 
223
   60         CONTINUE
 
224
          ELSE
 
225
              DO 80 J = 1,N
 
226
                  IF ((X(JX).NE.ZERO) .OR. (Y(JY).NE.ZERO)) THEN
 
227
                      TEMP1 = ALPHA*DCONJG(Y(JY))
 
228
                      TEMP2 = DCONJG(ALPHA*X(JX))
 
229
                      AP(KK) = DBLE(AP(KK)) +
 
230
     +                         DBLE(X(JX)*TEMP1+Y(JY)*TEMP2)
 
231
                      IX = JX
 
232
                      IY = JY
 
233
                      DO 70 K = KK + 1,KK + N - J
 
234
                          IX = IX + INCX
 
235
                          IY = IY + INCY
 
236
                          AP(K) = AP(K) + X(IX)*TEMP1 + Y(IY)*TEMP2
 
237
   70                 CONTINUE
 
238
                  ELSE
 
239
                      AP(KK) = DBLE(AP(KK))
 
240
                  END IF
 
241
                  JX = JX + INCX
 
242
                  JY = JY + INCY
 
243
                  KK = KK + N - J + 1
 
244
   80         CONTINUE
 
245
          END IF
 
246
      END IF
 
247
*
 
248
      RETURN
 
249
*
 
250
*     End of ZHPR2 .
 
251
*
 
252
      END