LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ zlahilb()

subroutine zlahilb ( integer  N,
integer  NRHS,
complex*16, dimension(lda,n)  A,
integer  LDA,
complex*16, dimension(ldx, nrhs)  X,
integer  LDX,
complex*16, dimension(ldb, nrhs)  B,
integer  LDB,
double precision, dimension(n)  WORK,
integer  INFO,
character*3  PATH 
)

ZLAHILB

Purpose:
 ZLAHILB generates an N by N scaled Hilbert matrix in A along with
 NRHS right-hand sides in B and solutions in X such that A*X=B.

 The Hilbert matrix is scaled by M = LCM(1, 2, ..., 2*N-1) so that all
 entries are integers.  The right-hand sides are the first NRHS
 columns of M * the identity matrix, and the solutions are the
 first NRHS columns of the inverse Hilbert matrix.

 The condition number of the Hilbert matrix grows exponentially with
 its size, roughly as O(e ** (3.5*N)).  Additionally, the inverse
 Hilbert matrices beyond a relatively small dimension cannot be
 generated exactly without extra precision.  Precision is exhausted
 when the largest entry in the inverse Hilbert matrix is greater than
 2 to the power of the number of bits in the fraction of the data type
 used plus one, which is 24 for single precision.

 In single, the generated solution is exact for N <= 6 and has
 small componentwise error for 7 <= N <= 11.
Parameters
[in]N
          N is INTEGER
          The dimension of the matrix A.
[in]NRHS
          NRHS is INTEGER
          The requested number of right-hand sides.
[out]A
          A is COMPLEX array, dimension (LDA, N)
          The generated scaled Hilbert matrix.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= N.
[out]X
          X is COMPLEX array, dimension (LDX, NRHS)
          The generated exact solutions.  Currently, the first NRHS
          columns of the inverse Hilbert matrix.
[in]LDX
          LDX is INTEGER
          The leading dimension of the array X.  LDX >= N.
[out]B
          B is REAL array, dimension (LDB, NRHS)
          The generated right-hand sides.  Currently, the first NRHS
          columns of LCM(1, 2, ..., 2*N-1) * the identity matrix.
[in]LDB
          LDB is INTEGER
          The leading dimension of the array B.  LDB >= N.
[out]WORK
          WORK is REAL array, dimension (N)
[out]INFO
          INFO is INTEGER
          = 0: successful exit
          = 1: N is too large; the data is still generated but may not
               be not exact.
          < 0: if INFO = -i, the i-th argument had an illegal value
[in]PATH
          PATH is CHARACTER*3
          The LAPACK path name.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 132 of file zlahilb.f.

134 *
135 * -- LAPACK test routine --
136 * -- LAPACK is a software package provided by Univ. of Tennessee, --
137 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
138 *
139 * .. Scalar Arguments ..
140  INTEGER N, NRHS, LDA, LDX, LDB, INFO
141 * .. Array Arguments ..
142  DOUBLE PRECISION WORK(N)
143  COMPLEX*16 A(LDA,N), X(LDX, NRHS), B(LDB, NRHS)
144  CHARACTER*3 PATH
145 * ..
146 *
147 * =====================================================================
148 * .. Local Scalars ..
149  INTEGER TM, TI, R
150  INTEGER M
151  INTEGER I, J
152  COMPLEX*16 TMP
153  CHARACTER*2 C2
154 * ..
155 * .. Parameters ..
156 * NMAX_EXACT the largest dimension where the generated data is
157 * exact.
158 * NMAX_APPROX the largest dimension where the generated data has
159 * a small componentwise relative error.
160 * ??? complex uses how many bits ???
161  INTEGER NMAX_EXACT, NMAX_APPROX, SIZE_D
162  parameter(nmax_exact = 6, nmax_approx = 11, size_d = 8)
163 *
164 * d's are generated from random permutation of those eight elements.
165  COMPLEX*16 d1(8), d2(8), invd1(8), invd2(8)
166  DATA d1 /(-1,0),(0,1),(-1,-1),(0,-1),(1,0),(-1,1),(1,1),(1,-1)/
167  DATA d2 /(-1,0),(0,-1),(-1,1),(0,1),(1,0),(-1,-1),(1,-1),(1,1)/
168 
169  DATA invd1 /(-1,0),(0,-1),(-.5,.5),(0,1),(1,0),
170  $ (-.5,-.5),(.5,-.5),(.5,.5)/
171  DATA invd2 /(-1,0),(0,1),(-.5,-.5),(0,-1),(1,0),
172  $ (-.5,.5),(.5,.5),(.5,-.5)/
173 * ..
174 * .. External Functions
175  EXTERNAL zlaset, lsamen
176  INTRINSIC dble
177  LOGICAL LSAMEN
178 * ..
179 * .. Executable Statements ..
180  c2 = path( 2: 3 )
181 *
182 * Test the input arguments
183 *
184  info = 0
185  IF (n .LT. 0 .OR. n .GT. nmax_approx) THEN
186  info = -1
187  ELSE IF (nrhs .LT. 0) THEN
188  info = -2
189  ELSE IF (lda .LT. n) THEN
190  info = -4
191  ELSE IF (ldx .LT. n) THEN
192  info = -6
193  ELSE IF (ldb .LT. n) THEN
194  info = -8
195  END IF
196  IF (info .LT. 0) THEN
197  CALL xerbla('ZLAHILB', -info)
198  RETURN
199  END IF
200  IF (n .GT. nmax_exact) THEN
201  info = 1
202  END IF
203 *
204 * Compute M = the LCM of the integers [1, 2*N-1]. The largest
205 * reasonable N is small enough that integers suffice (up to N = 11).
206  m = 1
207  DO i = 2, (2*n-1)
208  tm = m
209  ti = i
210  r = mod(tm, ti)
211  DO WHILE (r .NE. 0)
212  tm = ti
213  ti = r
214  r = mod(tm, ti)
215  END DO
216  m = (m / ti) * i
217  END DO
218 *
219 * Generate the scaled Hilbert matrix in A
220 * If we are testing SY routines,
221 * take D1_i = D2_i, else, D1_i = D2_i*
222  IF ( lsamen( 2, c2, 'SY' ) ) THEN
223  DO j = 1, n
224  DO i = 1, n
225  a(i, j) = d1(mod(j,size_d)+1) * (dble(m) / (i + j - 1))
226  $ * d1(mod(i,size_d)+1)
227  END DO
228  END DO
229  ELSE
230  DO j = 1, n
231  DO i = 1, n
232  a(i, j) = d1(mod(j,size_d)+1) * (dble(m) / (i + j - 1))
233  $ * d2(mod(i,size_d)+1)
234  END DO
235  END DO
236  END IF
237 *
238 * Generate matrix B as simply the first NRHS columns of M * the
239 * identity.
240  tmp = dble(m)
241  CALL zlaset('Full', n, nrhs, (0.0d+0,0.0d+0), tmp, b, ldb)
242 *
243 * Generate the true solutions in X. Because B = the first NRHS
244 * columns of M*I, the true solutions are just the first NRHS columns
245 * of the inverse Hilbert matrix.
246  work(1) = n
247  DO j = 2, n
248  work(j) = ( ( (work(j-1)/(j-1)) * (j-1 - n) ) /(j-1) )
249  $ * (n +j -1)
250  END DO
251 
252 * If we are testing SY routines,
253 * take D1_i = D2_i, else, D1_i = D2_i*
254  IF ( lsamen( 2, c2, 'SY' ) ) THEN
255  DO j = 1, nrhs
256  DO i = 1, n
257  x(i, j) = invd1(mod(j,size_d)+1) *
258  $ ((work(i)*work(j)) / (i + j - 1))
259  $ * invd1(mod(i,size_d)+1)
260  END DO
261  END DO
262  ELSE
263  DO j = 1, nrhs
264  DO i = 1, n
265  x(i, j) = invd2(mod(j,size_d)+1) *
266  $ ((work(i)*work(j)) / (i + j - 1))
267  $ * invd1(mod(i,size_d)+1)
268  END DO
269  END DO
270  END IF
logical function lsamen(N, CA, CB)
LSAMEN
Definition: lsamen.f:74
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine zlaset(UPLO, M, N, ALPHA, BETA, A, LDA)
ZLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition: zlaset.f:106
Here is the call graph for this function:
Here is the caller graph for this function: