00001 /* dget04.f -- translated by f2c (version 20061008). 00002 You must link the resulting object file with libf2c: 00003 on Microsoft Windows system, link with libf2c.lib; 00004 on Linux or Unix systems, link with .../path/to/libf2c.a -lm 00005 or, if you install libf2c.a in a standard place, with -lf2c -lm 00006 -- in that order, at the end of the command line, as in 00007 cc *.o -lf2c -lm 00008 Source for libf2c is in /netlib/f2c/libf2c.zip, e.g., 00009 00010 http://www.netlib.org/f2c/libf2c.zip 00011 */ 00012 00013 #include "f2c.h" 00014 #include "blaswrap.h" 00015 00016 /* Table of constant values */ 00017 00018 static integer c__1 = 1; 00019 00020 /* Subroutine */ int dget04_(integer *n, integer *nrhs, doublereal *x, 00021 integer *ldx, doublereal *xact, integer *ldxact, doublereal *rcond, 00022 doublereal *resid) 00023 { 00024 /* System generated locals */ 00025 integer x_dim1, x_offset, xact_dim1, xact_offset, i__1, i__2; 00026 doublereal d__1, d__2, d__3; 00027 00028 /* Local variables */ 00029 integer i__, j, ix; 00030 doublereal eps, xnorm; 00031 extern doublereal dlamch_(char *); 00032 doublereal diffnm; 00033 extern integer idamax_(integer *, doublereal *, integer *); 00034 00035 00036 /* -- LAPACK test routine (version 3.1) -- */ 00037 /* Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */ 00038 /* November 2006 */ 00039 00040 /* .. Scalar Arguments .. */ 00041 /* .. */ 00042 /* .. Array Arguments .. */ 00043 /* .. */ 00044 00045 /* Purpose */ 00046 /* ======= */ 00047 00048 /* DGET04 computes the difference between a computed solution and the */ 00049 /* true solution to a system of linear equations. */ 00050 00051 /* RESID = ( norm(X-XACT) * RCOND ) / ( norm(XACT) * EPS ), */ 00052 /* where RCOND is the reciprocal of the condition number and EPS is the */ 00053 /* machine epsilon. */ 00054 00055 /* Arguments */ 00056 /* ========= */ 00057 00058 /* N (input) INTEGER */ 00059 /* The number of rows of the matrices X and XACT. N >= 0. */ 00060 00061 /* NRHS (input) INTEGER */ 00062 /* The number of columns of the matrices X and XACT. NRHS >= 0. */ 00063 00064 /* X (input) DOUBLE PRECISION array, dimension (LDX,NRHS) */ 00065 /* The computed solution vectors. Each vector is stored as a */ 00066 /* column of the matrix X. */ 00067 00068 /* LDX (input) INTEGER */ 00069 /* The leading dimension of the array X. LDX >= max(1,N). */ 00070 00071 /* XACT (input) DOUBLE PRECISION array, dimension( LDX, NRHS ) */ 00072 /* The exact solution vectors. Each vector is stored as a */ 00073 /* column of the matrix XACT. */ 00074 00075 /* LDXACT (input) INTEGER */ 00076 /* The leading dimension of the array XACT. LDXACT >= max(1,N). */ 00077 00078 /* RCOND (input) DOUBLE PRECISION */ 00079 /* The reciprocal of the condition number of the coefficient */ 00080 /* matrix in the system of equations. */ 00081 00082 /* RESID (output) DOUBLE PRECISION */ 00083 /* The maximum over the NRHS solution vectors of */ 00084 /* ( norm(X-XACT) * RCOND ) / ( norm(XACT) * EPS ) */ 00085 00086 /* ===================================================================== */ 00087 00088 /* .. Parameters .. */ 00089 /* .. */ 00090 /* .. Local Scalars .. */ 00091 /* .. */ 00092 /* .. External Functions .. */ 00093 /* .. */ 00094 /* .. Intrinsic Functions .. */ 00095 /* .. */ 00096 /* .. Executable Statements .. */ 00097 00098 /* Quick exit if N = 0 or NRHS = 0. */ 00099 00100 /* Parameter adjustments */ 00101 x_dim1 = *ldx; 00102 x_offset = 1 + x_dim1; 00103 x -= x_offset; 00104 xact_dim1 = *ldxact; 00105 xact_offset = 1 + xact_dim1; 00106 xact -= xact_offset; 00107 00108 /* Function Body */ 00109 if (*n <= 0 || *nrhs <= 0) { 00110 *resid = 0.; 00111 return 0; 00112 } 00113 00114 /* Exit with RESID = 1/EPS if RCOND is invalid. */ 00115 00116 eps = dlamch_("Epsilon"); 00117 if (*rcond < 0.) { 00118 *resid = 1. / eps; 00119 return 0; 00120 } 00121 00122 /* Compute the maximum of */ 00123 /* norm(X - XACT) / ( norm(XACT) * EPS ) */ 00124 /* over all the vectors X and XACT . */ 00125 00126 *resid = 0.; 00127 i__1 = *nrhs; 00128 for (j = 1; j <= i__1; ++j) { 00129 ix = idamax_(n, &xact[j * xact_dim1 + 1], &c__1); 00130 xnorm = (d__1 = xact[ix + j * xact_dim1], abs(d__1)); 00131 diffnm = 0.; 00132 i__2 = *n; 00133 for (i__ = 1; i__ <= i__2; ++i__) { 00134 /* Computing MAX */ 00135 d__2 = diffnm, d__3 = (d__1 = x[i__ + j * x_dim1] - xact[i__ + j * 00136 xact_dim1], abs(d__1)); 00137 diffnm = max(d__2,d__3); 00138 /* L10: */ 00139 } 00140 if (xnorm <= 0.) { 00141 if (diffnm > 0.) { 00142 *resid = 1. / eps; 00143 } 00144 } else { 00145 /* Computing MAX */ 00146 d__1 = *resid, d__2 = diffnm / xnorm * *rcond; 00147 *resid = max(d__1,d__2); 00148 } 00149 /* L20: */ 00150 } 00151 if (*resid * eps < 1.) { 00152 *resid /= eps; 00153 } 00154 00155 return 0; 00156 00157 /* End of DGET04 */ 00158 00159 } /* dget04_ */