LAPACK  3.6.0
LAPACK: Linear Algebra PACKage
Collaboration diagram for real:

Functions

real function sasum (N, SX, INCX)
 SASUM More...
 
subroutine saxpy (N, SA, SX, INCX, SY, INCY)
 SAXPY More...
 
real function scabs1 (Z)
 SCABS1 More...
 
real function scasum (N, CX, INCX)
 SCASUM More...
 
real function scnrm2 (N, X, INCX)
 SCNRM2 More...
 
subroutine scopy (N, SX, INCX, SY, INCY)
 SCOPY More...
 
real function sdot (N, SX, INCX, SY, INCY)
 SDOT More...
 
real function sdsdot (N, SB, SX, INCX, SY, INCY)
 SDSDOT More...
 
real function snrm2 (N, X, INCX)
 SNRM2 More...
 
subroutine srot (N, SX, INCX, SY, INCY, C, S)
 SROT More...
 
subroutine srotg (SA, SB, C, S)
 SROTG More...
 
subroutine srotm (N, SX, INCX, SY, INCY, SPARAM)
 SROTM More...
 
subroutine srotmg (SD1, SD2, SX1, SY1, SPARAM)
 SROTMG More...
 
subroutine sscal (N, SA, SX, INCX)
 SSCAL More...
 
subroutine sswap (N, SX, INCX, SY, INCY)
 SSWAP More...
 

Detailed Description

This is the group of real LEVEL 1 BLAS routines.

Function Documentation

real function sasum ( integer  N,
real, dimension(*)  SX,
integer  INCX 
)

SASUM

Purpose:
    SASUM takes the sum of the absolute values.
    uses unrolled loops for increment equal to one.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011
Further Details:
     jack dongarra, linpack, 3/11/78.
     modified 3/93 to return if incx .le. 0.
     modified 12/3/93, array(1) declarations changed to array(*)

Definition at line 54 of file sasum.f.

54 *
55 * -- Reference BLAS level1 routine (version 3.4.0) --
56 * -- Reference BLAS is a software package provided by Univ. of Tennessee, --
57 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
58 * November 2011
59 *
60 * .. Scalar Arguments ..
61  INTEGER incx,n
62 * ..
63 * .. Array Arguments ..
64  REAL sx(*)
65 * ..
66 *
67 * =====================================================================
68 *
69 * .. Local Scalars ..
70  REAL stemp
71  INTEGER i,m,mp1,nincx
72 * ..
73 * .. Intrinsic Functions ..
74  INTRINSIC abs,mod
75 * ..
76  sasum = 0.0e0
77  stemp = 0.0e0
78  IF (n.LE.0 .OR. incx.LE.0) RETURN
79  IF (incx.EQ.1) THEN
80 * code for increment equal to 1
81 *
82 *
83 * clean-up loop
84 *
85  m = mod(n,6)
86  IF (m.NE.0) THEN
87  DO i = 1,m
88  stemp = stemp + abs(sx(i))
89  END DO
90  IF (n.LT.6) THEN
91  sasum = stemp
92  RETURN
93  END IF
94  END IF
95  mp1 = m + 1
96  DO i = mp1,n,6
97  stemp = stemp + abs(sx(i)) + abs(sx(i+1)) +
98  $ abs(sx(i+2)) + abs(sx(i+3)) +
99  $ abs(sx(i+4)) + abs(sx(i+5))
100  END DO
101  ELSE
102 *
103 * code for increment not equal to 1
104 *
105  nincx = n*incx
106  DO i = 1,nincx,incx
107  stemp = stemp + abs(sx(i))
108  END DO
109  END IF
110  sasum = stemp
111  RETURN
real function sasum(N, SX, INCX)
SASUM
Definition: sasum.f:54

Here is the caller graph for this function:

subroutine saxpy ( integer  N,
real  SA,
real, dimension(*)  SX,
integer  INCX,
real, dimension(*)  SY,
integer  INCY 
)

SAXPY

Purpose:
    SAXPY constant times a vector plus a vector.
    uses unrolled loops for increments equal to one.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011
Further Details:
     jack dongarra, linpack, 3/11/78.
     modified 12/3/93, array(1) declarations changed to array(*)

Definition at line 54 of file saxpy.f.

54 *
55 * -- Reference BLAS level1 routine (version 3.4.0) --
56 * -- Reference BLAS is a software package provided by Univ. of Tennessee, --
57 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
58 * November 2011
59 *
60 * .. Scalar Arguments ..
61  REAL sa
62  INTEGER incx,incy,n
63 * ..
64 * .. Array Arguments ..
65  REAL sx(*),sy(*)
66 * ..
67 *
68 * =====================================================================
69 *
70 * .. Local Scalars ..
71  INTEGER i,ix,iy,m,mp1
72 * ..
73 * .. Intrinsic Functions ..
74  INTRINSIC mod
75 * ..
76  IF (n.LE.0) RETURN
77  IF (sa.EQ.0.0) RETURN
78  IF (incx.EQ.1 .AND. incy.EQ.1) THEN
79 *
80 * code for both increments equal to 1
81 *
82 *
83 * clean-up loop
84 *
85  m = mod(n,4)
86  IF (m.NE.0) THEN
87  DO i = 1,m
88  sy(i) = sy(i) + sa*sx(i)
89  END DO
90  END IF
91  IF (n.LT.4) RETURN
92  mp1 = m + 1
93  DO i = mp1,n,4
94  sy(i) = sy(i) + sa*sx(i)
95  sy(i+1) = sy(i+1) + sa*sx(i+1)
96  sy(i+2) = sy(i+2) + sa*sx(i+2)
97  sy(i+3) = sy(i+3) + sa*sx(i+3)
98  END DO
99  ELSE
100 *
101 * code for unequal increments or equal increments
102 * not equal to 1
103 *
104  ix = 1
105  iy = 1
106  IF (incx.LT.0) ix = (-n+1)*incx + 1
107  IF (incy.LT.0) iy = (-n+1)*incy + 1
108  DO i = 1,n
109  sy(iy) = sy(iy) + sa*sx(ix)
110  ix = ix + incx
111  iy = iy + incy
112  END DO
113  END IF
114  RETURN
real function scabs1 ( complex  Z)

SCABS1

Purpose:
 SCABS1 computes |Re(.)| + |Im(.)| of a complex number
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2015

Definition at line 40 of file scabs1.f.

40 *
41 * -- Reference BLAS level1 routine (version 3.6.0) --
42 * -- Reference BLAS is a software package provided by Univ. of Tennessee, --
43 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
44 * November 2015
45 *
46 * .. Scalar Arguments ..
47  COMPLEX z
48 * ..
49 *
50 * =====================================================================
51 *
52 * .. Intrinsic Functions ..
53  INTRINSIC abs,aimag,real
54 * ..
55  scabs1 = abs(REAL(z)) + abs(aimag(z))
56  RETURN
real function scabs1(Z)
SCABS1
Definition: scabs1.f:40
real function scasum ( integer  N,
complex, dimension(*)  CX,
integer  INCX 
)

SCASUM

Purpose:
    SCASUM takes the sum of the (|Re(.)| + |Im(.)|)'s of a complex vector and
    returns a single precision result.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2015
Further Details:
     jack dongarra, linpack, 3/11/78.
     modified 3/93 to return if incx .le. 0.
     modified 12/3/93, array(1) declarations changed to array(*)

Definition at line 54 of file scasum.f.

54 *
55 * -- Reference BLAS level1 routine (version 3.6.0) --
56 * -- Reference BLAS is a software package provided by Univ. of Tennessee, --
57 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
58 * November 2015
59 *
60 * .. Scalar Arguments ..
61  INTEGER incx,n
62 * ..
63 * .. Array Arguments ..
64  COMPLEX cx(*)
65 * ..
66 *
67 * =====================================================================
68 *
69 * .. Local Scalars ..
70  REAL stemp
71  INTEGER i,nincx
72 * ..
73 * .. Intrinsic Functions ..
74  INTRINSIC abs,aimag,real
75 * ..
76  scasum = 0.0e0
77  stemp = 0.0e0
78  IF (n.LE.0 .OR. incx.LE.0) RETURN
79  IF (incx.EQ.1) THEN
80 *
81 * code for increment equal to 1
82 *
83  DO i = 1,n
84  stemp = stemp + abs(REAL(cx(i))) + abs(aimag(cx(i)))
85  END DO
86  ELSE
87 *
88 * code for increment not equal to 1
89 *
90  nincx = n*incx
91  DO i = 1,nincx,incx
92  stemp = stemp + abs(REAL(cx(i))) + abs(aimag(cx(i)))
93  END DO
94  END IF
95  scasum = stemp
96  RETURN
real function scasum(N, CX, INCX)
SCASUM
Definition: scasum.f:54

Here is the caller graph for this function:

real function scnrm2 ( integer  N,
complex, dimension(*)  X,
integer  INCX 
)

SCNRM2

Purpose:
 SCNRM2 returns the euclidean norm of a vector via the function
 name, so that

    SCNRM2 := sqrt( x**H*x )
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011
Further Details:
  -- This version written on 25-October-1982.
     Modified on 14-October-1993 to inline the call to CLASSQ.
     Sven Hammarling, Nag Ltd.

Definition at line 56 of file scnrm2.f.

56 *
57 * -- Reference BLAS level1 routine (version 3.4.0) --
58 * -- Reference BLAS is a software package provided by Univ. of Tennessee, --
59 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
60 * November 2011
61 *
62 * .. Scalar Arguments ..
63  INTEGER incx,n
64 * ..
65 * .. Array Arguments ..
66  COMPLEX x(*)
67 * ..
68 *
69 * =====================================================================
70 *
71 * .. Parameters ..
72  REAL one,zero
73  parameter(one=1.0e+0,zero=0.0e+0)
74 * ..
75 * .. Local Scalars ..
76  REAL norm,scale,ssq,temp
77  INTEGER ix
78 * ..
79 * .. Intrinsic Functions ..
80  INTRINSIC abs,aimag,REAL,sqrt
81 * ..
82  IF (n.LT.1 .OR. incx.LT.1) THEN
83  norm = zero
84  ELSE
85  scale = zero
86  ssq = one
87 * The following loop is equivalent to this call to the LAPACK
88 * auxiliary routine:
89 * CALL CLASSQ( N, X, INCX, SCALE, SSQ )
90 *
91  DO 10 ix = 1,1 + (n-1)*incx,incx
92  IF (REAL(x(ix)).NE.zero) then
93  temp = abs(REAL(x(ix)))
94  IF (scale.LT.temp) THEN
95  ssq = one + ssq* (scale/temp)**2
96  scale = temp
97  ELSE
98  ssq = ssq + (temp/scale)**2
99  END IF
100  END IF
101  IF (aimag(x(ix)).NE.zero) THEN
102  temp = abs(aimag(x(ix)))
103  IF (scale.LT.temp) THEN
104  ssq = one + ssq* (scale/temp)**2
105  scale = temp
106  ELSE
107  ssq = ssq + (temp/scale)**2
108  END IF
109  END IF
110  10 CONTINUE
111  norm = scale*sqrt(ssq)
112  END IF
113 *
114  scnrm2 = norm
115  RETURN
116 *
117 * End of SCNRM2.
118 *
real function scnrm2(N, X, INCX)
SCNRM2
Definition: scnrm2.f:56

Here is the caller graph for this function:

subroutine scopy ( integer  N,
real, dimension(*)  SX,
integer  INCX,
real, dimension(*)  SY,
integer  INCY 
)

SCOPY

Purpose:
    SCOPY copies a vector, x, to a vector, y.
    uses unrolled loops for increments equal to 1.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011
Further Details:
     jack dongarra, linpack, 3/11/78.
     modified 12/3/93, array(1) declarations changed to array(*)

Definition at line 53 of file scopy.f.

53 *
54 * -- Reference BLAS level1 routine (version 3.4.0) --
55 * -- Reference BLAS is a software package provided by Univ. of Tennessee, --
56 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
57 * November 2011
58 *
59 * .. Scalar Arguments ..
60  INTEGER incx,incy,n
61 * ..
62 * .. Array Arguments ..
63  REAL sx(*),sy(*)
64 * ..
65 *
66 * =====================================================================
67 *
68 * .. Local Scalars ..
69  INTEGER i,ix,iy,m,mp1
70 * ..
71 * .. Intrinsic Functions ..
72  INTRINSIC mod
73 * ..
74  IF (n.LE.0) RETURN
75  IF (incx.EQ.1 .AND. incy.EQ.1) THEN
76 *
77 * code for both increments equal to 1
78 *
79 *
80 * clean-up loop
81 *
82  m = mod(n,7)
83  IF (m.NE.0) THEN
84  DO i = 1,m
85  sy(i) = sx(i)
86  END DO
87  IF (n.LT.7) RETURN
88  END IF
89  mp1 = m + 1
90  DO i = mp1,n,7
91  sy(i) = sx(i)
92  sy(i+1) = sx(i+1)
93  sy(i+2) = sx(i+2)
94  sy(i+3) = sx(i+3)
95  sy(i+4) = sx(i+4)
96  sy(i+5) = sx(i+5)
97  sy(i+6) = sx(i+6)
98  END DO
99  ELSE
100 *
101 * code for unequal increments or equal increments
102 * not equal to 1
103 *
104  ix = 1
105  iy = 1
106  IF (incx.LT.0) ix = (-n+1)*incx + 1
107  IF (incy.LT.0) iy = (-n+1)*incy + 1
108  DO i = 1,n
109  sy(iy) = sx(ix)
110  ix = ix + incx
111  iy = iy + incy
112  END DO
113  END IF
114  RETURN
real function sdot ( integer  N,
real, dimension(*)  SX,
integer  INCX,
real, dimension(*)  SY,
integer  INCY 
)

SDOT

Purpose:
    SDOT forms the dot product of two vectors.
    uses unrolled loops for increments equal to one.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011
Further Details:
     jack dongarra, linpack, 3/11/78.
     modified 12/3/93, array(1) declarations changed to array(*)

Definition at line 53 of file sdot.f.

53 *
54 * -- Reference BLAS level1 routine (version 3.4.0) --
55 * -- Reference BLAS is a software package provided by Univ. of Tennessee, --
56 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
57 * November 2011
58 *
59 * .. Scalar Arguments ..
60  INTEGER incx,incy,n
61 * ..
62 * .. Array Arguments ..
63  REAL sx(*),sy(*)
64 * ..
65 *
66 * =====================================================================
67 *
68 * .. Local Scalars ..
69  REAL stemp
70  INTEGER i,ix,iy,m,mp1
71 * ..
72 * .. Intrinsic Functions ..
73  INTRINSIC mod
74 * ..
75  stemp = 0.0e0
76  sdot = 0.0e0
77  IF (n.LE.0) RETURN
78  IF (incx.EQ.1 .AND. incy.EQ.1) THEN
79 *
80 * code for both increments equal to 1
81 *
82 *
83 * clean-up loop
84 *
85  m = mod(n,5)
86  IF (m.NE.0) THEN
87  DO i = 1,m
88  stemp = stemp + sx(i)*sy(i)
89  END DO
90  IF (n.LT.5) THEN
91  sdot=stemp
92  RETURN
93  END IF
94  END IF
95  mp1 = m + 1
96  DO i = mp1,n,5
97  stemp = stemp + sx(i)*sy(i) + sx(i+1)*sy(i+1) +
98  $ sx(i+2)*sy(i+2) + sx(i+3)*sy(i+3) + sx(i+4)*sy(i+4)
99  END DO
100  ELSE
101 *
102 * code for unequal increments or equal increments
103 * not equal to 1
104 *
105  ix = 1
106  iy = 1
107  IF (incx.LT.0) ix = (-n+1)*incx + 1
108  IF (incy.LT.0) iy = (-n+1)*incy + 1
109  DO i = 1,n
110  stemp = stemp + sx(ix)*sy(iy)
111  ix = ix + incx
112  iy = iy + incy
113  END DO
114  END IF
115  sdot = stemp
116  RETURN
real function sdot(N, SX, INCX, SY, INCY)
SDOT
Definition: sdot.f:53

Here is the caller graph for this function:

real function sdsdot ( integer  N,
real  SB,
real, dimension(*)  SX,
integer  INCX,
real, dimension(*)  SY,
integer  INCY 
)

SDSDOT

Purpose:
 
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011

Definition at line 142 of file sdsdot.f.

142 *
143 * -- Reference BLAS level1 routine (version 3.4.0) --
144 * -- Reference BLAS is a software package provided by Univ. of Tennessee, --
145 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
146 * November 2011
147 *
148 * .. Scalar Arguments ..
149  REAL sb
150  INTEGER incx,incy,n
151 * ..
152 * .. Array Arguments ..
153  REAL sx(*),sy(*)
154 * ..
155 *
156 * PURPOSE
157 * =======
158 *
159 * Compute the inner product of two vectors with extended
160 * precision accumulation.
161 *
162 * Returns S.P. result with dot product accumulated in D.P.
163 * SDSDOT = SB + sum for I = 0 to N-1 of SX(LX+I*INCX)*SY(LY+I*INCY),
164 * where LX = 1 if INCX .GE. 0, else LX = 1+(1-N)*INCX, and LY is
165 * defined in a similar way using INCY.
166 *
167 * AUTHOR
168 * ======
169 * Lawson, C. L., (JPL), Hanson, R. J., (SNLA),
170 * Kincaid, D. R., (U. of Texas), Krogh, F. T., (JPL)
171 *
172 * ARGUMENTS
173 * =========
174 *
175 * N (input) INTEGER
176 * number of elements in input vector(s)
177 *
178 * SB (input) REAL
179 * single precision scalar to be added to inner product
180 *
181 * SX (input) REAL array, dimension (N)
182 * single precision vector with N elements
183 *
184 * INCX (input) INTEGER
185 * storage spacing between elements of SX
186 *
187 * SY (input) REAL array, dimension (N)
188 * single precision vector with N elements
189 *
190 * INCY (input) INTEGER
191 * storage spacing between elements of SY
192 *
193 * SDSDOT (output) REAL
194 * single precision dot product (SB if N .LE. 0)
195 *
196 * Further Details
197 * ===============
198 *
199 * REFERENCES
200 *
201 * C. L. Lawson, R. J. Hanson, D. R. Kincaid and F. T.
202 * Krogh, Basic linear algebra subprograms for Fortran
203 * usage, Algorithm No. 539, Transactions on Mathematical
204 * Software 5, 3 (September 1979), pp. 308-323.
205 *
206 * REVISION HISTORY (YYMMDD)
207 *
208 * 791001 DATE WRITTEN
209 * 890531 Changed all specific intrinsics to generic. (WRB)
210 * 890831 Modified array declarations. (WRB)
211 * 890831 REVISION DATE from Version 3.2
212 * 891214 Prologue converted to Version 4.0 format. (BAB)
213 * 920310 Corrected definition of LX in DESCRIPTION. (WRB)
214 * 920501 Reformatted the REFERENCES section. (WRB)
215 * 070118 Reformat to LAPACK coding style
216 *
217 * =====================================================================
218 *
219 * .. Local Scalars ..
220  DOUBLE PRECISION dsdot
221  INTEGER i,kx,ky,ns
222 * ..
223 * .. Intrinsic Functions ..
224  INTRINSIC dble
225 * ..
226  dsdot = sb
227  IF (n.LE.0) THEN
228  sdsdot = dsdot
229  RETURN
230  END IF
231  IF (incx.EQ.incy .AND. incx.GT.0) THEN
232 *
233 * Code for equal and positive increments.
234 *
235  ns = n*incx
236  DO i = 1,ns,incx
237  dsdot = dsdot + dble(sx(i))*dble(sy(i))
238  END DO
239  ELSE
240 *
241 * Code for unequal or nonpositive increments.
242 *
243  kx = 1
244  ky = 1
245  IF (incx.LT.0) kx = 1 + (1-n)*incx
246  IF (incy.LT.0) ky = 1 + (1-n)*incy
247  DO i = 1,n
248  dsdot = dsdot + dble(sx(kx))*dble(sy(ky))
249  kx = kx + incx
250  ky = ky + incy
251  END DO
252  END IF
253  sdsdot = dsdot
254  RETURN
double precision function dsdot(N, SX, INCX, SY, INCY)
DSDOT
Definition: dsdot.f:121
real function sdsdot(N, SB, SX, INCX, SY, INCY)
SDSDOT
Definition: sdsdot.f:142

Here is the caller graph for this function:

real function snrm2 ( integer  N,
real, dimension(*)  X,
integer  INCX 
)

SNRM2

Purpose:
 SNRM2 returns the euclidean norm of a vector via the function
 name, so that

    SNRM2 := sqrt( x'*x ).
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011
Further Details:
  -- This version written on 25-October-1982.
     Modified on 14-October-1993 to inline the call to SLASSQ.
     Sven Hammarling, Nag Ltd.

Definition at line 56 of file snrm2.f.

56 *
57 * -- Reference BLAS level1 routine (version 3.4.0) --
58 * -- Reference BLAS is a software package provided by Univ. of Tennessee, --
59 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
60 * November 2011
61 *
62 * .. Scalar Arguments ..
63  INTEGER incx,n
64 * ..
65 * .. Array Arguments ..
66  REAL x(*)
67 * ..
68 *
69 * =====================================================================
70 *
71 * .. Parameters ..
72  REAL one,zero
73  parameter(one=1.0e+0,zero=0.0e+0)
74 * ..
75 * .. Local Scalars ..
76  REAL absxi,norm,scale,ssq
77  INTEGER ix
78 * ..
79 * .. Intrinsic Functions ..
80  INTRINSIC abs,sqrt
81 * ..
82  IF (n.LT.1 .OR. incx.LT.1) THEN
83  norm = zero
84  ELSE IF (n.EQ.1) THEN
85  norm = abs(x(1))
86  ELSE
87  scale = zero
88  ssq = one
89 * The following loop is equivalent to this call to the LAPACK
90 * auxiliary routine:
91 * CALL SLASSQ( N, X, INCX, SCALE, SSQ )
92 *
93  DO 10 ix = 1,1 + (n-1)*incx,incx
94  IF (x(ix).NE.zero) THEN
95  absxi = abs(x(ix))
96  IF (scale.LT.absxi) THEN
97  ssq = one + ssq* (scale/absxi)**2
98  scale = absxi
99  ELSE
100  ssq = ssq + (absxi/scale)**2
101  END IF
102  END IF
103  10 CONTINUE
104  norm = scale*sqrt(ssq)
105  END IF
106 *
107  snrm2 = norm
108  RETURN
109 *
110 * End of SNRM2.
111 *
real function snrm2(N, X, INCX)
SNRM2
Definition: snrm2.f:56

Here is the caller graph for this function:

subroutine srot ( integer  N,
real, dimension(*)  SX,
integer  INCX,
real, dimension(*)  SY,
integer  INCY,
real  C,
real  S 
)

SROT

Purpose:
    applies a plane rotation.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011
Further Details:
     jack dongarra, linpack, 3/11/78.
     modified 12/3/93, array(1) declarations changed to array(*)

Definition at line 53 of file srot.f.

53 *
54 * -- Reference BLAS level1 routine (version 3.4.0) --
55 * -- Reference BLAS is a software package provided by Univ. of Tennessee, --
56 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
57 * November 2011
58 *
59 * .. Scalar Arguments ..
60  REAL c,s
61  INTEGER incx,incy,n
62 * ..
63 * .. Array Arguments ..
64  REAL sx(*),sy(*)
65 * ..
66 *
67 * =====================================================================
68 *
69 * .. Local Scalars ..
70  REAL stemp
71  INTEGER i,ix,iy
72 * ..
73  IF (n.LE.0) RETURN
74  IF (incx.EQ.1 .AND. incy.EQ.1) THEN
75 *
76 * code for both increments equal to 1
77 *
78  DO i = 1,n
79  stemp = c*sx(i) + s*sy(i)
80  sy(i) = c*sy(i) - s*sx(i)
81  sx(i) = stemp
82  END DO
83  ELSE
84 *
85 * code for unequal increments or equal increments not equal
86 * to 1
87 *
88  ix = 1
89  iy = 1
90  IF (incx.LT.0) ix = (-n+1)*incx + 1
91  IF (incy.LT.0) iy = (-n+1)*incy + 1
92  DO i = 1,n
93  stemp = c*sx(ix) + s*sy(iy)
94  sy(iy) = c*sy(iy) - s*sx(ix)
95  sx(ix) = stemp
96  ix = ix + incx
97  iy = iy + incy
98  END DO
99  END IF
100  RETURN

Here is the caller graph for this function:

subroutine srotg ( real  SA,
real  SB,
real  C,
real  S 
)

SROTG

Purpose:
    SROTG construct givens plane rotation.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011
Further Details:
     jack dongarra, linpack, 3/11/78.

Definition at line 48 of file srotg.f.

48 *
49 * -- Reference BLAS level1 routine (version 3.4.0) --
50 * -- Reference BLAS is a software package provided by Univ. of Tennessee, --
51 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
52 * November 2011
53 *
54 * .. Scalar Arguments ..
55  REAL c,s,sa,sb
56 * ..
57 *
58 * =====================================================================
59 *
60 * .. Local Scalars ..
61  REAL r,roe,scale,z
62 * ..
63 * .. Intrinsic Functions ..
64  INTRINSIC abs,sign,sqrt
65 * ..
66  roe = sb
67  IF (abs(sa).GT.abs(sb)) roe = sa
68  scale = abs(sa) + abs(sb)
69  IF (scale.EQ.0.0) THEN
70  c = 1.0
71  s = 0.0
72  r = 0.0
73  z = 0.0
74  ELSE
75  r = scale*sqrt((sa/scale)**2+ (sb/scale)**2)
76  r = sign(1.0,roe)*r
77  c = sa/r
78  s = sb/r
79  z = 1.0
80  IF (abs(sa).GT.abs(sb)) z = s
81  IF (abs(sb).GE.abs(sa) .AND. c.NE.0.0) z = 1.0/c
82  END IF
83  sa = r
84  sb = z
85  RETURN

Here is the caller graph for this function:

subroutine srotm ( integer  N,
real, dimension(*)  SX,
integer  INCX,
real, dimension(*)  SY,
integer  INCY,
real, dimension(5)  SPARAM 
)

SROTM

Purpose:
    APPLY THE MODIFIED GIVENS TRANSFORMATION, H, TO THE 2 BY N MATRIX

    (SX**T) , WHERE **T INDICATES TRANSPOSE. THE ELEMENTS OF SX ARE IN
    (SX**T)

    SX(LX+I*INCX), I = 0 TO N-1, WHERE LX = 1 IF INCX .GE. 0, ELSE
    LX = (-INCX)*N, AND SIMILARLY FOR SY USING USING LY AND INCY.
    WITH SPARAM(1)=SFLAG, H HAS ONE OF THE FOLLOWING FORMS..

    SFLAG=-1.E0     SFLAG=0.E0        SFLAG=1.E0     SFLAG=-2.E0

      (SH11  SH12)    (1.E0  SH12)    (SH11  1.E0)    (1.E0  0.E0)
    H=(          )    (          )    (          )    (          )
      (SH21  SH22),   (SH21  1.E0),   (-1.E0 SH22),   (0.E0  1.E0).
    SEE  SROTMG FOR A DESCRIPTION OF DATA STORAGE IN SPARAM.
Parameters
[in]N
          N is INTEGER
         number of elements in input vector(s)
[in,out]SX
          SX is REAL array, dimension N
         double precision vector with N elements
[in]INCX
          INCX is INTEGER
         storage spacing between elements of SX
[in,out]SY
          SY is REAL array, dimension N
         double precision vector with N elements
[in]INCY
          INCY is INTEGER
         storage spacing between elements of SY
[in,out]SPARAM
          SPARAM is REAL array, dimension 5
     SPARAM(1)=SFLAG
     SPARAM(2)=SH11
     SPARAM(3)=SH21
     SPARAM(4)=SH12
     SPARAM(5)=SH22
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011

Definition at line 101 of file srotm.f.

101 *
102 * -- Reference BLAS level1 routine (version 3.4.0) --
103 * -- Reference BLAS is a software package provided by Univ. of Tennessee, --
104 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
105 * November 2011
106 *
107 * .. Scalar Arguments ..
108  INTEGER incx,incy,n
109 * ..
110 * .. Array Arguments ..
111  REAL sparam(5),sx(*),sy(*)
112 * ..
113 *
114 * =====================================================================
115 *
116 * .. Local Scalars ..
117  REAL sflag,sh11,sh12,sh21,sh22,two,w,z,zero
118  INTEGER i,kx,ky,nsteps
119 * ..
120 * .. Data statements ..
121  DATA zero,two/0.e0,2.e0/
122 * ..
123 *
124  sflag = sparam(1)
125  IF (n.LE.0 .OR. (sflag+two.EQ.zero)) RETURN
126  IF (incx.EQ.incy.AND.incx.GT.0) THEN
127 *
128  nsteps = n*incx
129  IF (sflag.LT.zero) THEN
130  sh11 = sparam(2)
131  sh12 = sparam(4)
132  sh21 = sparam(3)
133  sh22 = sparam(5)
134  DO i = 1,nsteps,incx
135  w = sx(i)
136  z = sy(i)
137  sx(i) = w*sh11 + z*sh12
138  sy(i) = w*sh21 + z*sh22
139  END DO
140  ELSE IF (sflag.EQ.zero) THEN
141  sh12 = sparam(4)
142  sh21 = sparam(3)
143  DO i = 1,nsteps,incx
144  w = sx(i)
145  z = sy(i)
146  sx(i) = w + z*sh12
147  sy(i) = w*sh21 + z
148  END DO
149  ELSE
150  sh11 = sparam(2)
151  sh22 = sparam(5)
152  DO i = 1,nsteps,incx
153  w = sx(i)
154  z = sy(i)
155  sx(i) = w*sh11 + z
156  sy(i) = -w + sh22*z
157  END DO
158  END IF
159  ELSE
160  kx = 1
161  ky = 1
162  IF (incx.LT.0) kx = 1 + (1-n)*incx
163  IF (incy.LT.0) ky = 1 + (1-n)*incy
164 *
165  IF (sflag.LT.zero) THEN
166  sh11 = sparam(2)
167  sh12 = sparam(4)
168  sh21 = sparam(3)
169  sh22 = sparam(5)
170  DO i = 1,n
171  w = sx(kx)
172  z = sy(ky)
173  sx(kx) = w*sh11 + z*sh12
174  sy(ky) = w*sh21 + z*sh22
175  kx = kx + incx
176  ky = ky + incy
177  END DO
178  ELSE IF (sflag.EQ.zero) THEN
179  sh12 = sparam(4)
180  sh21 = sparam(3)
181  DO i = 1,n
182  w = sx(kx)
183  z = sy(ky)
184  sx(kx) = w + z*sh12
185  sy(ky) = w*sh21 + z
186  kx = kx + incx
187  ky = ky + incy
188  END DO
189  ELSE
190  sh11 = sparam(2)
191  sh22 = sparam(5)
192  DO i = 1,n
193  w = sx(kx)
194  z = sy(ky)
195  sx(kx) = w*sh11 + z
196  sy(ky) = -w + sh22*z
197  kx = kx + incx
198  ky = ky + incy
199  END DO
200  END IF
201  END IF
202  RETURN

Here is the caller graph for this function:

subroutine srotmg ( real  SD1,
real  SD2,
real  SX1,
real  SY1,
real, dimension(5)  SPARAM 
)

SROTMG

Purpose:
    CONSTRUCT THE MODIFIED GIVENS TRANSFORMATION MATRIX H WHICH ZEROS
    THE SECOND COMPONENT OF THE 2-VECTOR  (SQRT(SD1)*SX1,SQRT(SD2)*>    SY2)**T.
    WITH SPARAM(1)=SFLAG, H HAS ONE OF THE FOLLOWING FORMS..

    SFLAG=-1.E0     SFLAG=0.E0        SFLAG=1.E0     SFLAG=-2.E0

      (SH11  SH12)    (1.E0  SH12)    (SH11  1.E0)    (1.E0  0.E0)
    H=(          )    (          )    (          )    (          )
      (SH21  SH22),   (SH21  1.E0),   (-1.E0 SH22),   (0.E0  1.E0).
    LOCATIONS 2-4 OF SPARAM CONTAIN SH11,SH21,SH12, AND SH22
    RESPECTIVELY. (VALUES OF 1.E0, -1.E0, OR 0.E0 IMPLIED BY THE
    VALUE OF SPARAM(1) ARE NOT STORED IN SPARAM.)

    THE VALUES OF GAMSQ AND RGAMSQ SET IN THE DATA STATEMENT MAY BE
    INEXACT.  THIS IS OK AS THEY ARE ONLY USED FOR TESTING THE SIZE
    OF SD1 AND SD2.  ALL ACTUAL SCALING OF DATA IS DONE USING GAM.
Parameters
[in,out]SD1
          SD1 is REAL
[in,out]SD2
          SD2 is REAL
[in,out]SX1
          SX1 is REAL
[in]SY1
          SY1 is REAL
[in,out]SPARAM
          SPARAM is REAL array, dimension 5
     SPARAM(1)=SFLAG
     SPARAM(2)=SH11
     SPARAM(3)=SH21
     SPARAM(4)=SH12
     SPARAM(5)=SH22
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011

Definition at line 92 of file srotmg.f.

92 *
93 * -- Reference BLAS level1 routine (version 3.4.0) --
94 * -- Reference BLAS is a software package provided by Univ. of Tennessee, --
95 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
96 * November 2011
97 *
98 * .. Scalar Arguments ..
99  REAL sd1,sd2,sx1,sy1
100 * ..
101 * .. Array Arguments ..
102  REAL sparam(5)
103 * ..
104 *
105 * =====================================================================
106 *
107 * .. Local Scalars ..
108  REAL gam,gamsq,one,rgamsq,sflag,sh11,sh12,sh21,sh22,sp1,sp2,sq1,
109  $ sq2,stemp,su,two,zero
110 * ..
111 * .. Intrinsic Functions ..
112  INTRINSIC abs
113 * ..
114 * .. Data statements ..
115 *
116  DATA zero,one,two/0.e0,1.e0,2.e0/
117  DATA gam,gamsq,rgamsq/4096.e0,1.67772e7,5.96046e-8/
118 * ..
119 
120  IF (sd1.LT.zero) THEN
121 * GO ZERO-H-D-AND-SX1..
122  sflag = -one
123  sh11 = zero
124  sh12 = zero
125  sh21 = zero
126  sh22 = zero
127 *
128  sd1 = zero
129  sd2 = zero
130  sx1 = zero
131  ELSE
132 * CASE-SD1-NONNEGATIVE
133  sp2 = sd2*sy1
134  IF (sp2.EQ.zero) THEN
135  sflag = -two
136  sparam(1) = sflag
137  RETURN
138  END IF
139 * REGULAR-CASE..
140  sp1 = sd1*sx1
141  sq2 = sp2*sy1
142  sq1 = sp1*sx1
143 *
144  IF (abs(sq1).GT.abs(sq2)) THEN
145  sh21 = -sy1/sx1
146  sh12 = sp2/sp1
147 *
148  su = one - sh12*sh21
149 *
150  IF (su.GT.zero) THEN
151  sflag = zero
152  sd1 = sd1/su
153  sd2 = sd2/su
154  sx1 = sx1*su
155  END IF
156  ELSE
157 
158  IF (sq2.LT.zero) THEN
159 * GO ZERO-H-D-AND-SX1..
160  sflag = -one
161  sh11 = zero
162  sh12 = zero
163  sh21 = zero
164  sh22 = zero
165 *
166  sd1 = zero
167  sd2 = zero
168  sx1 = zero
169  ELSE
170  sflag = one
171  sh11 = sp1/sp2
172  sh22 = sx1/sy1
173  su = one + sh11*sh22
174  stemp = sd2/su
175  sd2 = sd1/su
176  sd1 = stemp
177  sx1 = sy1*su
178  END IF
179  END IF
180 
181 * PROCESURE..SCALE-CHECK
182  IF (sd1.NE.zero) THEN
183  DO WHILE ((sd1.LE.rgamsq) .OR. (sd1.GE.gamsq))
184  IF (sflag.EQ.zero) THEN
185  sh11 = one
186  sh22 = one
187  sflag = -one
188  ELSE
189  sh21 = -one
190  sh12 = one
191  sflag = -one
192  END IF
193  IF (sd1.LE.rgamsq) THEN
194  sd1 = sd1*gam**2
195  sx1 = sx1/gam
196  sh11 = sh11/gam
197  sh12 = sh12/gam
198  ELSE
199  sd1 = sd1/gam**2
200  sx1 = sx1*gam
201  sh11 = sh11*gam
202  sh12 = sh12*gam
203  END IF
204  ENDDO
205  END IF
206 
207  IF (sd2.NE.zero) THEN
208  DO WHILE ( (abs(sd2).LE.rgamsq) .OR. (abs(sd2).GE.gamsq) )
209  IF (sflag.EQ.zero) THEN
210  sh11 = one
211  sh22 = one
212  sflag = -one
213  ELSE
214  sh21 = -one
215  sh12 = one
216  sflag = -one
217  END IF
218  IF (abs(sd2).LE.rgamsq) THEN
219  sd2 = sd2*gam**2
220  sh21 = sh21/gam
221  sh22 = sh22/gam
222  ELSE
223  sd2 = sd2/gam**2
224  sh21 = sh21*gam
225  sh22 = sh22*gam
226  END IF
227  END DO
228  END IF
229 
230  END IF
231 
232  IF (sflag.LT.zero) THEN
233  sparam(2) = sh11
234  sparam(3) = sh21
235  sparam(4) = sh12
236  sparam(5) = sh22
237  ELSE IF (sflag.EQ.zero) THEN
238  sparam(3) = sh21
239  sparam(4) = sh12
240  ELSE
241  sparam(2) = sh11
242  sparam(5) = sh22
243  END IF
244 
245  sparam(1) = sflag
246  RETURN

Here is the caller graph for this function:

subroutine sscal ( integer  N,
real  SA,
real, dimension(*)  SX,
integer  INCX 
)

SSCAL

Purpose:
    scales a vector by a constant.
    uses unrolled loops for increment equal to 1.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011
Further Details:
     jack dongarra, linpack, 3/11/78.
     modified 3/93 to return if incx .le. 0.
     modified 12/3/93, array(1) declarations changed to array(*)

Definition at line 55 of file sscal.f.

55 *
56 * -- Reference BLAS level1 routine (version 3.4.0) --
57 * -- Reference BLAS is a software package provided by Univ. of Tennessee, --
58 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
59 * November 2011
60 *
61 * .. Scalar Arguments ..
62  REAL sa
63  INTEGER incx,n
64 * ..
65 * .. Array Arguments ..
66  REAL sx(*)
67 * ..
68 *
69 * =====================================================================
70 *
71 * .. Local Scalars ..
72  INTEGER i,m,mp1,nincx
73 * ..
74 * .. Intrinsic Functions ..
75  INTRINSIC mod
76 * ..
77  IF (n.LE.0 .OR. incx.LE.0) RETURN
78  IF (incx.EQ.1) THEN
79 *
80 * code for increment equal to 1
81 *
82 *
83 * clean-up loop
84 *
85  m = mod(n,5)
86  IF (m.NE.0) THEN
87  DO i = 1,m
88  sx(i) = sa*sx(i)
89  END DO
90  IF (n.LT.5) RETURN
91  END IF
92  mp1 = m + 1
93  DO i = mp1,n,5
94  sx(i) = sa*sx(i)
95  sx(i+1) = sa*sx(i+1)
96  sx(i+2) = sa*sx(i+2)
97  sx(i+3) = sa*sx(i+3)
98  sx(i+4) = sa*sx(i+4)
99  END DO
100  ELSE
101 *
102 * code for increment not equal to 1
103 *
104  nincx = n*incx
105  DO i = 1,nincx,incx
106  sx(i) = sa*sx(i)
107  END DO
108  END IF
109  RETURN
subroutine sswap ( integer  N,
real, dimension(*)  SX,
integer  INCX,
real, dimension(*)  SY,
integer  INCY 
)

SSWAP

Purpose:
    interchanges two vectors.
    uses unrolled loops for increments equal to 1.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011
Further Details:
     jack dongarra, linpack, 3/11/78.
     modified 12/3/93, array(1) declarations changed to array(*)

Definition at line 53 of file sswap.f.

53 *
54 * -- Reference BLAS level1 routine (version 3.4.0) --
55 * -- Reference BLAS is a software package provided by Univ. of Tennessee, --
56 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
57 * November 2011
58 *
59 * .. Scalar Arguments ..
60  INTEGER incx,incy,n
61 * ..
62 * .. Array Arguments ..
63  REAL sx(*),sy(*)
64 * ..
65 *
66 * =====================================================================
67 *
68 * .. Local Scalars ..
69  REAL stemp
70  INTEGER i,ix,iy,m,mp1
71 * ..
72 * .. Intrinsic Functions ..
73  INTRINSIC mod
74 * ..
75  IF (n.LE.0) RETURN
76  IF (incx.EQ.1 .AND. incy.EQ.1) THEN
77 *
78 * code for both increments equal to 1
79 *
80 *
81 * clean-up loop
82 *
83  m = mod(n,3)
84  IF (m.NE.0) THEN
85  DO i = 1,m
86  stemp = sx(i)
87  sx(i) = sy(i)
88  sy(i) = stemp
89  END DO
90  IF (n.LT.3) RETURN
91  END IF
92  mp1 = m + 1
93  DO i = mp1,n,3
94  stemp = sx(i)
95  sx(i) = sy(i)
96  sy(i) = stemp
97  stemp = sx(i+1)
98  sx(i+1) = sy(i+1)
99  sy(i+1) = stemp
100  stemp = sx(i+2)
101  sx(i+2) = sy(i+2)
102  sy(i+2) = stemp
103  END DO
104  ELSE
105 *
106 * code for unequal increments or equal increments not equal
107 * to 1
108 *
109  ix = 1
110  iy = 1
111  IF (incx.LT.0) ix = (-n+1)*incx + 1
112  IF (incy.LT.0) iy = (-n+1)*incy + 1
113  DO i = 1,n
114  stemp = sx(ix)
115  sx(ix) = sy(iy)
116  sy(iy) = stemp
117  ix = ix + incx
118  iy = iy + incy
119  END DO
120  END IF
121  RETURN