1*> \brief \b DPPT02
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 DPPT02( UPLO, N, NRHS, A, X, LDX, B, LDB, RWORK,
12*                          RESID )
13*
14*       .. Scalar Arguments ..
15*       CHARACTER          UPLO
16*       INTEGER            LDB, LDX, N, NRHS
17*       DOUBLE PRECISION   RESID
18*       ..
19*       .. Array Arguments ..
20*       DOUBLE PRECISION   A( * ), B( LDB, * ), RWORK( * ), X( LDX, * )
21*       ..
22*
23*
24*> \par Purpose:
25*  =============
26*>
27*> \verbatim
28*>
29*> DPPT02 computes the residual in the solution of a symmetric system
30*> of linear equations  A*x = b  when packed storage is used for the
31*> coefficient matrix.  The ratio computed is
32*>
33*>    RESID = norm(B - A*X) / ( norm(A) * norm(X) * EPS),
34*>
35*> where EPS is the machine precision.
36*> \endverbatim
37*
38*  Arguments:
39*  ==========
40*
41*> \param[in] UPLO
42*> \verbatim
43*>          UPLO is CHARACTER*1
44*>          Specifies whether the upper or lower triangular part of the
45*>          symmetric matrix A is stored:
46*>          = 'U':  Upper triangular
47*>          = 'L':  Lower triangular
48*> \endverbatim
49*>
50*> \param[in] N
51*> \verbatim
52*>          N is INTEGER
53*>          The number of rows and columns of the matrix A.  N >= 0.
54*> \endverbatim
55*>
56*> \param[in] NRHS
57*> \verbatim
58*>          NRHS is INTEGER
59*>          The number of columns of B, the matrix of right hand sides.
60*>          NRHS >= 0.
61*> \endverbatim
62*>
63*> \param[in] A
64*> \verbatim
65*>          A is DOUBLE PRECISION array, dimension (N*(N+1)/2)
66*>          The original symmetric matrix A, stored as a packed
67*>          triangular matrix.
68*> \endverbatim
69*>
70*> \param[in] X
71*> \verbatim
72*>          X is DOUBLE PRECISION array, dimension (LDX,NRHS)
73*>          The computed solution vectors for the system of linear
74*>          equations.
75*> \endverbatim
76*>
77*> \param[in] LDX
78*> \verbatim
79*>          LDX is INTEGER
80*>          The leading dimension of the array X.   LDX >= max(1,N).
81*> \endverbatim
82*>
83*> \param[in,out] B
84*> \verbatim
85*>          B is DOUBLE PRECISION array, dimension (LDB,NRHS)
86*>          On entry, the right hand side vectors for the system of
87*>          linear equations.
88*>          On exit, B is overwritten with the difference B - A*X.
89*> \endverbatim
90*>
91*> \param[in] LDB
92*> \verbatim
93*>          LDB is INTEGER
94*>          The leading dimension of the array B.  LDB >= max(1,N).
95*> \endverbatim
96*>
97*> \param[out] RWORK
98*> \verbatim
99*>          RWORK is DOUBLE PRECISION array, dimension (N)
100*> \endverbatim
101*>
102*> \param[out] RESID
103*> \verbatim
104*>          RESID is DOUBLE PRECISION
105*>          The maximum over the number of right hand sides of
106*>          norm(B - A*X) / ( norm(A) * norm(X) * EPS ).
107*> \endverbatim
108*
109*  Authors:
110*  ========
111*
112*> \author Univ. of Tennessee
113*> \author Univ. of California Berkeley
114*> \author Univ. of Colorado Denver
115*> \author NAG Ltd.
116*
117*> \date November 2011
118*
119*> \ingroup double_lin
120*
121*  =====================================================================
122      SUBROUTINE DPPT02( UPLO, N, NRHS, A, X, LDX, B, LDB, RWORK,
123     $                   RESID )
124*
125*  -- LAPACK test routine (version 3.4.0) --
126*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
127*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
128*     November 2011
129*
130*     .. Scalar Arguments ..
131      CHARACTER          UPLO
132      INTEGER            LDB, LDX, N, NRHS
133      DOUBLE PRECISION   RESID
134*     ..
135*     .. Array Arguments ..
136      DOUBLE PRECISION   A( * ), B( LDB, * ), RWORK( * ), X( LDX, * )
137*     ..
138*
139*  =====================================================================
140*
141*     .. Parameters ..
142      DOUBLE PRECISION   ZERO, ONE
143      PARAMETER          ( ZERO = 0.0D+0, ONE = 1.0D+0 )
144*     ..
145*     .. Local Scalars ..
146      INTEGER            J
147      DOUBLE PRECISION   ANORM, BNORM, EPS, XNORM
148*     ..
149*     .. External Functions ..
150      DOUBLE PRECISION   DASUM, DLAMCH, DLANSP
151      EXTERNAL           DASUM, DLAMCH, DLANSP
152*     ..
153*     .. External Subroutines ..
154      EXTERNAL           DSPMV
155*     ..
156*     .. Intrinsic Functions ..
157      INTRINSIC          MAX
158*     ..
159*     .. Executable Statements ..
160*
161*     Quick exit if N = 0 or NRHS = 0.
162*
163      IF( N.LE.0 .OR. NRHS.LE.0 ) THEN
164         RESID = ZERO
165         RETURN
166      END IF
167*
168*     Exit with RESID = 1/EPS if ANORM = 0.
169*
170      EPS = DLAMCH( 'Epsilon' )
171      ANORM = DLANSP( '1', UPLO, N, A, RWORK )
172      IF( ANORM.LE.ZERO ) THEN
173         RESID = ONE / EPS
174         RETURN
175      END IF
176*
177*     Compute  B - A*X  for the matrix of right hand sides B.
178*
179      DO 10 J = 1, NRHS
180         CALL DSPMV( UPLO, N, -ONE, A, X( 1, J ), 1, ONE, B( 1, J ), 1 )
181   10 CONTINUE
182*
183*     Compute the maximum over the number of right hand sides of
184*        norm( B - A*X ) / ( norm(A) * norm(X) * EPS ) .
185*
186      RESID = ZERO
187      DO 20 J = 1, NRHS
188         BNORM = DASUM( N, B( 1, J ), 1 )
189         XNORM = DASUM( N, X( 1, J ), 1 )
190         IF( XNORM.LE.ZERO ) THEN
191            RESID = ONE / EPS
192         ELSE
193            RESID = MAX( RESID, ( ( BNORM / ANORM ) / XNORM ) / EPS )
194         END IF
195   20 CONTINUE
196*
197      RETURN
198*
199*     End of DPPT02
200*
201      END
202