开发者

Is it possible to solve a non-square under/over constrained matrix using Accelerate/LAPACK?

Is it possible to solve a non-square under/over constrained matrix using Accelerate/LAPACK? Such as the following two matrices. If any variables are under constrained they should equal 0 instead of being infinite.

So in the under constrained case: A, D & E would equal 0, while B, C & F equal -1.

In the over constrained case all variables wou开发者_如何学运维ld be equal to -1.

Under Constrained:

 ____                        ____
| (A) (B) (C) (D) (E) (F)        |
| -1   0   0   1   0   0   |  0  |
|  1   0   0   0  -1   0   |  0  |
|  0  -1   1   0   0   0   |  0  |
|  0   1   0   0   0  -1   |  0  |
|  0   1   0   0   0   0   | -1  |
|____                        ____|

Over Constrained:

 ____                        ____
|                                |
| -1   0   0   1   0   0   |  0  |
|  1   0   0   0  -1   0   |  0  |
|  0  -1   1   0   0   0   |  0  |
|  0   1   0   0   0  -1   |  0  |
|  0   1   0   0   0   0   | -1  |
|  0   0   1  -1   0   0   |  0  |
|  1  -1   0   0   0   0   |  0  |
|____                        ____|


Yes!

void SolveUnderdeterminedSystem() {

    __CLPK_integer m = 5;
    __CLPK_integer n = 6;
    __CLPK_integer nrhs = 1;
    double A[30] = {
        -1.0,  1.0,  0.0,  0.0,  0.0,
         0.0,  0.0, -1.0,  1.0,  1.0,
         0.0,  0.0,  1.0,  0.0,  0.0,
         1.0,  0.0,  0.0,  0.0,  0.0,
         0.0, -1.0,  0.0,  0.0,  0.0,
         0.0,  0.0,  0.0, -1.0,  0.0
    };
    __CLPK_integer lda = 5;
    double x[6] = { 0.0, 0.0, 0.0, 0.0, -1.0, 0.0 };
    __CLPK_integer ldb = 6;
    /* Need to allocate at least 2*min(m,n) workspace. */
    double work[12];
    __CLPK_integer workSize = 12;
    __CLPK_integer info;

    dgels_("N", &m, &n, &nrhs, A, &lda, x, &ldb, work, &workSize, &info);

    if (info)
        printf("Could not solve system; dgels exited with error %d\n", info);
    else
        printf("Solution is [%f, %f, %f, %f, %f, %f]\n",
               x[0], x[1], x[2], x[3], x[4], x[5]);
}

The same routine will also solve over-determined systems in the least-squares sense (the result will be a minimizer of the residual ||Ax - b||).

Note that dgels_ assumes that the matrix has full rank (i.e., rank(A) = min(m, n)). If this is not the case, you will need to use a different routine (dgelsd_) that uses an SVD factorization instead of QR.

You seem to be asking a lot of questions about LAPACK. It would be well worth your time to read the documentation.

0

上一篇:

下一篇:

精彩评论

暂无评论...
验证码 换一张
取 消

最新问答

问答排行榜