159
160
161
162
163
164
165 LOGICAL TSTERR
166 INTEGER NMAX, NN, NOUT, NRHS
167 REAL THRESH
168
169
170 LOGICAL DOTYPE( * )
171 INTEGER NVAL( * )
172 REAL RWORK( * ), S( * )
173 COMPLEX A( * ), AFAC( * ), ASAV( * ), B( * ),
174 $ BSAV( * ), WORK( * ), X( * ), XACT( * )
175
176
177
178
179
180 REAL ONE, ZERO
181 parameter( one = 1.0e+0, zero = 0.0e+0 )
182 INTEGER NTYPES
183 parameter( ntypes = 9 )
184 INTEGER NTESTS
185 parameter( ntests = 6 )
186
187
188 LOGICAL EQUIL, NOFACT, PREFAC, ZEROT
189 CHARACTER DIST, EQUED, FACT, TYPE, UPLO, XTYPE
190 CHARACTER*3 PATH
191 INTEGER I, IEQUED, IFACT, IMAT, IN, INFO, IOFF, IUPLO,
192 $ IZERO, K, K1, KL, KU, LDA, MODE, N, NB, NBMIN,
193 $ NERRS, NFACT, NFAIL, NIMAT, NRUN, NT
194 REAL AINVNM, AMAX, ANORM, CNDNUM, RCOND, RCONDC,
195 $ ROLDC, SCOND
196
197
198 CHARACTER EQUEDS( 2 ), FACTS( 3 ), UPLOS( 2 )
199 INTEGER ISEED( 4 ), ISEEDY( 4 )
200 REAL RESULT( NTESTS )
201
202
203 LOGICAL LSAME
204 REAL CLANHE, SGET06
206
207
212
213
214 LOGICAL LERR, OK
215 CHARACTER*32 SRNAMT
216 INTEGER INFOT, NUNIT
217
218
219 COMMON / infoc / infot, nunit, ok, lerr
220 COMMON / srnamc / srnamt
221
222
223 INTRINSIC cmplx, max
224
225
226 DATA iseedy / 1988, 1989, 1990, 1991 /
227 DATA uplos / 'U', 'L' /
228 DATA facts / 'F', 'N', 'E' /
229 DATA equeds / 'N', 'Y' /
230
231
232
233
234
235 path( 1: 1 ) = 'Complex precision'
236 path( 2: 3 ) = 'PO'
237 nrun = 0
238 nfail = 0
239 nerrs = 0
240 DO 10 i = 1, 4
241 iseed( i ) = iseedy( i )
242 10 CONTINUE
243
244
245
246 IF( tsterr )
247 $
CALL cerrvx( path, nout )
248 infot = 0
249
250
251
252 nb = 1
253 nbmin = 2
256
257
258
259 DO 130 in = 1, nn
260 n = nval( in )
261 lda = max( n, 1 )
262 xtype = 'N'
263 nimat = ntypes
264 IF( n.LE.0 )
265 $ nimat = 1
266
267 DO 120 imat = 1, nimat
268
269
270
271 IF( .NOT.dotype( imat ) )
272 $ GO TO 120
273
274
275
276 zerot = imat.GE.3 .AND. imat.LE.5
277 IF( zerot .AND. n.LT.imat-2 )
278 $ GO TO 120
279
280
281
282 DO 110 iuplo = 1, 2
283 uplo = uplos( iuplo )
284
285
286
287
288 CALL clatb4( path, imat, n, n,
TYPE, KL, KU, ANORM, MODE,
289 $ CNDNUM, DIST )
290
291 srnamt = 'CLATMS'
292 CALL clatms( n, n, dist, iseed,
TYPE, RWORK, MODE,
293 $ CNDNUM, ANORM, KL, KU, UPLO, A, LDA, WORK,
294 $ INFO )
295
296
297
298 IF( info.NE.0 ) THEN
299 CALL alaerh( path,
'CLATMS', info, 0, uplo, n, n, -1,
300 $ -1, -1, imat, nfail, nerrs, nout )
301 GO TO 110
302 END IF
303
304
305
306
307 IF( zerot ) THEN
308 IF( imat.EQ.3 ) THEN
309 izero = 1
310 ELSE IF( imat.EQ.4 ) THEN
311 izero = n
312 ELSE
313 izero = n / 2 + 1
314 END IF
315 ioff = ( izero-1 )*lda
316
317
318
319 IF( iuplo.EQ.1 ) THEN
320 DO 20 i = 1, izero - 1
321 a( ioff+i ) = zero
322 20 CONTINUE
323 ioff = ioff + izero
324 DO 30 i = izero, n
325 a( ioff ) = zero
326 ioff = ioff + lda
327 30 CONTINUE
328 ELSE
329 ioff = izero
330 DO 40 i = 1, izero - 1
331 a( ioff ) = zero
332 ioff = ioff + lda
333 40 CONTINUE
334 ioff = ioff - izero
335 DO 50 i = izero, n
336 a( ioff+i ) = zero
337 50 CONTINUE
338 END IF
339 ELSE
340 izero = 0
341 END IF
342
343
344
345 CALL claipd( n, a, lda+1, 0 )
346
347
348
349 CALL clacpy( uplo, n, n, a, lda, asav, lda )
350
351 DO 100 iequed = 1, 2
352 equed = equeds( iequed )
353 IF( iequed.EQ.1 ) THEN
354 nfact = 3
355 ELSE
356 nfact = 1
357 END IF
358
359 DO 90 ifact = 1, nfact
360 fact = facts( ifact )
361 prefac =
lsame( fact,
'F' )
362 nofact =
lsame( fact,
'N' )
363 equil =
lsame( fact,
'E' )
364
365 IF( zerot ) THEN
366 IF( prefac )
367 $ GO TO 90
368 rcondc = zero
369
370 ELSE IF( .NOT.
lsame( fact,
'N' ) )
THEN
371
372
373
374
375
376
377 CALL clacpy( uplo, n, n, asav, lda, afac, lda )
378 IF( equil .OR. iequed.GT.1 ) THEN
379
380
381
382
383 CALL cpoequ( n, afac, lda, s, scond, amax,
384 $ info )
385 IF( info.EQ.0 .AND. n.GT.0 ) THEN
386 IF( iequed.GT.1 )
387 $ scond = zero
388
389
390
391 CALL claqhe( uplo, n, afac, lda, s, scond,
392 $ amax, equed )
393 END IF
394 END IF
395
396
397
398
399 IF( equil )
400 $ roldc = rcondc
401
402
403
404 anorm =
clanhe(
'1', uplo, n, afac, lda, rwork )
405
406
407
408 CALL cpotrf( uplo, n, afac, lda, info )
409
410
411
412 CALL clacpy( uplo, n, n, afac, lda, a, lda )
413 CALL cpotri( uplo, n, a, lda, info )
414
415
416
417 ainvnm =
clanhe(
'1', uplo, n, a, lda, rwork )
418 IF( anorm.LE.zero .OR. ainvnm.LE.zero ) THEN
419 rcondc = one
420 ELSE
421 rcondc = ( one / anorm ) / ainvnm
422 END IF
423 END IF
424
425
426
427 CALL clacpy( uplo, n, n, asav, lda, a, lda )
428
429
430
431 srnamt = 'CLARHS'
432 CALL clarhs( path, xtype, uplo,
' ', n, n, kl, ku,
433 $ nrhs, a, lda, xact, lda, b, lda,
434 $ iseed, info )
435 xtype = 'C'
436 CALL clacpy(
'Full', n, nrhs, b, lda, bsav, lda )
437
438 IF( nofact ) THEN
439
440
441
442
443
444
445 CALL clacpy( uplo, n, n, a, lda, afac, lda )
446 CALL clacpy(
'Full', n, nrhs, b, lda, x, lda )
447
448 srnamt = 'CPOSV '
449 CALL cposv( uplo, n, nrhs, afac, lda, x, lda,
450 $ info )
451
452
453
454 IF( info.NE.izero ) THEN
455 CALL alaerh( path,
'CPOSV ', info, izero,
456 $ uplo, n, n, -1, -1, nrhs, imat,
457 $ nfail, nerrs, nout )
458 GO TO 70
459 ELSE IF( info.NE.0 ) THEN
460 GO TO 70
461 END IF
462
463
464
465
466 CALL cpot01( uplo, n, a, lda, afac, lda, rwork,
467 $ result( 1 ) )
468
469
470
471 CALL clacpy(
'Full', n, nrhs, b, lda, work,
472 $ lda )
473 CALL cpot02( uplo, n, nrhs, a, lda, x, lda,
474 $ work, lda, rwork, result( 2 ) )
475
476
477
478 CALL cget04( n, nrhs, x, lda, xact, lda, rcondc,
479 $ result( 3 ) )
480 nt = 3
481
482
483
484
485 DO 60 k = 1, nt
486 IF( result( k ).GE.thresh ) THEN
487 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
488 $
CALL aladhd( nout, path )
489 WRITE( nout, fmt = 9999 )'CPOSV ', uplo,
490 $ n, imat, k, result( k )
491 nfail = nfail + 1
492 END IF
493 60 CONTINUE
494 nrun = nrun + nt
495 70 CONTINUE
496 END IF
497
498
499
500 IF( .NOT.prefac )
501 $
CALL claset( uplo, n, n, cmplx( zero ),
502 $ cmplx( zero ), afac, lda )
503 CALL claset(
'Full', n, nrhs, cmplx( zero ),
504 $ cmplx( zero ), x, lda )
505 IF( iequed.GT.1 .AND. n.GT.0 ) THEN
506
507
508
509
510 CALL claqhe( uplo, n, a, lda, s, scond, amax,
511 $ equed )
512 END IF
513
514
515
516
517 srnamt = 'CPOSVX'
518 CALL cposvx( fact, uplo, n, nrhs, a, lda, afac,
519 $ lda, equed, s, b, lda, x, lda, rcond,
520 $ rwork, rwork( nrhs+1 ), work,
521 $ rwork( 2*nrhs+1 ), info )
522
523
524
525 IF( info.NE.izero ) THEN
526 CALL alaerh( path,
'CPOSVX', info, izero,
527 $ fact // uplo, n, n, -1, -1, nrhs,
528 $ imat, nfail, nerrs, nout )
529 GO TO 90
530 END IF
531
532 IF( info.EQ.0 ) THEN
533 IF( .NOT.prefac ) THEN
534
535
536
537
538 CALL cpot01( uplo, n, a, lda, afac, lda,
539 $ rwork( 2*nrhs+1 ), result( 1 ) )
540 k1 = 1
541 ELSE
542 k1 = 2
543 END IF
544
545
546
547 CALL clacpy(
'Full', n, nrhs, bsav, lda, work,
548 $ lda )
549 CALL cpot02( uplo, n, nrhs, asav, lda, x, lda,
550 $ work, lda, rwork( 2*nrhs+1 ),
551 $ result( 2 ) )
552
553
554
555 IF( nofact .OR. ( prefac .AND.
lsame( equed,
556 $ 'N' ) ) ) THEN
557 CALL cget04( n, nrhs, x, lda, xact, lda,
558 $ rcondc, result( 3 ) )
559 ELSE
560 CALL cget04( n, nrhs, x, lda, xact, lda,
561 $ roldc, result( 3 ) )
562 END IF
563
564
565
566
567 CALL cpot05( uplo, n, nrhs, asav, lda, b, lda,
568 $ x, lda, xact, lda, rwork,
569 $ rwork( nrhs+1 ), result( 4 ) )
570 ELSE
571 k1 = 6
572 END IF
573
574
575
576
577 result( 6 ) =
sget06( rcond, rcondc )
578
579
580
581
582 DO 80 k = k1, 6
583 IF( result( k ).GE.thresh ) THEN
584 IF( nfail.EQ.0 .AND. nerrs.EQ.0 )
585 $
CALL aladhd( nout, path )
586 IF( prefac ) THEN
587 WRITE( nout, fmt = 9997 )'CPOSVX', fact,
588 $ uplo, n, equed, imat, k, result( k )
589 ELSE
590 WRITE( nout, fmt = 9998 )'CPOSVX', fact,
591 $ uplo, n, imat, k, result( k )
592 END IF
593 nfail = nfail + 1
594 END IF
595 80 CONTINUE
596 nrun = nrun + 7 - k1
597 90 CONTINUE
598 100 CONTINUE
599 110 CONTINUE
600 120 CONTINUE
601 130 CONTINUE
602
603
604
605 CALL alasvm( path, nout, nfail, nrun, nerrs )
606
607 9999 FORMAT( 1x, a, ', UPLO=''', a1, ''', N =', i5, ', type ', i1,
608 $ ', test(', i1, ')=', g12.5 )
609 9998 FORMAT( 1x, a, ', FACT=''', a1, ''', UPLO=''', a1, ''', N=', i5,
610 $ ', type ', i1, ', test(', i1, ')=', g12.5 )
611 9997 FORMAT( 1x, a, ', FACT=''', a1, ''', UPLO=''', a1, ''', N=', i5,
612 $ ', EQUED=''', a1, ''', type ', i1, ', test(', i1, ') =',
613 $ g12.5 )
614 RETURN
615
616
617
subroutine alasvm(type, nout, nfail, nrun, nerrs)
ALASVM
subroutine clarhs(path, xtype, uplo, trans, m, n, kl, ku, nrhs, a, lda, x, ldx, b, ldb, iseed, info)
CLARHS
subroutine xlaenv(ispec, nvalue)
XLAENV
subroutine aladhd(iounit, path)
ALADHD
subroutine alaerh(path, subnam, info, infoe, opts, m, n, kl, ku, n5, imat, nfail, nerrs, nout)
ALAERH
subroutine cerrvx(path, nunit)
CERRVX
subroutine cget04(n, nrhs, x, ldx, xact, ldxact, rcond, resid)
CGET04
subroutine claipd(n, a, inda, vinda)
CLAIPD
subroutine clatb4(path, imat, m, n, type, kl, ku, anorm, mode, cndnum, dist)
CLATB4
subroutine clatms(m, n, dist, iseed, sym, d, mode, cond, dmax, kl, ku, pack, a, lda, work, info)
CLATMS
subroutine cpot01(uplo, n, a, lda, afac, ldafac, rwork, resid)
CPOT01
subroutine cpot02(uplo, n, nrhs, a, lda, x, ldx, b, ldb, rwork, resid)
CPOT02
subroutine cpot05(uplo, n, nrhs, a, lda, b, ldb, x, ldx, xact, ldxact, ferr, berr, reslts)
CPOT05
subroutine clacpy(uplo, m, n, a, lda, b, ldb)
CLACPY copies all or part of one two-dimensional array to another.
real function clanhe(norm, uplo, n, a, lda, work)
CLANHE returns the value of the 1-norm, or the Frobenius norm, or the infinity norm,...
subroutine claqhe(uplo, n, a, lda, s, scond, amax, equed)
CLAQHE scales a Hermitian matrix.
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.
logical function lsame(ca, cb)
LSAME
subroutine cpoequ(n, a, lda, s, scond, amax, info)
CPOEQU
subroutine cposv(uplo, n, nrhs, a, lda, b, ldb, info)
CPOSV computes the solution to system of linear equations A * X = B for PO matrices
subroutine cposvx(fact, uplo, n, nrhs, a, lda, af, ldaf, equed, s, b, ldb, x, ldx, rcond, ferr, berr, work, rwork, info)
CPOSVX computes the solution to system of linear equations A * X = B for PO matrices
subroutine cpotrf(uplo, n, a, lda, info)
CPOTRF
subroutine cpotri(uplo, n, a, lda, info)
CPOTRI
real function sget06(rcond, rcondc)
SGET06