1

我正在尝试将 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;
}
4

0 回答 0