Использование ускорения (CLAPACK) для решения расширенной матрицы? - PullRequest
4 голосов
/ 30 ноября 2010

Кто-нибудь знает, какую функцию / метод использовать в Accelerate (CLAPACK) для решения расширенной матрицы, такой как приведенная ниже?Ищем любой пример кода, ссылки на образцы, подсказки о том, как решить матрицу.Я просматривал документацию, но почти все связано с более сложными графическими системами, и существуют сотни, казалось бы, похожих методов.

 ____            ____
|                    |
|  4   3  -1   |  2  |
| -2   3   8   |  0  |
|  0   2   6   | -1  |
|____            ____|

1 Ответ

10 голосов
/ 30 ноября 2010
#include <Accelerate/Accelerate.h>
#include <stdio.h>

int main(int argc, char *argv[]) {

    /* Dimension of the matrix */
    __CLPK_integer n = 3;

    /* Number of right-hand side vectors to solve for */
    __CLPK_integer nrhs = 1;

    /* Note the ordering of the entries in A.  LAPACK uses "column major"
       ordering as follows:

        0 3 6
        1 4 7
        2 5 8

       This is a Fortran-ism that persists in CLAPACK.                    */
    double A[9] = {4.0, -2.0, 0.0, 3.0, 3.0, 2.0, -1.0, 8.0, 6.0 };

    /* "Leading dimension" of the matrix; most of the time, this is just the
       matrix dimension (but not always; you'll learn about this by the time
       you need to use it. */
    __CLPK_integer lda = 3;

    /* Integer array to hold information about the matrix factorization */
    __CLPK_integer ipiv[3];

    /* Right hand side to solve for */
    double x[3] = { 2.0, 0.0, -1.0 };

    /* Leading dimension of the right hand side vector */
    __CLPK_integer ldb = 3;

    /* Status variable */
    __CLPK_integer info;

    /* Solve the augmented system with a call to dgesv_.  Note that this
       routine will overwrite the contents of the array A with a factored
       form of the matrix.  If you need the original matrix, you need to
       copy it before calling dgesv_.  Note that all the scalar arguments
       are passed as pointers; this too is a Fortran-ism.                 */
    dgesv_(&n, &nrhs, A, &lda, ipiv, x, &ldb, &info);

    /* Handle error conditions */
    if (info)
        printf("Could not solve system; dgesv exited with error %d\n", info);

    /* If no error, print the result */
    else
        printf("Solution is [%f, %f, %f]\n", x[0], x[1], x[2]);

    return 0;
}

Скомпилируйте и запустите:

 scanon$ gcc solve.c -framework Accelerate
 scanon$ ./a.out 
 Solution is [-0.479167, 1.125000, -0.541667]
...