LAPACK  3.8.0
LAPACK: Linear Algebra PACKage
slahilb.f
Go to the documentation of this file.
1 *> \brief \b SLAHILB
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 * Definition:
9 * ===========
10 *
11 * SUBROUTINE SLAHILB( N, NRHS, A, LDA, X, LDX, B, LDB, WORK, INFO)
12 *
13 * .. Scalar Arguments ..
14 * INTEGER N, NRHS, LDA, LDX, LDB, INFO
15 * .. Array Arguments ..
16 * REAL A(LDA, N), X(LDX, NRHS), B(LDB, NRHS), WORK(N)
17 * ..
18 *
19 *
20 *> \par Purpose:
21 * =============
22 *>
23 *> \verbatim
24 *>
25 *> SLAHILB generates an N by N scaled Hilbert matrix in A along with
26 *> NRHS right-hand sides in B and solutions in X such that A*X=B.
27 *>
28 *> The Hilbert matrix is scaled by M = LCM(1, 2, ..., 2*N-1) so that all
29 *> entries are integers. The right-hand sides are the first NRHS
30 *> columns of M * the identity matrix, and the solutions are the
31 *> first NRHS columns of the inverse Hilbert matrix.
32 *>
33 *> The condition number of the Hilbert matrix grows exponentially with
34 *> its size, roughly as O(e ** (3.5*N)). Additionally, the inverse
35 *> Hilbert matrices beyond a relatively small dimension cannot be
36 *> generated exactly without extra precision. Precision is exhausted
37 *> when the largest entry in the inverse Hilbert matrix is greater than
38 *> 2 to the power of the number of bits in the fraction of the data type
39 *> used plus one, which is 24 for single precision.
40 *>
41 *> In single, the generated solution is exact for N <= 6 and has
42 *> small componentwise error for 7 <= N <= 11.
43 *> \endverbatim
44 *
45 * Arguments:
46 * ==========
47 *
48 *> \param[in] N
49 *> \verbatim
50 *> N is INTEGER
51 *> The dimension of the matrix A.
52 *> \endverbatim
53 *>
54 *> \param[in] NRHS
55 *> \verbatim
56 *> NRHS is INTEGER
57 *> The requested number of right-hand sides.
58 *> \endverbatim
59 *>
60 *> \param[out] A
61 *> \verbatim
62 *> A is REAL array, dimension (LDA, N)
63 *> The generated scaled Hilbert matrix.
64 *> \endverbatim
65 *>
66 *> \param[in] LDA
67 *> \verbatim
68 *> LDA is INTEGER
69 *> The leading dimension of the array A. LDA >= N.
70 *> \endverbatim
71 *>
72 *> \param[out] X
73 *> \verbatim
74 *> X is REAL array, dimension (LDX, NRHS)
75 *> The generated exact solutions. Currently, the first NRHS
76 *> columns of the inverse Hilbert matrix.
77 *> \endverbatim
78 *>
79 *> \param[in] LDX
80 *> \verbatim
81 *> LDX is INTEGER
82 *> The leading dimension of the array X. LDX >= N.
83 *> \endverbatim
84 *>
85 *> \param[out] B
86 *> \verbatim
87 *> B is REAL array, dimension (LDB, NRHS)
88 *> The generated right-hand sides. Currently, the first NRHS
89 *> columns of LCM(1, 2, ..., 2*N-1) * the identity matrix.
90 *> \endverbatim
91 *>
92 *> \param[in] LDB
93 *> \verbatim
94 *> LDB is INTEGER
95 *> The leading dimension of the array B. LDB >= N.
96 *> \endverbatim
97 *>
98 *> \param[out] WORK
99 *> \verbatim
100 *> WORK is REAL array, dimension (N)
101 *> \endverbatim
102 *>
103 *> \param[out] INFO
104 *> \verbatim
105 *> INFO is INTEGER
106 *> = 0: successful exit
107 *> = 1: N is too large; the data is still generated but may not
108 *> be not exact.
109 *> < 0: if INFO = -i, the i-th argument had an illegal value
110 *> \endverbatim
111 *
112 * Authors:
113 * ========
114 *
115 *> \author Univ. of Tennessee
116 *> \author Univ. of California Berkeley
117 *> \author Univ. of Colorado Denver
118 *> \author NAG Ltd.
119 *
120 *> \date November 2017
121 *
122 *> \ingroup real_matgen
123 *
124 * =====================================================================
125  SUBROUTINE slahilb( N, NRHS, A, LDA, X, LDX, B, LDB, WORK, INFO)
126 *
127 * -- LAPACK test routine (version 3.8.0) --
128 * -- LAPACK is a software package provided by Univ. of Tennessee, --
129 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
130 * November 2017
131 *
132 * .. Scalar Arguments ..
133  INTEGER N, NRHS, LDA, LDX, LDB, INFO
134 * .. Array Arguments ..
135  REAL A(lda, n), X(ldx, nrhs), B(ldb, nrhs), WORK(n)
136 * ..
137 *
138 * =====================================================================
139 * .. Local Scalars ..
140  INTEGER TM, TI, R
141  INTEGER M
142  INTEGER I, J
143 * ..
144 * .. Parameters ..
145 * NMAX_EXACT the largest dimension where the generated data is
146 * exact.
147 * NMAX_APPROX the largest dimension where the generated data has
148 * a small componentwise relative error.
149  INTEGER NMAX_EXACT, NMAX_APPROX
150  parameter(nmax_exact = 6, nmax_approx = 11)
151 * ..
152 * .. External Subroutines ..
153  EXTERNAL xerbla
154 * ..
155 * .. External Functions
156  EXTERNAL slaset
157  INTRINSIC real
158 * ..
159 * .. Executable Statements ..
160 *
161 * Test the input arguments
162 *
163  info = 0
164  IF (n .LT. 0 .OR. n .GT. nmax_approx) THEN
165  info = -1
166  ELSE IF (nrhs .LT. 0) THEN
167  info = -2
168  ELSE IF (lda .LT. n) THEN
169  info = -4
170  ELSE IF (ldx .LT. n) THEN
171  info = -6
172  ELSE IF (ldb .LT. n) THEN
173  info = -8
174  END IF
175  IF (info .LT. 0) THEN
176  CALL xerbla('SLAHILB', -info)
177  RETURN
178  END IF
179  IF (n .GT. nmax_exact) THEN
180  info = 1
181  END IF
182 *
183 * Compute M = the LCM of the integers [1, 2*N-1]. The largest
184 * reasonable N is small enough that integers suffice (up to N = 11).
185  m = 1
186  DO i = 2, (2*n-1)
187  tm = m
188  ti = i
189  r = mod(tm, ti)
190  DO WHILE (r .NE. 0)
191  tm = ti
192  ti = r
193  r = mod(tm, ti)
194  END DO
195  m = (m / ti) * i
196  END DO
197 *
198 * Generate the scaled Hilbert matrix in A
199  DO j = 1, n
200  DO i = 1, n
201  a(i, j) = REAL(M) / (i + j - 1)
202  END DO
203  END DO
204 *
205 * Generate matrix B as simply the first NRHS columns of M * the
206 * identity.
207  CALL slaset('Full', n, nrhs, 0.0, REAL(M), B, LDB)
208 *
209 * Generate the true solutions in X. Because B = the first NRHS
210 * columns of M*I, the true solutions are just the first NRHS columns
211 * of the inverse Hilbert matrix.
212  work(1) = n
213  DO j = 2, n
214  work(j) = ( ( (work(j-1)/(j-1)) * (j-1 - n) ) /(j-1) )
215  $ * (n +j -1)
216  END DO
217 *
218  DO j = 1, nrhs
219  DO i = 1, n
220  x(i, j) = (work(i)*work(j)) / (i + j - 1)
221  END DO
222  END DO
223 *
224  END
225 
subroutine slahilb(N, NRHS, A, LDA, X, LDX, B, LDB, WORK, INFO)
SLAHILB
Definition: slahilb.f:126
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine slaset(UPLO, M, N, ALPHA, BETA, A, LDA)
SLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values...
Definition: slaset.f:112