LAPACK 3.12.1
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ slatm4()

subroutine slatm4 ( integer itype,
integer n,
integer nz1,
integer nz2,
integer isign,
real amagn,
real rcond,
real triang,
integer idist,
integer, dimension( 4 ) iseed,
real, dimension( lda, * ) a,
integer lda )

SLATM4

Purpose:
!>
!> SLATM4 generates basic square matrices, which may later be
!> multiplied by others in order to produce test matrices.  It is
!> intended mainly to be used to test the generalized eigenvalue
!> routines.
!>
!> It first generates the diagonal and (possibly) subdiagonal,
!> according to the value of ITYPE, NZ1, NZ2, ISIGN, AMAGN, and RCOND.
!> It then fills in the upper triangle with random numbers, if TRIANG is
!> non-zero.
!> 
Parameters
[in]ITYPE
!>          ITYPE is INTEGER
!>          The  of matrix on the diagonal and sub-diagonal.
!>          If ITYPE < 0, then type abs(ITYPE) is generated and then
!>             swapped end for end (A(I,J) := A'(N-J,N-I).)  See also
!>             the description of AMAGN and ISIGN.
!>
!>          Special types:
!>          = 0:  the zero matrix.
!>          = 1:  the identity.
!>          = 2:  a transposed Jordan block.
!>          = 3:  If N is odd, then a k+1 x k+1 transposed Jordan block
!>                followed by a k x k identity block, where k=(N-1)/2.
!>                If N is even, then k=(N-2)/2, and a zero diagonal entry
!>                is tacked onto the end.
!>
!>          Diagonal types.  The diagonal consists of NZ1 zeros, then
!>             k=N-NZ1-NZ2 nonzeros.  The subdiagonal is zero.  ITYPE
!>             specifies the nonzero diagonal entries as follows:
!>          = 4:  1, ..., k
!>          = 5:  1, RCOND, ..., RCOND
!>          = 6:  1, ..., 1, RCOND
!>          = 7:  1, a, a^2, ..., a^(k-1)=RCOND
!>          = 8:  1, 1-d, 1-2*d, ..., 1-(k-1)*d=RCOND
!>          = 9:  random numbers chosen from (RCOND,1)
!>          = 10: random numbers with distribution IDIST (see SLARND.)
!> 
[in]N
!>          N is INTEGER
!>          The order of the matrix.
!> 
[in]NZ1
!>          NZ1 is INTEGER
!>          If abs(ITYPE) > 3, then the first NZ1 diagonal entries will
!>          be zero.
!> 
[in]NZ2
!>          NZ2 is INTEGER
!>          If abs(ITYPE) > 3, then the last NZ2 diagonal entries will
!>          be zero.
!> 
[in]ISIGN
!>          ISIGN is INTEGER
!>          = 0: The sign of the diagonal and subdiagonal entries will
!>               be left unchanged.
!>          = 1: The diagonal and subdiagonal entries will have their
!>               sign changed at random.
!>          = 2: If ITYPE is 2 or 3, then the same as ISIGN=1.
!>               Otherwise, with probability 0.5, odd-even pairs of
!>               diagonal entries A(2*j-1,2*j-1), A(2*j,2*j) will be
!>               converted to a 2x2 block by pre- and post-multiplying
!>               by distinct random orthogonal rotations.  The remaining
!>               diagonal entries will have their sign changed at random.
!> 
[in]AMAGN
!>          AMAGN is REAL
!>          The diagonal and subdiagonal entries will be multiplied by
!>          AMAGN.
!> 
[in]RCOND
!>          RCOND is REAL
!>          If abs(ITYPE) > 4, then the smallest diagonal entry will be
!>          entry will be RCOND.  RCOND must be between 0 and 1.
!> 
[in]TRIANG
!>          TRIANG is REAL
!>          The entries above the diagonal will be random numbers with
!>          magnitude bounded by TRIANG (i.e., random numbers multiplied
!>          by TRIANG.)
!> 
[in]IDIST
!>          IDIST is INTEGER
!>          Specifies the type of distribution to be used to generate a
!>          random matrix.
!>          = 1:  UNIFORM( 0, 1 )
!>          = 2:  UNIFORM( -1, 1 )
!>          = 3:  NORMAL ( 0, 1 )
!> 
[in,out]ISEED
!>          ISEED is INTEGER array, dimension (4)
!>          On entry ISEED specifies the seed of the random number
!>          generator.  The values of ISEED are changed on exit, and can
!>          be used in the next call to SLATM4 to continue the same
!>          random number sequence.
!>          Note: ISEED(4) should be odd, for the random number generator
!>          used at present.
!> 
[out]A
!>          A is REAL array, dimension (LDA, N)
!>          Array to be computed.
!> 
[in]LDA
!>          LDA is INTEGER
!>          Leading dimension of A.  Must be at least 1 and at least N.
!> 
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 173 of file slatm4.f.

175*
176* -- LAPACK test routine --
177* -- LAPACK is a software package provided by Univ. of Tennessee, --
178* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
179*
180* .. Scalar Arguments ..
181 INTEGER IDIST, ISIGN, ITYPE, LDA, N, NZ1, NZ2
182 REAL AMAGN, RCOND, TRIANG
183* ..
184* .. Array Arguments ..
185 INTEGER ISEED( 4 )
186 REAL A( LDA, * )
187* ..
188*
189* =====================================================================
190*
191* .. Parameters ..
192 REAL ZERO, ONE, TWO
193 parameter( zero = 0.0e0, one = 1.0e0, two = 2.0e0 )
194 REAL HALF
195 parameter( half = one / two )
196* ..
197* .. Local Scalars ..
198 INTEGER I, IOFF, ISDB, ISDE, JC, JD, JR, K, KBEG, KEND,
199 $ KLEN
200 REAL ALPHA, CL, CR, SAFMIN, SL, SR, SV1, SV2, TEMP
201* ..
202* .. External Functions ..
203 REAL SLAMCH, SLARAN, SLARND
204 EXTERNAL slamch, slaran, slarnd
205* ..
206* .. External Subroutines ..
207 EXTERNAL slaset
208* ..
209* .. Intrinsic Functions ..
210 INTRINSIC abs, exp, log, max, min, mod, real, sqrt
211* ..
212* .. Executable Statements ..
213*
214 IF( n.LE.0 )
215 $ RETURN
216 CALL slaset( 'Full', n, n, zero, zero, a, lda )
217*
218* Insure a correct ISEED
219*
220 IF( mod( iseed( 4 ), 2 ).NE.1 )
221 $ iseed( 4 ) = iseed( 4 ) + 1
222*
223* Compute diagonal and subdiagonal according to ITYPE, NZ1, NZ2,
224* and RCOND
225*
226 IF( itype.NE.0 ) THEN
227 IF( abs( itype ).GE.4 ) THEN
228 kbeg = max( 1, min( n, nz1+1 ) )
229 kend = max( kbeg, min( n, n-nz2 ) )
230 klen = kend + 1 - kbeg
231 ELSE
232 kbeg = 1
233 kend = n
234 klen = n
235 END IF
236 isdb = 1
237 isde = 0
238 GO TO ( 10, 30, 50, 80, 100, 120, 140, 160,
239 $ 180, 200 )abs( itype )
240*
241* abs(ITYPE) = 1: Identity
242*
243 10 CONTINUE
244 DO 20 jd = 1, n
245 a( jd, jd ) = one
246 20 CONTINUE
247 GO TO 220
248*
249* abs(ITYPE) = 2: Transposed Jordan block
250*
251 30 CONTINUE
252 DO 40 jd = 1, n - 1
253 a( jd+1, jd ) = one
254 40 CONTINUE
255 isdb = 1
256 isde = n - 1
257 GO TO 220
258*
259* abs(ITYPE) = 3: Transposed Jordan block, followed by the
260* identity.
261*
262 50 CONTINUE
263 k = ( n-1 ) / 2
264 DO 60 jd = 1, k
265 a( jd+1, jd ) = one
266 60 CONTINUE
267 isdb = 1
268 isde = k
269 DO 70 jd = k + 2, 2*k + 1
270 a( jd, jd ) = one
271 70 CONTINUE
272 GO TO 220
273*
274* abs(ITYPE) = 4: 1,...,k
275*
276 80 CONTINUE
277 DO 90 jd = kbeg, kend
278 a( jd, jd ) = real( jd-nz1 )
279 90 CONTINUE
280 GO TO 220
281*
282* abs(ITYPE) = 5: One large D value:
283*
284 100 CONTINUE
285 DO 110 jd = kbeg + 1, kend
286 a( jd, jd ) = rcond
287 110 CONTINUE
288 a( kbeg, kbeg ) = one
289 GO TO 220
290*
291* abs(ITYPE) = 6: One small D value:
292*
293 120 CONTINUE
294 DO 130 jd = kbeg, kend - 1
295 a( jd, jd ) = one
296 130 CONTINUE
297 a( kend, kend ) = rcond
298 GO TO 220
299*
300* abs(ITYPE) = 7: Exponentially distributed D values:
301*
302 140 CONTINUE
303 a( kbeg, kbeg ) = one
304 IF( klen.GT.1 ) THEN
305 alpha = rcond**( one / real( klen-1 ) )
306 DO 150 i = 2, klen
307 a( nz1+i, nz1+i ) = alpha**real( i-1 )
308 150 CONTINUE
309 END IF
310 GO TO 220
311*
312* abs(ITYPE) = 8: Arithmetically distributed D values:
313*
314 160 CONTINUE
315 a( kbeg, kbeg ) = one
316 IF( klen.GT.1 ) THEN
317 alpha = ( one-rcond ) / real( klen-1 )
318 DO 170 i = 2, klen
319 a( nz1+i, nz1+i ) = real( klen-i )*alpha + rcond
320 170 CONTINUE
321 END IF
322 GO TO 220
323*
324* abs(ITYPE) = 9: Randomly distributed D values on ( RCOND, 1):
325*
326 180 CONTINUE
327 alpha = log( rcond )
328 DO 190 jd = kbeg, kend
329 a( jd, jd ) = exp( alpha*slaran( iseed ) )
330 190 CONTINUE
331 GO TO 220
332*
333* abs(ITYPE) = 10: Randomly distributed D values from DIST
334*
335 200 CONTINUE
336 DO 210 jd = kbeg, kend
337 a( jd, jd ) = slarnd( idist, iseed )
338 210 CONTINUE
339*
340 220 CONTINUE
341*
342* Scale by AMAGN
343*
344 DO 230 jd = kbeg, kend
345 a( jd, jd ) = amagn*real( a( jd, jd ) )
346 230 CONTINUE
347 DO 240 jd = isdb, isde
348 a( jd+1, jd ) = amagn*real( a( jd+1, jd ) )
349 240 CONTINUE
350*
351* If ISIGN = 1 or 2, assign random signs to diagonal and
352* subdiagonal
353*
354 IF( isign.GT.0 ) THEN
355 DO 250 jd = kbeg, kend
356 IF( real( a( jd, jd ) ).NE.zero ) THEN
357 IF( slaran( iseed ).GT.half )
358 $ a( jd, jd ) = -a( jd, jd )
359 END IF
360 250 CONTINUE
361 DO 260 jd = isdb, isde
362 IF( real( a( jd+1, jd ) ).NE.zero ) THEN
363 IF( slaran( iseed ).GT.half )
364 $ a( jd+1, jd ) = -a( jd+1, jd )
365 END IF
366 260 CONTINUE
367 END IF
368*
369* Reverse if ITYPE < 0
370*
371 IF( itype.LT.0 ) THEN
372 DO 270 jd = kbeg, ( kbeg+kend-1 ) / 2
373 temp = a( jd, jd )
374 a( jd, jd ) = a( kbeg+kend-jd, kbeg+kend-jd )
375 a( kbeg+kend-jd, kbeg+kend-jd ) = temp
376 270 CONTINUE
377 DO 280 jd = 1, ( n-1 ) / 2
378 temp = a( jd+1, jd )
379 a( jd+1, jd ) = a( n+1-jd, n-jd )
380 a( n+1-jd, n-jd ) = temp
381 280 CONTINUE
382 END IF
383*
384* If ISIGN = 2, and no subdiagonals already, then apply
385* random rotations to make 2x2 blocks.
386*
387 IF( isign.EQ.2 .AND. itype.NE.2 .AND. itype.NE.3 ) THEN
388 safmin = slamch( 'S' )
389 DO 290 jd = kbeg, kend - 1, 2
390 IF( slaran( iseed ).GT.half ) THEN
391*
392* Rotation on left.
393*
394 cl = two*slaran( iseed ) - one
395 sl = two*slaran( iseed ) - one
396 temp = one / max( safmin, sqrt( cl**2+sl**2 ) )
397 cl = cl*temp
398 sl = sl*temp
399*
400* Rotation on right.
401*
402 cr = two*slaran( iseed ) - one
403 sr = two*slaran( iseed ) - one
404 temp = one / max( safmin, sqrt( cr**2+sr**2 ) )
405 cr = cr*temp
406 sr = sr*temp
407*
408* Apply
409*
410 sv1 = a( jd, jd )
411 sv2 = a( jd+1, jd+1 )
412 a( jd, jd ) = cl*cr*sv1 + sl*sr*sv2
413 a( jd+1, jd ) = -sl*cr*sv1 + cl*sr*sv2
414 a( jd, jd+1 ) = -cl*sr*sv1 + sl*cr*sv2
415 a( jd+1, jd+1 ) = sl*sr*sv1 + cl*cr*sv2
416 END IF
417 290 CONTINUE
418 END IF
419*
420 END IF
421*
422* Fill in upper triangle (except for 2x2 blocks)
423*
424 IF( triang.NE.zero ) THEN
425 IF( isign.NE.2 .OR. itype.EQ.2 .OR. itype.EQ.3 ) THEN
426 ioff = 1
427 ELSE
428 ioff = 2
429 DO 300 jr = 1, n - 1
430 IF( a( jr+1, jr ).EQ.zero )
431 $ a( jr, jr+1 ) = triang*slarnd( idist, iseed )
432 300 CONTINUE
433 END IF
434*
435 DO 320 jc = 2, n
436 DO 310 jr = 1, jc - ioff
437 a( jr, jc ) = triang*slarnd( idist, iseed )
438 310 CONTINUE
439 320 CONTINUE
440 END IF
441*
442 RETURN
443*
444* End of SLATM4
445*
real function slamch(cmach)
SLAMCH
Definition slamch.f:68
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:108
real function slaran(iseed)
SLARAN
Definition slaran.f:67
real function slarnd(idist, iseed)
SLARND
Definition slarnd.f:73
Here is the call graph for this function:
Here is the caller graph for this function: