267 SUBROUTINE dlatm5( PRTYPE, M, N, A, LDA, B, LDB, C, LDC, D, LDD,
268 $ e,
lde, f, ldf, r, ldr, l, ldl, alpha, qblcka,
277 INTEGER lda, ldb, ldc, ldd,
lde, ldf, ldl, ldr, m, n,
278 $ prtype, qblcka, qblckb
279 DOUBLE PRECISION alpha
282 DOUBLE PRECISION a( lda, * ), b( ldb, * ), c( ldc, * ),
283 $ d( ldd, * ), e(
lde, * ), f( ldf, * ),
284 $ l( ldl, * ), r( ldr, * )
290 DOUBLE PRECISION one, zero, twenty, half, two
291 parameter( one = 1.0d+0, zero = 0.0d+0, twenty = 2.0d+1,
292 $ half = 0.5d+0, two = 2.0d+0 )
296 DOUBLE PRECISION imeps, reeps
299 INTRINSIC dble, mod, sin
306 IF( prtype.EQ.1 )
THEN
312 ELSE IF( i.EQ.j-1 )
THEN
325 b( i, j ) = one - alpha
327 ELSE IF( i.EQ.j-1 )
THEN
339 r( i, j ) = ( half-sin( dble( i / j ) ) )*twenty
340 l( i, j ) = r( i, j )
344 ELSE IF( prtype.EQ.2 .OR. prtype.EQ.3 )
THEN
348 a( i, j ) = ( half-sin( dble( i ) ) )*two
349 d( i, j ) = ( half-sin( dble( i*j ) ) )*two
360 b( i, j ) = ( half-sin( dble( i+j ) ) )*two
361 e( i, j ) = ( half-sin( dble( j ) ) )*two
371 r( i, j ) = ( half-sin( dble( i*j ) ) )*twenty
372 l( i, j ) = ( half-sin( dble( i+j ) ) )*twenty
376 IF( prtype.EQ.3 )
THEN
379 DO 130 k = 1, m - 1, qblcka
380 a( k+1, k+1 ) = a( k, k )
381 a( k+1, k ) = -sin( a( k, k+1 ) )
386 DO 140 k = 1, n - 1, qblckb
387 b( k+1, k+1 ) = b( k, k )
388 b( k+1, k ) = -sin( b( k, k+1 ) )
392 ELSE IF( prtype.EQ.4 )
THEN
395 a( i, j ) = ( half-sin( dble( i*j ) ) )*twenty
396 d( i, j ) = ( half-sin( dble( i+j ) ) )*two
402 b( i, j ) = ( half-sin( dble( i+j ) ) )*twenty
403 e( i, j ) = ( half-sin( dble( i*j ) ) )*two
409 r( i, j ) = ( half-sin( dble( j / i ) ) )*twenty
410 l( i, j ) = ( half-sin( dble( i*j ) ) )*two
414 ELSE IF( prtype.GE.5 )
THEN
415 reeps = half*two*twenty / alpha
416 imeps = ( half-two ) / alpha
419 r( i, j ) = ( half-sin( dble( i*j ) ) )*alpha / twenty
420 l( i, j ) = ( half-sin( dble( i+j ) ) )*alpha / twenty
432 $ a( i, i ) = one + reeps
433 IF( mod( i, 2 ).NE.0 .AND. i.LT.m )
THEN
435 ELSE IF( i.GT.1 )
THEN
438 ELSE IF( i.LE.8 )
THEN
444 IF( mod( i, 2 ).NE.0 .AND. i.LT.m )
THEN
446 ELSE IF( i.GT.1 )
THEN
451 IF( mod( i, 2 ).NE.0 .AND. i.LT.m )
THEN
452 a( i, i+1 ) = imeps*2
453 ELSE IF( i.GT.1 )
THEN
454 a( i, i-1 ) = -imeps*2
464 $ b( i, i ) = one - reeps
465 IF( mod( i, 2 ).NE.0 .AND. i.LT.n )
THEN
467 ELSE IF( i.GT.1 )
THEN
470 ELSE IF( i.LE.8 )
THEN
476 IF( mod( i, 2 ).NE.0 .AND. i.LT.n )
THEN
477 b( i, i+1 ) = one + imeps
478 ELSE IF( i.GT.1 )
THEN
479 b( i, i-1 ) = -one - imeps
482 b( i, i ) = one - reeps
483 IF( mod( i, 2 ).NE.0 .AND. i.LT.n )
THEN
484 b( i, i+1 ) = imeps*2
485 ELSE IF( i.GT.1 )
THEN
486 b( i, i-1 ) = -imeps*2
494 CALL
dgemm(
'N',
'N', m, n, m, one, a, lda, r, ldr, zero, c, ldc )
495 CALL
dgemm(
'N',
'N', m, n, n, -one, l, ldl, b, ldb, one, c, ldc )
496 CALL
dgemm(
'N',
'N', m, n, m, one, d, ldd, r, ldr, zero, f, ldf )
497 CALL
dgemm(
'N',
'N', m, n, n, -one, l, ldl, e,
lde, one, f, ldf )