 LAPACK  3.10.0 LAPACK: Linear Algebra PACKage

◆ clahilb()

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

CLAHILB

Purpose:
CLAHILB 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.

Definition at line 132 of file clahilb.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  REAL WORK(N)
143  COMPLEX 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 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 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 Subroutines ..
175  EXTERNAL xerbla
176 * ..
177 * .. External Functions
178  EXTERNAL claset, lsamen
179  INTRINSIC real
180  LOGICAL LSAMEN
181 * ..
182 * .. Executable Statements ..
183  c2 = path( 2: 3 )
184 *
185 * Test the input arguments
186 *
187  info = 0
188  IF (n .LT. 0 .OR. n .GT. nmax_approx) THEN
189  info = -1
190  ELSE IF (nrhs .LT. 0) THEN
191  info = -2
192  ELSE IF (lda .LT. n) THEN
193  info = -4
194  ELSE IF (ldx .LT. n) THEN
195  info = -6
196  ELSE IF (ldb .LT. n) THEN
197  info = -8
198  END IF
199  IF (info .LT. 0) THEN
200  CALL xerbla('CLAHILB', -info)
201  RETURN
202  END IF
203  IF (n .GT. nmax_exact) THEN
204  info = 1
205  END IF
206 *
207 * Compute M = the LCM of the integers [1, 2*N-1]. The largest
208 * reasonable N is small enough that integers suffice (up to N = 11).
209  m = 1
210  DO i = 2, (2*n-1)
211  tm = m
212  ti = i
213  r = mod(tm, ti)
214  DO WHILE (r .NE. 0)
215  tm = ti
216  ti = r
217  r = mod(tm, ti)
218  END DO
219  m = (m / ti) * i
220  END DO
221 *
222 * Generate the scaled Hilbert matrix in A
223 * If we are testing SY routines, take
224 * D1_i = D2_i, else, D1_i = D2_i*
225  IF ( lsamen( 2, c2, 'SY' ) ) THEN
226  DO j = 1, n
227  DO i = 1, n
228  a(i, j) = d1(mod(j,size_d)+1) * (real(m) / (i + j - 1))
229  \$ * d1(mod(i,size_d)+1)
230  END DO
231  END DO
232  ELSE
233  DO j = 1, n
234  DO i = 1, n
235  a(i, j) = d1(mod(j,size_d)+1) * (real(m) / (i + j - 1))
236  \$ * d2(mod(i,size_d)+1)
237  END DO
238  END DO
239  END IF
240 *
241 * Generate matrix B as simply the first NRHS columns of M * the
242 * identity.
243  tmp = real(m)
244  CALL claset('Full', n, nrhs, (0.0,0.0), tmp, b, ldb)
245 *
246 * Generate the true solutions in X. Because B = the first NRHS
247 * columns of M*I, the true solutions are just the first NRHS columns
248 * of the inverse Hilbert matrix.
249  work(1) = n
250  DO j = 2, n
251  work(j) = ( ( (work(j-1)/(j-1)) * (j-1 - n) ) /(j-1) )
252  \$ * (n +j -1)
253  END DO
254
255 * If we are testing SY routines,
256 * take D1_i = D2_i, else, D1_i = D2_i*
257  IF ( lsamen( 2, c2, 'SY' ) ) THEN
258  DO j = 1, nrhs
259  DO i = 1, n
260  x(i, j) =
261  \$ invd1(mod(j,size_d)+1) *
262  \$ ((work(i)*work(j)) / (i + j - 1))
263  \$ * invd1(mod(i,size_d)+1)
264  END DO
265  END DO
266  ELSE
267  DO j = 1, nrhs
268  DO i = 1, n
269  x(i, j) =
270  \$ invd2(mod(j,size_d)+1) *
271  \$ ((work(i)*work(j)) / (i + j - 1))
272  \$ * invd1(mod(i,size_d)+1)
273  END DO
274  END DO
275  END IF
logical function lsamen(N, CA, CB)
LSAMEN
Definition: lsamen.f:74
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine claset(UPLO, M, N, ALPHA, BETA, A, LDA)
CLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition: claset.f:106
Here is the call graph for this function: