CODEWITHSUNDEEP
clinear-algebralapacklapacke
Dekakaruk
Dekakaruk

Reputation: 717

Solving linear equation with LAPACKE

I'm trying to solve some linear equation (which is symmetrical, tridiagonal and positive). I have to use LAPACKE. My code is as follows:

#include <lapacke.h>
#include <stdio.h>



void print_mtrx(double * mtrx, int n, int m)
{
    int i, j;

    for(i = 0; i < n; i++)
    {
        for(j = 0; j < m; j++)
        {
            printf("%f ", mtrx[i*m+j]);
        }
        printf("\n");
    }
    printf("\n");
}

int main()
{
    double matrix[5*5] = {
        2,  0,  0,  0,  0,
        0,  2,  0,  0,  0,
        0,  0,  2,  0,  0,
        0,  0,  0,  2,  0,
        0,  0,  0,  0,  2
    };

    double rozw[5] = {1,2,3,4,5};

    double matrix2[5*5] = {
        7,  0,  0,  0,  0,
        0,  7,  0,  0,  0,
        0,  0,  7,  0,  0,
        0,  0,  0,  7,  0,
        0,  0,  0,  0,  7
    };



    LAPACKE_dptsv(LAPACK_COL_MAJOR, 5, 5, matrix, matrix2, rozw, 5);

    print_mtrx(matrix, 5, 5);
    print_mtrx(matrix2, 5, 5);
    print_mtrx(rozw, 5, 1);

}

LAPACKE's function seems to do nothing, without any errors. The main problem is, I have no idea what do the function parameters stand for. I've searched long, but there is no real documentation. Here is what I managed to found or guess:

How can I find real meaning of these arguments? How can I make my code work?

Upvotes: 2

Views: 1572

Answers (2)

Rufflewind
Rufflewind

Reputation: 8966

When it comes to documentation for BLAS and/or LAPACK, Intel is probably the most comprehensive out there. You can look up the docs for ?ptsv, which explains what each parameter is for.

(Hint: when searching for a BLAS or LAPACK in Google, be sure to drop the s/d/c/z prefix.)

Here's the relevant snippet:

The routine solves for X the real or complex system of linear equations A*X = B, where A is an n-by-n symmetric/Hermitian positive-definite tridiagonal matrix, the columns of matrix B are individual right-hand sides, and the columns of X are the corresponding solutions.

A is factored as A = L*D*LT (real flavors) or A = L*D*LH (complex flavors), and the factored form of A is then used to solve the system of equations A*X = B.

Input Parameters

n: The order of matrix A; n ≥ 0.

nrhs: The number of right-hand sides, the number of columns in B; nrhs ≥ 0.

d: Array, dimension at least max(1, n). Contains the diagonal elements of the tridiagonal matrix A.

e, b: Arrays: e(n - 1), b(ldb,*). The array e contains the (n - 1) subdiagonal elements of A. The array b contains the matrix B whose columns are the right-hand sides for the systems of equations. The second dimension of b must be at least max(1,nrhs).

ldb: The leading dimension of b; ldb ≥ max(1, n).

Output Parameters

d: Overwritten by the n diagonal elements of the diagonal matrix D from the L*D*LT (real) / L*D*LH (complex) factorization of A.

e: Overwritten by the (n - 1) subdiagonal elements of the unit bidiagonal factor L from the factorization of A.

b: Overwritten by the solution matrix X.

info: If info = 0, the execution is successful. If info = -i, the i-th parameter had an illegal value. If info = i, the leading minor of order i (and therefore the matrix A itself) is not positive-definite, and the solution has not been computed. The factorization has not been completed unless i = n.

Upvotes: 1

Dekakaruk
Dekakaruk

Reputation: 717

So one must to look into documentation of pure LAPACK (http://www.netlib.org/lapack/explore-html/d0/dea/dptsv_8f.html#af1bd4c731915bd8755a4da8086fd79a8), and also ignore incorrect (in case of LAPACKE) remark that LDB be greater or equal max(1,N).

The correct program is as follows:

#include <lapacke.h>
#include <stdio.h>



void print_mtrx(double * mtrx, int n, int m)
{
    int i, j;

    for(i = 0; i < n; i++)
    {
        for(j = 0; j < m; j++)
        {
            printf("%f ", mtrx[i*m+j]);
        }
        printf("\n");
    }
    printf("\n");
}

int main()
{
    double diagonal[5] = {5,1,5,1,5};
    double subdiagonal[4] = {0,0,0,0};

    double solution[5] = {1,2,3,4,5};


    LAPACKE_dptsv(LAPACK_ROW_MAJOR, 5 /*size of matrix*/, 1 /*number of columns in solution*/,
                  diagonal, subdiagonal, solution, 1 /*leading dimension of solution vector*/);

    print_mtrx(solution, 5, 1);
}

Upvotes: 2

Related Questions