Millepede-II V04-17-03
minresqlpBlasModule.f90
Go to the documentation of this file.
1
3
4!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
5! File minresqlpBlasModule.f90
6!
7! This file contains the following BLAS subroutines
8! ddot, dnrm2
9! required by subroutine MINRESQLP.
10!
11! Contributors:
12! Sou-Cheng Choi <sctchoi@uchicago.edu>
13! Computation Institute (CI)
14! University of Chicago
15! Chicago, IL 60637, USA
16!
17! Michael Saunders <saunders@stanford.edu>
18! Systems Optimization Laboratory (SOL)
19! Stanford University
20! Stanford, CA 94305-4026, USA
21!
22! History:
23! 24 Sep 2007: All parameters declared with correct intent
24! to avoid compiler warnings.
25! 24 Oct 2007: Use real(8) instead of double precision or -r8.
26! 24 May 2011: Use a module to package the BLAS subroutines. Use real(dp)
27! instead of real(8), where dp is a constant defined in
28! minresqlpDataModule and used in other program units.
29! 12 Jul 2011: Created complex version zminresqlpBlasModule.f90
30! from real version minresqlpBlasModule.f90.
31! 03 Aug 2013: dp constants 0.d0 and 1.d0 defined with _dp.
32!+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
33
35
36 use minresqlpdatamodule, only : dp, ip, zero, one
37 implicit none
38
39 public :: ddot, dnrm2
40
41contains
42
43!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
44!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
45!
47!
48! Discussion:
49! This routine uses double precision real arithmetic.
50! This routine uses unrolled loops for increments equal to one.
51!
52! Modified:
53! 16 May 2005
54!
55! Author:
56! Jack Dongarra
57! Fortran90 translation by John Burkardt.
58!
59! Reference:
60! Jack Dongarra, Jim Bunch, Cleve Moler, Pete Stewart,
61! LINPACK User's Guide,
62! SIAM, 1979,
63! ISBN13: 978-0-898711-72-1,
64! LC: QA214.L56.
65!
66! Charles Lawson, Richard Hanson, David Kincaid, Fred Krogh,
67! Algorithm 539,
68! Basic Linear Algebra Subprograms for Fortran Usage,
69! ACM Transactions on Mathematical Software,
70! Volume 5, Number 3, September 1979, pages 308-323.
71!
72! Parameters:
73!
74! Input, integer N, the number of entries in the vectors.
75!
76! Input, real ( kind = dp ) DX(*), the first vector.
77!
78! Input, integer INCX, the increment between successive entries in DX.
79!
80! Input, real ( kind = dp ) DY(*), the second vector.
81!
82! Input, integer INCY, the increment between successive entries in DY.
83!
84! Output, real ( kind = dp ) DDOT, the sum of the product of the
85! corresponding entries of DX and DY.
86
87 real(dp) function ddot(n,dx,incx,dy,incy)
88
89 implicit none
90 integer(ip), intent(in) :: n,incx,incy
91 real(dp), intent(in) :: dx(*),dy(*)
92
93 real(dp) :: dtemp
94 integer(ip) :: i,ix,iy,m
95
96 ddot = zero
97 dtemp = zero
98 if ( n <= 0 ) then
99 return
100 end if
101
102! Code for unequal increments or equal increments
103! not equal to 1.
104
105 if ( incx /= 1 .or. incy /= 1 ) then
106
107 if ( 0 <= incx ) then
108 ix = 1
109 else
110 ix = ( - n + 1 ) * incx + 1
111 end if
112
113 if ( 0 <= incy ) then
114 iy = 1
115 else
116 iy = ( - n + 1 ) * incy + 1
117 end if
118
119 do i = 1, n
120 dtemp = dtemp + dx(ix) * dy(iy)
121 ix = ix + incx
122 iy = iy + incy
123 end do
124
125! Code for both increments equal to 1.
126
127 else
128
129 m = mod( n, 5 )
130
131 do i = 1, m
132 dtemp = dtemp + dx(i) * dy(i)
133 end do
134
135 do i = m+1, n, 5
136 dtemp = dtemp + dx(i)*dy(i) + dx(i+1)*dy(i+1) + dx(i+2)*dy(i+2) &
137 + dx(i+3)*dy(i+3) + dx(i+4)*dy(i+4)
138 end do
139
140 end if
141
142 ddot = dtemp
143 return
144end function ddot
145
146!*****************************************************************************
147!
149!
150! Discussion:
151! This routine uses real(dp) real arithmetic.
152! DNRM2 ( X ) = sqrt ( X' * X )
153!
154! Modified:
155! 16 May 2005
156!
157! Author:
158! Sven Hammarling
159! Fortran90 translation by John Burkardt.
160!
161! Reference:
162! Jack Dongarra, Jim Bunch, Cleve Moler, Pete Stewart,
163! LINPACK User's Guide,
164! SIAM, 1979,
165! ISBN13: 978-0-898711-72-1,
166! LC: QA214.L56.
167!
168! Charles Lawson, Richard Hanson, David Kincaid, Fred Krogh,
169! Algorithm 539,
170! Basic Linear Algebra Subprograms for Fortran Usage,
171! ACM Transactions on Mathematical Software,
172! Volume 5, Number 3, September 1979, pages 308-323.
173!
174! Parameters:
175!
176! Input, integer N, the number of entries in the vector.
177!
178! Input, real ( kind = dp ) X(*), the vector whose norm is to be computed.
179!
180! Input, integer INCX, the increment between successive entries of X.
181!
182! Output, real ( kind = dp ) DNRM2, the Euclidean norm of X.
183
184 real(dp) function dnrm2 ( n, x, incx )
185
186 implicit none
187 integer(ip), intent(in) :: n,incx
188 real(dp), intent(in) :: x(*)
189
190 integer(ip) :: ix
191 real(dp) :: ssq,absxi,norm,scale
192
193 if ( n < 1 .or. incx < 1 ) then
194 norm = zero
195 else if ( n == 1 ) then
196 norm = abs( x(1) )
197 else
198 scale = zero
199 ssq = one
200
201 do ix = 1, 1 + ( n - 1 )*incx, incx
202 if ( x(ix) /= zero ) then
203 absxi = abs( x(ix) )
204 if ( scale < absxi ) then
205 ssq = 1_dp + ssq * ( scale / absxi )**2
206 scale = absxi
207 else
208 ssq = ssq + ( absxi / scale )**2
209 end if
210 end if
211 end do
212 norm = scale * sqrt( ssq )
213 end if
214
215 dnrm2 = norm
216 return
217end function dnrm2
218
219end module minresqlpblasmodule
real(dp) function, public ddot(n, dx, incx, dy, incy)
DDOT forms the dot product of two vectors.
real(dp) function, public dnrm2(n, x, incx)
DNRM2 returns the euclidean norm of a vector.
Defines precision and range in real(kind=dp) and integer(kind=ip) for portability and a few constants...
real(dp), parameter, public one
real(dp), parameter, public zero
integer, parameter, public ip
integer, parameter, public dp