LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
dptt01.f
Go to the documentation of this file.
1*> \brief \b DPTT01
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 DPTT01( N, D, E, DF, EF, WORK, RESID )
12*
13* .. Scalar Arguments ..
14* INTEGER N
15* DOUBLE PRECISION RESID
16* ..
17* .. Array Arguments ..
18* DOUBLE PRECISION D( * ), DF( * ), E( * ), EF( * ), WORK( * )
19* ..
20*
21*
22*> \par Purpose:
23* =============
24*>
25*> \verbatim
26*>
27*> DPTT01 reconstructs a tridiagonal matrix A from its L*D*L'
28*> factorization and computes the residual
29*> norm(L*D*L' - A) / ( n * norm(A) * EPS ),
30*> where EPS is the machine epsilon.
31*> \endverbatim
32*
33* Arguments:
34* ==========
35*
36*> \param[in] N
37*> \verbatim
38*> N is INTEGER
39*> The order of the matrix A.
40*> \endverbatim
41*>
42*> \param[in] D
43*> \verbatim
44*> D is DOUBLE PRECISION array, dimension (N)
45*> The n diagonal elements of the tridiagonal matrix A.
46*> \endverbatim
47*>
48*> \param[in] E
49*> \verbatim
50*> E is DOUBLE PRECISION array, dimension (N-1)
51*> The (n-1) subdiagonal elements of the tridiagonal matrix A.
52*> \endverbatim
53*>
54*> \param[in] DF
55*> \verbatim
56*> DF is DOUBLE PRECISION array, dimension (N)
57*> The n diagonal elements of the factor L from the L*D*L'
58*> factorization of A.
59*> \endverbatim
60*>
61*> \param[in] EF
62*> \verbatim
63*> EF is DOUBLE PRECISION array, dimension (N-1)
64*> The (n-1) subdiagonal elements of the factor L from the
65*> L*D*L' factorization of A.
66*> \endverbatim
67*>
68*> \param[out] WORK
69*> \verbatim
70*> WORK is DOUBLE PRECISION array, dimension (2*N)
71*> \endverbatim
72*>
73*> \param[out] RESID
74*> \verbatim
75*> RESID is DOUBLE PRECISION
76*> norm(L*D*L' - A) / (n * norm(A) * EPS)
77*> \endverbatim
78*
79* Authors:
80* ========
81*
82*> \author Univ. of Tennessee
83*> \author Univ. of California Berkeley
84*> \author Univ. of Colorado Denver
85*> \author NAG Ltd.
86*
87*> \ingroup double_lin
88*
89* =====================================================================
90 SUBROUTINE dptt01( N, D, E, DF, EF, WORK, RESID )
91*
92* -- LAPACK test routine --
93* -- LAPACK is a software package provided by Univ. of Tennessee, --
94* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
95*
96* .. Scalar Arguments ..
97 INTEGER N
98 DOUBLE PRECISION RESID
99* ..
100* .. Array Arguments ..
101 DOUBLE PRECISION D( * ), DF( * ), E( * ), EF( * ), WORK( * )
102* ..
103*
104* =====================================================================
105*
106* .. Parameters ..
107 DOUBLE PRECISION ONE, ZERO
108 parameter( one = 1.0d+0, zero = 0.0d+0 )
109* ..
110* .. Local Scalars ..
111 INTEGER I
112 DOUBLE PRECISION ANORM, DE, EPS
113* ..
114* .. External Functions ..
115 DOUBLE PRECISION DLAMCH
116 EXTERNAL dlamch
117* ..
118* .. Intrinsic Functions ..
119 INTRINSIC abs, dble, max
120* ..
121* .. Executable Statements ..
122*
123* Quick return if possible
124*
125 IF( n.LE.0 ) THEN
126 resid = zero
127 RETURN
128 END IF
129*
130 eps = dlamch( 'Epsilon' )
131*
132* Construct the difference L*D*L' - A.
133*
134 work( 1 ) = df( 1 ) - d( 1 )
135 DO 10 i = 1, n - 1
136 de = df( i )*ef( i )
137 work( n+i ) = de - e( i )
138 work( 1+i ) = de*ef( i ) + df( i+1 ) - d( i+1 )
139 10 CONTINUE
140*
141* Compute the 1-norms of the tridiagonal matrices A and WORK.
142*
143 IF( n.EQ.1 ) THEN
144 anorm = d( 1 )
145 resid = abs( work( 1 ) )
146 ELSE
147 anorm = max( d( 1 )+abs( e( 1 ) ), d( n )+abs( e( n-1 ) ) )
148 resid = max( abs( work( 1 ) )+abs( work( n+1 ) ),
149 $ abs( work( n ) )+abs( work( 2*n-1 ) ) )
150 DO 20 i = 2, n - 1
151 anorm = max( anorm, d( i )+abs( e( i ) )+abs( e( i-1 ) ) )
152 resid = max( resid, abs( work( i ) )+abs( work( n+i-1 ) )+
153 $ abs( work( n+i ) ) )
154 20 CONTINUE
155 END IF
156*
157* Compute norm(L*D*L' - A) / (n * norm(A) * EPS)
158*
159 IF( anorm.LE.zero ) THEN
160 IF( resid.NE.zero )
161 $ resid = one / eps
162 ELSE
163 resid = ( ( resid / dble( n ) ) / anorm ) / eps
164 END IF
165*
166 RETURN
167*
168* End of DPTT01
169*
170 END
subroutine dptt01(n, d, e, df, ef, work, resid)
DPTT01
Definition dptt01.f:91