![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b DPTT02 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 * Definition: 00009 * =========== 00010 * 00011 * SUBROUTINE DPTT02( N, NRHS, D, E, X, LDX, B, LDB, RESID ) 00012 * 00013 * .. Scalar Arguments .. 00014 * INTEGER LDB, LDX, N, NRHS 00015 * DOUBLE PRECISION RESID 00016 * .. 00017 * .. Array Arguments .. 00018 * DOUBLE PRECISION B( LDB, * ), D( * ), E( * ), X( LDX, * ) 00019 * .. 00020 * 00021 * 00022 *> \par Purpose: 00023 * ============= 00024 *> 00025 *> \verbatim 00026 *> 00027 *> DPTT02 computes the residual for the solution to a symmetric 00028 *> tridiagonal system of equations: 00029 *> RESID = norm(B - A*X) / (norm(A) * norm(X) * EPS), 00030 *> where EPS is the machine epsilon. 00031 *> \endverbatim 00032 * 00033 * Arguments: 00034 * ========== 00035 * 00036 *> \param[in] N 00037 *> \verbatim 00038 *> N is INTEGTER 00039 *> The order of the matrix A. 00040 *> \endverbatim 00041 *> 00042 *> \param[in] NRHS 00043 *> \verbatim 00044 *> NRHS is INTEGER 00045 *> The number of right hand sides, i.e., the number of columns 00046 *> of the matrices B and X. NRHS >= 0. 00047 *> \endverbatim 00048 *> 00049 *> \param[in] D 00050 *> \verbatim 00051 *> D is DOUBLE PRECISION array, dimension (N) 00052 *> The n diagonal elements of the tridiagonal matrix A. 00053 *> \endverbatim 00054 *> 00055 *> \param[in] E 00056 *> \verbatim 00057 *> E is DOUBLE PRECISION array, dimension (N-1) 00058 *> The (n-1) subdiagonal elements of the tridiagonal matrix A. 00059 *> \endverbatim 00060 *> 00061 *> \param[in] X 00062 *> \verbatim 00063 *> X is DOUBLE PRECISION array, dimension (LDX,NRHS) 00064 *> The n by nrhs matrix of solution vectors X. 00065 *> \endverbatim 00066 *> 00067 *> \param[in] LDX 00068 *> \verbatim 00069 *> LDX is INTEGER 00070 *> The leading dimension of the array X. LDX >= max(1,N). 00071 *> \endverbatim 00072 *> 00073 *> \param[in,out] B 00074 *> \verbatim 00075 *> B is DOUBLE PRECISION array, dimension (LDB,NRHS) 00076 *> On entry, the n by nrhs matrix of right hand side vectors B. 00077 *> On exit, B is overwritten with the difference B - A*X. 00078 *> \endverbatim 00079 *> 00080 *> \param[in] LDB 00081 *> \verbatim 00082 *> LDB is INTEGER 00083 *> The leading dimension of the array B. LDB >= max(1,N). 00084 *> \endverbatim 00085 *> 00086 *> \param[out] RESID 00087 *> \verbatim 00088 *> RESID is DOUBLE PRECISION 00089 *> norm(B - A*X) / (norm(A) * norm(X) * EPS) 00090 *> \endverbatim 00091 * 00092 * Authors: 00093 * ======== 00094 * 00095 *> \author Univ. of Tennessee 00096 *> \author Univ. of California Berkeley 00097 *> \author Univ. of Colorado Denver 00098 *> \author NAG Ltd. 00099 * 00100 *> \date November 2011 00101 * 00102 *> \ingroup double_lin 00103 * 00104 * ===================================================================== 00105 SUBROUTINE DPTT02( N, NRHS, D, E, X, LDX, B, LDB, RESID ) 00106 * 00107 * -- LAPACK test routine (version 3.4.0) -- 00108 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00109 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00110 * November 2011 00111 * 00112 * .. Scalar Arguments .. 00113 INTEGER LDB, LDX, N, NRHS 00114 DOUBLE PRECISION RESID 00115 * .. 00116 * .. Array Arguments .. 00117 DOUBLE PRECISION B( LDB, * ), D( * ), E( * ), X( LDX, * ) 00118 * .. 00119 * 00120 * ===================================================================== 00121 * 00122 * .. Parameters .. 00123 DOUBLE PRECISION ONE, ZERO 00124 PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 ) 00125 * .. 00126 * .. Local Scalars .. 00127 INTEGER J 00128 DOUBLE PRECISION ANORM, BNORM, EPS, XNORM 00129 * .. 00130 * .. External Functions .. 00131 DOUBLE PRECISION DASUM, DLAMCH, DLANST 00132 EXTERNAL DASUM, DLAMCH, DLANST 00133 * .. 00134 * .. Intrinsic Functions .. 00135 INTRINSIC MAX 00136 * .. 00137 * .. External Subroutines .. 00138 EXTERNAL DLAPTM 00139 * .. 00140 * .. Executable Statements .. 00141 * 00142 * Quick return if possible 00143 * 00144 IF( N.LE.0 ) THEN 00145 RESID = ZERO 00146 RETURN 00147 END IF 00148 * 00149 * Compute the 1-norm of the tridiagonal matrix A. 00150 * 00151 ANORM = DLANST( '1', N, D, E ) 00152 * 00153 * Exit with RESID = 1/EPS if ANORM = 0. 00154 * 00155 EPS = DLAMCH( 'Epsilon' ) 00156 IF( ANORM.LE.ZERO ) THEN 00157 RESID = ONE / EPS 00158 RETURN 00159 END IF 00160 * 00161 * Compute B - A*X. 00162 * 00163 CALL DLAPTM( N, NRHS, -ONE, D, E, X, LDX, ONE, B, LDB ) 00164 * 00165 * Compute the maximum over the number of right hand sides of 00166 * norm(B - A*X) / ( norm(A) * norm(X) * EPS ). 00167 * 00168 RESID = ZERO 00169 DO 10 J = 1, NRHS 00170 BNORM = DASUM( N, B( 1, J ), 1 ) 00171 XNORM = DASUM( N, X( 1, J ), 1 ) 00172 IF( XNORM.LE.ZERO ) THEN 00173 RESID = ONE / EPS 00174 ELSE 00175 RESID = MAX( RESID, ( ( BNORM / ANORM ) / XNORM ) / EPS ) 00176 END IF 00177 10 CONTINUE 00178 * 00179 RETURN 00180 * 00181 * End of DPTT02 00182 * 00183 END