00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013 #include "f2c.h"
00014 #include "blaswrap.h"
00015
00016
00017
00018 struct {
00019 char srnamt[32];
00020 } srnamc_;
00021
00022 #define srnamc_1 srnamc_
00023
00024
00025
00026 static real c_b4 = -1e10f;
00027 static real c_b9 = 0.f;
00028 static real c_b14 = -1.f;
00029 static real c_b15 = 1.f;
00030
00031 int sqrt02_(integer *m, integer *n, integer *k, real *a,
00032 real *af, real *q, real *r__, integer *lda, real *tau, real *work,
00033 integer *lwork, real *rwork, real *result)
00034 {
00035
00036 integer a_dim1, a_offset, af_dim1, af_offset, q_dim1, q_offset, r_dim1,
00037 r_offset, i__1;
00038
00039
00040 int s_copy(char *, char *, ftnlen, ftnlen);
00041
00042
00043 real eps;
00044 integer info;
00045 real resid;
00046 extern int sgemm_(char *, char *, integer *, integer *,
00047 integer *, real *, real *, integer *, real *, integer *, real *,
00048 real *, integer *);
00049 real anorm;
00050 extern int ssyrk_(char *, char *, integer *, integer *,
00051 real *, real *, integer *, real *, real *, integer *);
00052 extern doublereal slamch_(char *), slange_(char *, integer *,
00053 integer *, real *, integer *, real *);
00054 extern int slacpy_(char *, integer *, integer *, real *,
00055 integer *, real *, integer *), slaset_(char *, integer *,
00056 integer *, real *, real *, real *, integer *);
00057 extern doublereal slansy_(char *, char *, integer *, real *, integer *,
00058 real *);
00059 extern int sorgqr_(integer *, integer *, integer *, real
00060 *, integer *, real *, real *, integer *, integer *);
00061
00062
00063
00064
00065
00066
00067
00068
00069
00070
00071
00072
00073
00074
00075
00076
00077
00078
00079
00080
00081
00082
00083
00084
00085
00086
00087
00088
00089
00090
00091
00092
00093
00094
00095
00096
00097
00098
00099
00100
00101
00102
00103
00104
00105
00106
00107
00108
00109
00110
00111
00112
00113
00114
00115
00116
00117
00118
00119
00120
00121
00122
00123
00124
00125
00126
00127
00128
00129
00130
00131
00132
00133
00134
00135
00136
00137
00138
00139
00140
00141
00142
00143
00144
00145
00146
00147 r_dim1 = *lda;
00148 r_offset = 1 + r_dim1;
00149 r__ -= r_offset;
00150 q_dim1 = *lda;
00151 q_offset = 1 + q_dim1;
00152 q -= q_offset;
00153 af_dim1 = *lda;
00154 af_offset = 1 + af_dim1;
00155 af -= af_offset;
00156 a_dim1 = *lda;
00157 a_offset = 1 + a_dim1;
00158 a -= a_offset;
00159 --tau;
00160 --work;
00161 --rwork;
00162 --result;
00163
00164
00165 eps = slamch_("Epsilon");
00166
00167
00168
00169 slaset_("Full", m, n, &c_b4, &c_b4, &q[q_offset], lda);
00170 i__1 = *m - 1;
00171 slacpy_("Lower", &i__1, k, &af[af_dim1 + 2], lda, &q[q_dim1 + 2], lda);
00172
00173
00174
00175 s_copy(srnamc_1.srnamt, "SORGQR", (ftnlen)32, (ftnlen)6);
00176 sorgqr_(m, n, k, &q[q_offset], lda, &tau[1], &work[1], lwork, &info);
00177
00178
00179
00180 slaset_("Full", n, k, &c_b9, &c_b9, &r__[r_offset], lda);
00181 slacpy_("Upper", n, k, &af[af_offset], lda, &r__[r_offset], lda);
00182
00183
00184
00185 sgemm_("Transpose", "No transpose", n, k, m, &c_b14, &q[q_offset], lda, &
00186 a[a_offset], lda, &c_b15, &r__[r_offset], lda);
00187
00188
00189
00190 anorm = slange_("1", m, k, &a[a_offset], lda, &rwork[1]);
00191 resid = slange_("1", n, k, &r__[r_offset], lda, &rwork[1]);
00192 if (anorm > 0.f) {
00193 result[1] = resid / (real) max(1,*m) / anorm / eps;
00194 } else {
00195 result[1] = 0.f;
00196 }
00197
00198
00199
00200 slaset_("Full", n, n, &c_b9, &c_b15, &r__[r_offset], lda);
00201 ssyrk_("Upper", "Transpose", n, m, &c_b14, &q[q_offset], lda, &c_b15, &
00202 r__[r_offset], lda);
00203
00204
00205
00206 resid = slansy_("1", "Upper", n, &r__[r_offset], lda, &rwork[1]);
00207
00208 result[2] = resid / (real) max(1,*m) / eps;
00209
00210 return 0;
00211
00212
00213
00214 }