我正在尝试将 QR 分解从“C 中的数值食谱”转换为使用瘦(高)矩阵。有谁知道如何做到这一点?我认为问题在于户主如何与 A 相乘,但我无法弄清楚。
void qrdcmp(float **a, int n, float *c, float *d, int *sing)
{
int i,j,k;
float scale,sigma,sum,tau;
*sing=0;
for (k = 1; k < n; k++) {
scale = 0.0;
for (i = k; i <= n; i++) scale = FMAX(scale, fabs(a[i][k]));
if (scale == 0.0) { // Singular case.
*sing = 1;
c[k] = d[k] = 0.0;
} else { // Form Qk and Qk · A.
for (i = k; i <= n; i++) a[i][k] /= scale;
for (sum = 0.0, i = k; i <= n; i++) sum += SQR(a[i][k]);
sigma = SIGN(sqrt(sum), a[k][k]);
a[k][k] += sigma;
c[k] = sigma * a[k][k];
d[k] = -scale * sigma;
for (j = k + 1; j <= n; j++) {
for (sum = 0.0, i = k; i <= n; i++) sum += a[i][k] * a[i][j];
tau = sum / c[k];
for (i = k; i <= n; i++) a[i][j] -= tau * a[i][k];
}
}
}
d[n] = a[n][n];
if (d[n] == 0.0) *sing=1;
}