LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
dlarhs.f
Go to the documentation of this file.
1*> \brief \b DLARHS
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 DLARHS( PATH, XTYPE, UPLO, TRANS, M, N, KL, KU, NRHS,
12* A, LDA, X, LDX, B, LDB, ISEED, INFO )
13*
14* .. Scalar Arguments ..
15* CHARACTER TRANS, UPLO, XTYPE
16* CHARACTER*3 PATH
17* INTEGER INFO, KL, KU, LDA, LDB, LDX, M, N, NRHS
18* ..
19* .. Array Arguments ..
20* INTEGER ISEED( 4 )
21* DOUBLE PRECISION A( LDA, * ), B( LDB, * ), X( LDX, * )
22* ..
23*
24*
25*> \par Purpose:
26* =============
27*>
28*> \verbatim
29*>
30*> DLARHS chooses a set of NRHS random solution vectors and sets
31*> up the right hand sides for the linear system
32*> op(A) * X = B,
33*> where op(A) = A or A**T, depending on TRANS.
34*> \endverbatim
35*
36* Arguments:
37* ==========
38*
39*> \param[in] PATH
40*> \verbatim
41*> PATH is CHARACTER*3
42*> The type of the real matrix A. PATH may be given in any
43*> combination of upper and lower case. Valid types include
44*> xGE: General m x n matrix
45*> xGB: General banded matrix
46*> xPO: Symmetric positive definite, 2-D storage
47*> xPP: Symmetric positive definite packed
48*> xPB: Symmetric positive definite banded
49*> xSY: Symmetric indefinite, 2-D storage
50*> xSP: Symmetric indefinite packed
51*> xSB: Symmetric indefinite banded
52*> xTR: Triangular
53*> xTP: Triangular packed
54*> xTB: Triangular banded
55*> xQR: General m x n matrix
56*> xLQ: General m x n matrix
57*> xQL: General m x n matrix
58*> xRQ: General m x n matrix
59*> where the leading character indicates the precision.
60*> \endverbatim
61*>
62*> \param[in] XTYPE
63*> \verbatim
64*> XTYPE is CHARACTER*1
65*> Specifies how the exact solution X will be determined:
66*> = 'N': New solution; generate a random X.
67*> = 'C': Computed; use value of X on entry.
68*> \endverbatim
69*>
70*> \param[in] UPLO
71*> \verbatim
72*> UPLO is CHARACTER*1
73*> Specifies whether the upper or lower triangular part of the
74*> matrix A is stored, if A is symmetric.
75*> = 'U': Upper triangular
76*> = 'L': Lower triangular
77*> \endverbatim
78*>
79*> \param[in] TRANS
80*> \verbatim
81*> TRANS is CHARACTER*1
82*> Used only if A is nonsymmetric; specifies the operation
83*> applied to the matrix A.
84*> = 'N': B := A * X (No transpose)
85*> = 'T': B := A**T * X (Transpose)
86*> = 'C': B := A**H * X (Conjugate transpose = Transpose)
87*> \endverbatim
88*>
89*> \param[in] M
90*> \verbatim
91*> M is INTEGER
92*> The number or rows of the matrix A. M >= 0.
93*> \endverbatim
94*>
95*> \param[in] N
96*> \verbatim
97*> N is INTEGER
98*> The number of columns of the matrix A. N >= 0.
99*> \endverbatim
100*>
101*> \param[in] KL
102*> \verbatim
103*> KL is INTEGER
104*> Used only if A is a band matrix; specifies the number of
105*> subdiagonals of A if A is a general band matrix or if A is
106*> symmetric or triangular and UPLO = 'L'; specifies the number
107*> of superdiagonals of A if A is symmetric or triangular and
108*> UPLO = 'U'. 0 <= KL <= M-1.
109*> \endverbatim
110*>
111*> \param[in] KU
112*> \verbatim
113*> KU is INTEGER
114*> Used only if A is a general band matrix or if A is
115*> triangular.
116*>
117*> If PATH = xGB, specifies the number of superdiagonals of A,
118*> and 0 <= KU <= N-1.
119*>
120*> If PATH = xTR, xTP, or xTB, specifies whether or not the
121*> matrix has unit diagonal:
122*> = 1: matrix has non-unit diagonal (default)
123*> = 2: matrix has unit diagonal
124*> \endverbatim
125*>
126*> \param[in] NRHS
127*> \verbatim
128*> NRHS is INTEGER
129*> The number of right hand side vectors in the system A*X = B.
130*> \endverbatim
131*>
132*> \param[in] A
133*> \verbatim
134*> A is DOUBLE PRECISION array, dimension (LDA,N)
135*> The test matrix whose type is given by PATH.
136*> \endverbatim
137*>
138*> \param[in] LDA
139*> \verbatim
140*> LDA is INTEGER
141*> The leading dimension of the array A.
142*> If PATH = xGB, LDA >= KL+KU+1.
143*> If PATH = xPB, xSB, xHB, or xTB, LDA >= KL+1.
144*> Otherwise, LDA >= max(1,M).
145*> \endverbatim
146*>
147*> \param[in,out] X
148*> \verbatim
149*> X is or output) DOUBLE PRECISION array, dimension(LDX,NRHS)
150*> On entry, if XTYPE = 'C' (for 'Computed'), then X contains
151*> the exact solution to the system of linear equations.
152*> On exit, if XTYPE = 'N' (for 'New'), then X is initialized
153*> with random values.
154*> \endverbatim
155*>
156*> \param[in] LDX
157*> \verbatim
158*> LDX is INTEGER
159*> The leading dimension of the array X. If TRANS = 'N',
160*> LDX >= max(1,N); if TRANS = 'T', LDX >= max(1,M).
161*> \endverbatim
162*>
163*> \param[out] B
164*> \verbatim
165*> B is DOUBLE PRECISION array, dimension (LDB,NRHS)
166*> The right hand side vector(s) for the system of equations,
167*> computed from B = op(A) * X, where op(A) is determined by
168*> TRANS.
169*> \endverbatim
170*>
171*> \param[in] LDB
172*> \verbatim
173*> LDB is INTEGER
174*> The leading dimension of the array B. If TRANS = 'N',
175*> LDB >= max(1,M); if TRANS = 'T', LDB >= max(1,N).
176*> \endverbatim
177*>
178*> \param[in,out] ISEED
179*> \verbatim
180*> ISEED is INTEGER array, dimension (4)
181*> The seed vector for the random number generator (used in
182*> DLATMS). Modified on exit.
183*> \endverbatim
184*>
185*> \param[out] INFO
186*> \verbatim
187*> INFO is INTEGER
188*> = 0: successful exit
189*> < 0: if INFO = -i, the i-th argument had an illegal value
190*> \endverbatim
191*
192* Authors:
193* ========
194*
195*> \author Univ. of Tennessee
196*> \author Univ. of California Berkeley
197*> \author Univ. of Colorado Denver
198*> \author NAG Ltd.
199*
200*> \ingroup double_eig
201*
202* =====================================================================
203 SUBROUTINE dlarhs( PATH, XTYPE, UPLO, TRANS, M, N, KL, KU, NRHS,
204 $ A, LDA, X, LDX, B, LDB, ISEED, INFO )
205*
206* -- LAPACK test routine --
207* -- LAPACK is a software package provided by Univ. of Tennessee, --
208* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
209*
210* .. Scalar Arguments ..
211 CHARACTER TRANS, UPLO, XTYPE
212 CHARACTER*3 PATH
213 INTEGER INFO, KL, KU, LDA, LDB, LDX, M, N, NRHS
214* ..
215* .. Array Arguments ..
216 INTEGER ISEED( 4 )
217 DOUBLE PRECISION A( LDA, * ), B( LDB, * ), X( LDX, * )
218* ..
219*
220* =====================================================================
221*
222* .. Parameters ..
223 DOUBLE PRECISION ONE, ZERO
224 parameter( one = 1.0d+0, zero = 0.0d+0 )
225* ..
226* .. Local Scalars ..
227 LOGICAL BAND, GEN, NOTRAN, QRS, SYM, TRAN, TRI
228 CHARACTER C1, DIAG
229 CHARACTER*2 C2
230 INTEGER J, MB, NX
231* ..
232* .. External Functions ..
233 LOGICAL LSAME, LSAMEN
234 EXTERNAL lsame, lsamen
235* ..
236* .. External Subroutines ..
237 EXTERNAL dgbmv, dgemm, dlacpy, dlarnv, dsbmv, dspmv,
239* ..
240* .. Intrinsic Functions ..
241 INTRINSIC max
242* ..
243* .. Executable Statements ..
244*
245* Test the input parameters.
246*
247 info = 0
248 c1 = path( 1: 1 )
249 c2 = path( 2: 3 )
250 tran = lsame( trans, 'T' ) .OR. lsame( trans, 'C' )
251 notran = .NOT.tran
252 gen = lsame( path( 2: 2 ), 'G' )
253 qrs = lsame( path( 2: 2 ), 'Q' ) .OR. lsame( path( 3: 3 ), 'Q' )
254 sym = lsame( path( 2: 2 ), 'P' ) .OR. lsame( path( 2: 2 ), 'S' )
255 tri = lsame( path( 2: 2 ), 'T' )
256 band = lsame( path( 3: 3 ), 'B' )
257 IF( .NOT.lsame( c1, 'Double precision' ) ) THEN
258 info = -1
259 ELSE IF( .NOT.( lsame( xtype, 'N' ) .OR. lsame( xtype, 'C' ) ) )
260 $ THEN
261 info = -2
262 ELSE IF( ( sym .OR. tri ) .AND. .NOT.
263 $ ( lsame( uplo, 'U' ) .OR. lsame( uplo, 'L' ) ) ) THEN
264 info = -3
265 ELSE IF( ( gen .OR. qrs ) .AND. .NOT.
266 $ ( tran .OR. lsame( trans, 'N' ) ) ) THEN
267 info = -4
268 ELSE IF( m.LT.0 ) THEN
269 info = -5
270 ELSE IF( n.LT.0 ) THEN
271 info = -6
272 ELSE IF( band .AND. kl.LT.0 ) THEN
273 info = -7
274 ELSE IF( band .AND. ku.LT.0 ) THEN
275 info = -8
276 ELSE IF( nrhs.LT.0 ) THEN
277 info = -9
278 ELSE IF( ( .NOT.band .AND. lda.LT.max( 1, m ) ) .OR.
279 $ ( band .AND. ( sym .OR. tri ) .AND. lda.LT.kl+1 ) .OR.
280 $ ( band .AND. gen .AND. lda.LT.kl+ku+1 ) ) THEN
281 info = -11
282 ELSE IF( ( notran .AND. ldx.LT.max( 1, n ) ) .OR.
283 $ ( tran .AND. ldx.LT.max( 1, m ) ) ) THEN
284 info = -13
285 ELSE IF( ( notran .AND. ldb.LT.max( 1, m ) ) .OR.
286 $ ( tran .AND. ldb.LT.max( 1, n ) ) ) THEN
287 info = -15
288 END IF
289 IF( info.NE.0 ) THEN
290 CALL xerbla( 'DLARHS', -info )
291 RETURN
292 END IF
293*
294* Initialize X to NRHS random vectors unless XTYPE = 'C'.
295*
296 IF( tran ) THEN
297 nx = m
298 mb = n
299 ELSE
300 nx = n
301 mb = m
302 END IF
303 IF( .NOT.lsame( xtype, 'C' ) ) THEN
304 DO 10 j = 1, nrhs
305 CALL dlarnv( 2, iseed, n, x( 1, j ) )
306 10 CONTINUE
307 END IF
308*
309* Multiply X by op(A) using an appropriate
310* matrix multiply routine.
311*
312 IF( lsamen( 2, c2, 'GE' ) .OR. lsamen( 2, c2, 'QR' ) .OR.
313 $ lsamen( 2, c2, 'LQ' ) .OR. lsamen( 2, c2, 'QL' ) .OR.
314 $ lsamen( 2, c2, 'RQ' ) ) THEN
315*
316* General matrix
317*
318 CALL dgemm( trans, 'N', mb, nrhs, nx, one, a, lda, x, ldx,
319 $ zero, b, ldb )
320*
321 ELSE IF( lsamen( 2, c2, 'PO' ) .OR. lsamen( 2, c2, 'SY' ) ) THEN
322*
323* Symmetric matrix, 2-D storage
324*
325 CALL dsymm( 'Left', uplo, n, nrhs, one, a, lda, x, ldx, zero,
326 $ b, ldb )
327*
328 ELSE IF( lsamen( 2, c2, 'GB' ) ) THEN
329*
330* General matrix, band storage
331*
332 DO 20 j = 1, nrhs
333 CALL dgbmv( trans, mb, nx, kl, ku, one, a, lda, x( 1, j ),
334 $ 1, zero, b( 1, j ), 1 )
335 20 CONTINUE
336*
337 ELSE IF( lsamen( 2, c2, 'PB' ) ) THEN
338*
339* Symmetric matrix, band storage
340*
341 DO 30 j = 1, nrhs
342 CALL dsbmv( uplo, n, kl, one, a, lda, x( 1, j ), 1, zero,
343 $ b( 1, j ), 1 )
344 30 CONTINUE
345*
346 ELSE IF( lsamen( 2, c2, 'PP' ) .OR. lsamen( 2, c2, 'SP' ) ) THEN
347*
348* Symmetric matrix, packed storage
349*
350 DO 40 j = 1, nrhs
351 CALL dspmv( uplo, n, one, a, x( 1, j ), 1, zero, b( 1, j ),
352 $ 1 )
353 40 CONTINUE
354*
355 ELSE IF( lsamen( 2, c2, 'TR' ) ) THEN
356*
357* Triangular matrix. Note that for triangular matrices,
358* KU = 1 => non-unit triangular
359* KU = 2 => unit triangular
360*
361 CALL dlacpy( 'Full', n, nrhs, x, ldx, b, ldb )
362 IF( ku.EQ.2 ) THEN
363 diag = 'U'
364 ELSE
365 diag = 'N'
366 END IF
367 CALL dtrmm( 'Left', uplo, trans, diag, n, nrhs, one, a, lda, b,
368 $ ldb )
369*
370 ELSE IF( lsamen( 2, c2, 'TP' ) ) THEN
371*
372* Triangular matrix, packed storage
373*
374 CALL dlacpy( 'Full', n, nrhs, x, ldx, b, ldb )
375 IF( ku.EQ.2 ) THEN
376 diag = 'U'
377 ELSE
378 diag = 'N'
379 END IF
380 DO 50 j = 1, nrhs
381 CALL dtpmv( uplo, trans, diag, n, a, b( 1, j ), 1 )
382 50 CONTINUE
383*
384 ELSE IF( lsamen( 2, c2, 'TB' ) ) THEN
385*
386* Triangular matrix, banded storage
387*
388 CALL dlacpy( 'Full', n, nrhs, x, ldx, b, ldb )
389 IF( ku.EQ.2 ) THEN
390 diag = 'U'
391 ELSE
392 diag = 'N'
393 END IF
394 DO 60 j = 1, nrhs
395 CALL dtbmv( uplo, trans, diag, n, kl, a, lda, b( 1, j ), 1 )
396 60 CONTINUE
397*
398 ELSE
399*
400* If PATH is none of the above, return with an error code.
401*
402 info = -1
403 CALL xerbla( 'DLARHS', -info )
404 END IF
405*
406 RETURN
407*
408* End of DLARHS
409*
410 END
subroutine dlarhs(path, xtype, uplo, trans, m, n, kl, ku, nrhs, a, lda, x, ldx, b, ldb, iseed, info)
DLARHS
Definition dlarhs.f:205
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine dgbmv(trans, m, n, kl, ku, alpha, a, lda, x, incx, beta, y, incy)
DGBMV
Definition dgbmv.f:188
subroutine dgemm(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
DGEMM
Definition dgemm.f:188
subroutine dsbmv(uplo, n, k, alpha, a, lda, x, incx, beta, y, incy)
DSBMV
Definition dsbmv.f:184
subroutine dsymm(side, uplo, m, n, alpha, a, lda, b, ldb, beta, c, ldc)
DSYMM
Definition dsymm.f:189
subroutine dspmv(uplo, n, alpha, ap, x, incx, beta, y, incy)
DSPMV
Definition dspmv.f:147
subroutine dlacpy(uplo, m, n, a, lda, b, ldb)
DLACPY copies all or part of one two-dimensional array to another.
Definition dlacpy.f:103
subroutine dlarnv(idist, iseed, n, x)
DLARNV returns a vector of random numbers from a uniform or normal distribution.
Definition dlarnv.f:97
subroutine dtbmv(uplo, trans, diag, n, k, a, lda, x, incx)
DTBMV
Definition dtbmv.f:186
subroutine dtpmv(uplo, trans, diag, n, ap, x, incx)
DTPMV
Definition dtpmv.f:142
subroutine dtrmm(side, uplo, transa, diag, m, n, alpha, a, lda, b, ldb)
DTRMM
Definition dtrmm.f:177