1       SUBROUTINE SPTT01( N, D, E, DF, EF, WORK, RESID )
2 *
3 *  -- LAPACK test routine (version 3.1) --
4 *     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
5 *     November 2006
6 *
7 *     .. Scalar Arguments ..
8       INTEGER            N
9       REAL               RESID
10 *     ..
11 *     .. Array Arguments ..
12       REAL               D( * ), DF( * ), E( * ), EF( * ), WORK( * )
13 *     ..
14 *
15 *  Purpose
16 *  =======
17 *
18 *  SPTT01 reconstructs a tridiagonal matrix A from its L*D*L'
19 *  factorization and computes the residual
20 *     norm(L*D*L' - A) / ( n * norm(A) * EPS ),
21 *  where EPS is the machine epsilon.
22 *
23 *  Arguments
24 *  =========
25 *
26 *  N       (input) INTEGTER
27 *          The order of the matrix A.
28 *
29 *  D       (input) REAL array, dimension (N)
30 *          The n diagonal elements of the tridiagonal matrix A.
31 *
32 *  E       (input) REAL array, dimension (N-1)
33 *          The (n-1) subdiagonal elements of the tridiagonal matrix A.
34 *
35 *  DF      (input) REAL array, dimension (N)
36 *          The n diagonal elements of the factor L from the L*D*L'
37 *          factorization of A.
38 *
39 *  EF      (input) REAL array, dimension (N-1)
40 *          The (n-1) subdiagonal elements of the factor L from the
41 *          L*D*L' factorization of A.
42 *
43 *  WORK    (workspace) REAL array, dimension (2*N)
44 *
45 *  RESID   (output) REAL
46 *          norm(L*D*L' - A) / (n * norm(A) * EPS)
47 *
48 *  =====================================================================
49 *
50 *     .. Parameters ..
51       REAL               ONE, ZERO
52       PARAMETER          ( ONE = 1.0E+0, ZERO = 0.0E+0 )
53 *     ..
54 *     .. Local Scalars ..
55       INTEGER            I
56       REAL               ANORM, DE, EPS
57 *     ..
58 *     .. External Functions ..
59       REAL               SLAMCH
60       EXTERNAL           SLAMCH
61 *     ..
62 *     .. Intrinsic Functions ..
63       INTRINSIC          ABSMAX, REAL
64 *     ..
65 *     .. Executable Statements ..
66 *
67 *     Quick return if possible
68 *
69       IF( N.LE.0 ) THEN
70          RESID = ZERO
71          RETURN
72       END IF
73 *
74       EPS = SLAMCH( 'Epsilon' )
75 *
76 *     Construct the difference L*D*L' - A.
77 *
78       WORK( 1 ) = DF( 1 ) - D( 1 )
79       DO 10 I = 1, N - 1
80          DE = DF( I )*EF( I )
81          WORK( N+I ) = DE - E( I )
82          WORK( 1+I ) = DE*EF( I ) + DF( I+1 ) - D( I+1 )
83    10 CONTINUE
84 *
85 *     Compute the 1-norms of the tridiagonal matrices A and WORK.
86 *
87       IF( N.EQ.1 ) THEN
88          ANORM = D( 1 )
89          RESID = ABS( WORK( 1 ) )
90       ELSE
91          ANORM = MAX( D( 1 )+ABS( E( 1 ) ), D( N )+ABS( E( N-1 ) ) )
92          RESID = MAXABS( WORK( 1 ) )+ABS( WORK( N+1 ) ),
93      $ABS( WORK( N ) )+ABS( WORK( 2*N-1 ) ) ) 94 DO 20 I = 2, N - 1 95 ANORM = MAX( ANORM, D( I )+ABS( E( I ) )+ABS( E( I-1 ) ) ) 96 RESID = MAX( RESID, ABS( WORK( I ) )+ABS( WORK( N+I-1 ) )+ 97$              ABS( WORK( N+I ) ) )
98    20    CONTINUE
99       END IF
100 *
101 *     Compute norm(L*D*L' - A) / (n * norm(A) * EPS)
102 *
103       IF( ANORM.LE.ZERO ) THEN
104          IF( RESID.NE.ZERO )
105      \$      RESID = ONE / EPS
106       ELSE
107          RESID = ( ( RESID / REAL( N ) ) / ANORM ) / EPS
108       END IF
109 *
110       RETURN
111 *
112 *     End of SPTT01
113 *
114       END