Dplasma
edu.emory.mathcs.jplasma.tdouble

Class Dplasma



  • public class Dplasmaextends Object
    User's API.
    • Method Summary

      Methods 
      Modifier and TypeMethod and Description
      static Stringlapack_const(int plasma_const)
      Returns LAPACK string constant that corresponds to PLASMA integer constant.
      static int[]plasma_Allocate_IPIV(int M, int N)
      This routine allocates the memory needed for the pivot array IPIV used in LU factorization.
      static double[]plasma_Allocate_L(int M, int N)
      This routine allocates the memory needed for the lower factor L used in LU factorization.
      static double[]plasma_Allocate_T(int M, int N)
      This routine allocates the memory needed for the triangular factor T used in QR and LQ factorization.
      static intplasma_DGELS(int trans, int M, int N, int NRHS, double[] A, int A_offset, int LDA, double[] T, int T_offset, double[] B, int B_offset, int LDB)
      Solves overdetermined real linear systems involving an M-by-N matrix A using a QR factorization of A.
      static intplasma_DGEQRF(int M, int N, double[] A, int A_offset, int LDA, double[] T, int T_offset)
      Computes a QR factorization of a real M-by-N matrix A: A = Q * R.
      static intplasma_DGESV(int N, int NRHS, double[] A, int A_offset, int LDA, double[] L, int L_offset, int[] IPIV, int IPIV_offset, double[] B, int B_offset, int LDB)
      Computes the solution to a real system of linear equations A * X = B, where A is an N-by-N matrix and X and B are N-by-NRHS matrices.
      static intplasma_DGETRF(int M, int N, double[] A, int A_offset, int LDA, double[] L, int L_offset, int[] IPIV, int IPIV_offset)
      Computes an LU factorization of a general M-by-N matrix A using partial pivoting with row interchanges.
      static intplasma_DGETRS(int M, int NRHS, int N, double[] A, int A_offset, int LDA, double[] L, int L_offset, int[] IPIV, int IPIV_offset, double[] B, int B_offset, int LDB)
      Solves a system of linear equations A * X = B with a general M-by-N matrix A and a M-by-NRHS B matrix using the LU factorization computed by plasma_DGETRF.
      static intplasma_DORMQR(int side, int trans, int M, int NRHS, int N, double[] A, int A_offset, int LDA, double[] T, int T_offset, double[] B, int B_offset, int LDB)
      Overwrites the general real M-by-N matrix B with
      static intplasma_DPOSV(int uplo, int N, int NRHS, double[] A, int A_offset, int LDA, double[] B, int B_offset, int LDB)
      Computes the solution to a real system of linear equations A * X = B, where A is an N-by-N symmetric positive definite matrix and X and B are N-by-NRHS matrices.
      static intplasma_DPOTRF(int uplo, int N, double[] A, int A_offset, int LDA)
      Computes the Cholesky factorization of a real symmetric positive definite matrix A.
      static intplasma_DPOTRS(int uplo, int N, int NRHS, double[] A, int A_offset, int LDA, double[] B, int B_offset, int LDB)
      Solves a system of linear equations A*X = B with a symmetric positive definite matrix A using the Cholesky factorization A = U**T*U or A = L*L**T computed by DPOTRF.
      static intplasma_DTRSM(int side, int uplo, int transA, int diag, int N, int NRHS, double[] A, int A_offset, int LDA, double[] B, int B_offset, int LDB)
      Solves one of the matrix equations op( A )*X = alpha*B, or X*op( A ) = alpha*B, where alpha is a scalar, X and B are M-by-N matrices, A is a unit, or non-unit, upper or lower triangular matrix and op( A ) is one of op( A ) = A or op( A ) = A'.
      static intplasma_DTRSMPL(int M, int NRHS, int N, double[] A, int A_offset, int LDA, double[] L, int L_offset, int[] IPIV, int IPIV_offset, double[] B, int B_offset, int LDB)
      Applies the factor L with the pivot IPIV from plasma_DGETRF to solve P*L*X = B, where L is an M-by-N matrix and X and B are N-by-NRHS matrices.
      static intplasma_Finalize()
      PLASMA completion.
      static intplasma_get_int(int param)
      Get PLASMA integer parameter
      static intplasma_Init(int M, int N, int NRHS)
      PLASMA initialization.
      static intplasma_set_int(int param, int value)
      Set PLASMA integer parameter
    • Method Detail

      • lapack_const

        public static String lapack_const(int plasma_const)
        Returns LAPACK string constant that corresponds to PLASMA integer constant.
        Parameters:
        plasma_const - PLASMA constant.
        Returns:
        LAPACK constant.
      • plasma_set_int

        public static int plasma_set_int(int param,                 int value)
        Set PLASMA integer parameter
        Parameters:
        param - PLASMA parameter
        value - the value of the parameter.
        Returns:
        Success or error code.
      • plasma_get_int

        public static int plasma_get_int(int param)
        Get PLASMA integer parameter
        Parameters:
        param - PLASMA parameter
        Returns:
        the value of the parameter or the error code
      • plasma_Init

        public static int plasma_Init(int M,              int N,              int NRHS)
        PLASMA initialization. This routine checks internal hardware constraints and arranges PLASMA internal structures to fit the hardware.
        Parameters:
        M - The number of rows.
        N - The number of columns.
        NRHS - Number of right hand sides.
        Returns:
        Success or error code.
      • plasma_Finalize

        public static int plasma_Finalize()
        PLASMA completion. This routine ends the parallel environment by joining and destroy- ing the threads. Also, it releases any internal memory allocation needed during the execution.
        Returns:
        Success or error code.
      • plasma_Allocate_T

        public static double[] plasma_Allocate_T(int M,                         int N)
        This routine allocates the memory needed for the triangular factor T used in QR and LQ factorization.
        Parameters:
        M - The number of rows.
        N - The number of columns.
        Returns:
        User's storage for T.
      • plasma_Allocate_L

        public static double[] plasma_Allocate_L(int M,                         int N)
        This routine allocates the memory needed for the lower factor L used in LU factorization.
        Parameters:
        M - The number of rows.
        N - The number of columns.
        Returns:
        User's storage for L.
      • plasma_Allocate_IPIV

        public static int[] plasma_Allocate_IPIV(int M,                         int N)
        This routine allocates the memory needed for the pivot array IPIV used in LU factorization.
        Parameters:
        M - The number of rows.
        N - The number of columns.
        Returns:
        User's storage for IPIV.
      • plasma_DPOSV

        public static int plasma_DPOSV(int uplo,               int N,               int NRHS,               double[] A,               int A_offset,               int LDA,               double[] B,               int B_offset,               int LDB)
        Computes the solution to a real system of linear equations A * X = B, where A is an N-by-N symmetric positive definite matrix and X and B are N-by-NRHS matrices. The Cholesky decomposition is used to factor A as A = U**T* U, if uplo = PlasmaUpper, or A = L * L**T, if uplo = PlasmaLower, where U is an upper triangular matrix and L is a lower triangular matrix. The factored form of A is then used to solve the system of equations A * X = B.
        Parameters:
        uplo - = PlasmaUpper: upper triangle of A is stored; = PlasmaLower: lower triangle of A is stored.
        N - The number of linear equations, i.e., the order of the matrix A. N >= 0.
        NRHS - The number of right hand sides, i.e., the number of columns of the matrix B. NRHS >= 0.
        A - An array of dimension LDA-by-N. On entry, the symmetric matrix A. If if uplo = PlasmaUpper, the leading N-by-N upper triangular part of A contains the upper triangular part of the matrix A, and the strictly lower triangular part of A is not referenced. If if uplo = PlasmaLower, the leading N-by-N lower triangular part of A contains the lower triangular part of the matrix A, and the strictly upper triangular part of A is not referenced. On exit, if return value = PLASMA_SUCCESS, the factor U or L from the Cholesky factorization A = U**T*U or A = L*L**T.
        A_offset - The index of the first element in the array A.
        LDA - The leading dimension of the array A. LDA >= max(1,N).
        B - An array of dimension LDB-by-NRHS. On entry, the N-by-NRHS right hand side matrix B. On exit, if return value = PLASMA_SUCCESS, the N-by-NRHS solution matrix X.
        B_offset - The index of the first element in the array B.
        LDB - The leading dimension of the array B. LDB >= max(1,N).
        Returns:
        Success or error code.
      • plasma_DPOTRF

        public static int plasma_DPOTRF(int uplo,                int N,                double[] A,                int A_offset,                int LDA)
        Computes the Cholesky factorization of a real symmetric positive definite matrix A. The factorization has the form A = U**T * U, if UPLO = 'U', or A = L * L**T, if UPLO = 'L', where U is an upper triangular matrix and L is lower triangular.
        Parameters:
        uplo - = PlasmaUpper: upper triangle of A is stored; = PlasmaLower: lower triangle of A is stored.
        N - The order of the matrix A. N >= 0.
        A - An array of dimension LDA-by-N. On entry, the symmetric matrix A. If uplo = PlasmaUpper, the leading N-by-N upper triangular part of A contains the upper triangular part of the matrix A, and the strictly lower triangular part of A is not referenced. If uplo = PlasmaLower, the leading N-by-N lower triangular part of A contains the lower triangular part of the matrix A, and the strictly upper triangular part of A is not referenced. On exit, if return value = PLASMA_SUCCESS, the factor U or L from the Cholesky factorization A = U**T*U or A = L*L**T.
        A_offset - The index of the first element in the array A.
        LDA - The leading dimension of the array A. LDA >= max(1,N).
        Returns:
        Success or error code.
      • plasma_DPOTRS

        public static int plasma_DPOTRS(int uplo,                int N,                int NRHS,                double[] A,                int A_offset,                int LDA,                double[] B,                int B_offset,                int LDB)
        Solves a system of linear equations A*X = B with a symmetric positive definite matrix A using the Cholesky factorization A = U**T*U or A = L*L**T computed by DPOTRF.
        Parameters:
        uplo - = PlasmaUpper: upper triangle of A is stored; = PlasmaLower: lower triangle of A is stored.
        N - The order of the matrix A. N >= 0.
        NRHS - The number of right hand sides, i.e., the number of columns of the matrix B. NRHS >= 0.
        A - An array of dimension LDA-by-N. The triangular factor U or L from the Cholesky factorization A = U**T*U or A = L*L**T, as computed by plasma_DPOTRF.
        A_offset - The index of the first element in the array A.
        LDA - The leading dimension of the array A. LDA >= max(1,N).
        B - An array of dimension LDB-by-NRHS. On entry, the right hand side matrix B. On exit, the solution matrix X.
        B_offset - The index of the first element in the array B.
        LDB - The leading dimension of the array B. LDB >= max(1,N).
        Returns:
        Success or error code.
      • plasma_DTRSM

        public static int plasma_DTRSM(int side,               int uplo,               int transA,               int diag,               int N,               int NRHS,               double[] A,               int A_offset,               int LDA,               double[] B,               int B_offset,               int LDB)
        Solves one of the matrix equations op( A )*X = alpha*B, or X*op( A ) = alpha*B, where alpha is a scalar, X and B are M-by-N matrices, A is a unit, or non-unit, upper or lower triangular matrix and op( A ) is one of op( A ) = A or op( A ) = A'. The matrix X is overwritten on B. This routine works only for side=PlasmaLeft.
        Parameters:
        side - Specifies whether op( A ) appears on the left or right of X as follows: side = PlasmaLeft op( A )*X = alpha*B. side = PlasmaRight X*op( A ) = alpha*B.
        uplo - Specifies whether the matrix A is an upper or lower triangular matrix as follows: uplo = PlasmaUpper A is an upper triangular matrix. uplo = PlasmaLower A is a lower triangular matrix.
        transA - Specifies the form of op( A ) to be used in the matrix multiplication as follows: transA = PlasmaNoTrans op( A ) = A. transA = PlasmaTrans op( A ) = A'.
        diag - Specifies whether or not A is unit triangular as follows: diag = PlasmaUnit A is assumed to be unit triangular. diag = PlasmaNonUnit A is not assumed to be unit triangular.
        N - The number of linear equations, i.e., the order of the matrix A. N >= 0.
        NRHS - The number of right hand sides, i.e., the number of columns of the matrix B. NRHS >= 0.
        A - An array of dimension LDA-by-K, where K is M when side = PlasmaLeft and is N when side = PlasmaRight. Before entry with uplo = PlasmaUpper, the leading K by K upper triangular part of the array A must contain the upper triangular matrix and the strictly lower triangular part of A is not referenced. Before entry with uplo = PlasmaLower, the leading K by K lower triangular part of the array A must contain the lower triangular matrix and the strictly upper triangular part of A is not referenced. Note that when diag = PlasmaUnit, the diagonal elements of A are not referenced either, but are assumed to be unity.
        A_offset - The index of the first element in the array A.
        LDA - The leading dimension of the array A. LDA >= max(1,N).
        B - An array of DIMENSION LDB-by-N. Before entry, the leading M-by-N part of the array B must contain the right-hand side matrix B, and on exit is overwritten by the solution matrix X.
        B_offset - The index of the first element in the array B.
        LDB - The leading dimension of the array B. LDA >= max(1,M).
        Returns:
        Success or error code.
      • plasma_DGELS

        public static int plasma_DGELS(int trans,               int M,               int N,               int NRHS,               double[] A,               int A_offset,               int LDA,               double[] T,               int T_offset,               double[] B,               int B_offset,               int LDB)
        Solves overdetermined real linear systems involving an M-by-N matrix A using a QR factorization of A. It is assumed that A has full rank. The following options are provided:
        i. If M >= N: find the least squares solution of an overdetermined system, i.e., solve the least squares problem minimize || B - A*X ||. This routine works only for trans = PlasmaNoTrans.
        Parameters:
        trans - = PlasmaNoTrans: the linear system involves A; = PlasmaTrans: the linear system involves A**T.
        M - The number of rows of the matrix A. M >= 0.
        N - The number of columns of the matrix A. N >= 0.
        NRHS - The number of right hand sides, i.e., the number of columns of the matrices B and X. NRHS >=0.
        A - An array of dimension LDA-by-N On entry, the M-by-N matrix A. On exit, A is overwritten by details of its QR factorization as returned by plasma_DGEQRF;
        A_offset - The index of the first element in the array A.
        LDA - The index of the first element in the array A.
        T - The triangular factors. This array has to be allocated by plasma_Allocate_T.
        T_offset - The index of the first element in the array T_offset.
        B - An array of dimension LDB-by-NRHS On entry, the matrix B of right hand side vectors, stored columnwise; B is M-by-NRHS if trans = PlasmaNoTrans, or N-by-NRHS if trans = PlasmaTrans. On exit, B is overwritten by the solution vectors, stored columnwise: if trans = PlasmaNoTrans and M >= N, rows 1 to N of B contain the least squares solution vectors; the residual sum of squares for the solution in each column is given by the sum of squares of elements N+1 to M in that column; if trans = PlasmaNoTrans and M < N, rows 1 to N of B contain the minimum norm solution vectors; if trans = PlasmaTrans and M >= N, rows 1 to M of B contain the minimum norm solution vectors; if trans = PlasmaTrans and M < N, rows 1 to M of B contain the least squares solution vectors; the residual sum of squares for the solution in each column is given by the sum of squares of elements M+1 to N in that column.
        B_offset - The index of the first element in the array B.
        LDB - The leading dimension of the array B. LDA >= max(1,M).
        Returns:
        Success or error code.
      • plasma_DGEQRF

        public static int plasma_DGEQRF(int M,                int N,                double[] A,                int A_offset,                int LDA,                double[] T,                int T_offset)
        Computes a QR factorization of a real M-by-N matrix A: A = Q * R.
        Parameters:
        M - The number of rows of the matrix A. M >= 0.
        N - The number of columns of the matrix A. N >= 0.
        A - An array of dimension LDA,-by-N. On entry, the M-by-N matrix A. On exit, the elements on and above the diagonal of the array contain the min(M,N)-by-N upper trapezoidal matrix R (R is upper triangular if M >= N); the elements below the diagonal, with the array TAU, represent the orthogonal matrix Q as a product of min(M,N) elementary reflectors (see Further Details).
        A_offset - The index of the first element in the array A.
        LDA - The leading dimension of the array A. LDA >= max(1,M).
        T - An array allocated by plasma_Allocate_T. The scalar factors of the elementary reflectors (see Further Details).
        T_offset - The index of the first element in the array T.
        Returns:
        Success or error code.

        Further Details

        The matrix Q is represented as a product of elementary reflectors Q = H(1) H(2) . . . H(k), where k = min(M,N). Each H(i) has the form H(i) = I - tau * v * v' where tau is a real scalar, and v is a real vector with v(1:i-1) = 0 and v(i) = 1; v(i+1:m) is stored on exit in A(i+1:m,i), and tau in T(i).

      • plasma_DORMQR

        public static int plasma_DORMQR(int side,                int trans,                int M,                int NRHS,                int N,                double[] A,                int A_offset,                int LDA,                double[] T,                int T_offset,                double[] B,                int B_offset,                int LDB)
        Overwrites the general real M-by-N matrix B with

        side = PlasmaLeft side = PlasmaRight
        trans = PlasmaNoTrans Q * B B * Q
        trans = PlasmaTrans Q**T * B B * Q**T

        where Q is a real orthogonal matrix defined as the product of k elementary reflectors Q = H(1) H(2) . . . H(k) as returned by plasma_DGEQRF. Q is of order M if side = PlasmaLeft and of order N if side = PlasmaRight. This routine works only for side = PlasmaLeft.
        Parameters:
        side - = PlasmaLeft: apply Q or Q**T from the Left; = PlasmaRight: apply Q or Q**T from the Right.
        trans - = PlasmaNoTrans: No transpose, apply Q; = PlasmaTrans: Transpose, apply Q**T.
        M - The number of rows of the matrix B. M >= 0.
        NRHS - The number of right hand sides. NRHS >= 0.
        N - The number of columns of the matrix B. N >= 0.
        A - An array of dimension LDA-by-K. The i-th column must contain the vector which defines the elementary reflector H(i), for i = 1,2,...,k, as returned by plasma_DGEQRF in the first k columns of its array argument A. A is modified by the routine but restored on exit.
        A_offset - The index of the first element in the array A.
        LDA - The leading dimension of the array A. If side = PlasmaLeft, LDA >= max(1,M); if side = PlasmaRight, LDA >= max(1,N).
        T - T(i) must contain the scalar factor of the elementary reflector H(i), as returned by plasma_DGEQRF.
        T_offset - The index of the first element in the array T.
        B - An array of dimension LDB-by-N. On entry, the M-by-N matrix B. On exit, B is overwritten by Q*B or Q**T*B or B*Q**T or B*Q.
        B_offset - The index of the first element in the array B.
        LDB - The leading dimension of the array C. LDB >= max(1,M).
        Returns:
        Success or error code.
      • plasma_DGESV

        public static int plasma_DGESV(int N,               int NRHS,               double[] A,               int A_offset,               int LDA,               double[] L,               int L_offset,               int[] IPIV,               int IPIV_offset,               double[] B,               int B_offset,               int LDB)
        Computes the solution to a real system of linear equations A * X = B, where A is an N-by-N matrix and X and B are N-by-NRHS matrices. The LU decomposition with partial pivoting and row interchanges is used to factor A as A = P * L * U, where P is a permutation matrix, L is unit lower triangular, and U is upper triangular. The factored form of A is then used to solve the system of equations A * X = B.
        Parameters:
        N - The number of linear equations, i.e., the order of the matrix A. N >= 0.
        NRHS - The number of right hand sides, i.e., the number of columns of the matrix B. NRHS >= 0.
        A - An array of dimension LDA-by-N. On entry, the N-by-N coefficient matrix A. On exit, the factors L and U from the factorization A = P*L*U; the unit diagonal elements of L are not stored.
        A_offset - The index of the first element in the array A.
        LDA - The leading dimension of the array A. LDA >= max(1,N).
        L - The lower factor L. This array has to be allocated by plasma_Allocate_L.
        L_offset - The index of the first element in the array L.
        IPIV - The pivot indices that define the permutation matrix P; row i of the matrix was interchanged with row IPIV(i). This array has to be allocated by plasma_Allocate_IPIV.
        IPIV_offset - The index of the first element in the array IPIV.
        B - An array of dimension LDB-by-NRHS. On entry, the N-by-NRHS matrix of right hand side matrix B. On exit, if return value = PLASMA_SUCCESS, the N-by-NRHS solution matrix X.
        B_offset - The index of the first element in the array B.
        LDB - The leading dimension of the array B. LDB >= max(1,N).
        Returns:
        Success or error code.
      • plasma_DGETRF

        public static int plasma_DGETRF(int M,                int N,                double[] A,                int A_offset,                int LDA,                double[] L,                int L_offset,                int[] IPIV,                int IPIV_offset)
        Computes an LU factorization of a general M-by-N matrix A using partial pivoting with row interchanges. The factorization has the form A = P * L * U where P is a permutation matrix, L is lower triangular with unit diagonal elements (lower trapezoidal if M > N), and U is upper triangular (upper trapezoidal if M < N).
        Parameters:
        M - The number of rows of the matrix A. M >= 0.
        N - The number of columns of the matrix A. N >= 0.
        A - An array of dimension LDA-by-N. On entry, the M-by-N matrix to be factored. On exit, the factors L and U from the factorization A = P*L*U; the unit diagonal elements of L are not stored.
        A_offset - The index of the first element in the array A.
        LDA - The leading dimension of the array A. LDA >= max(1,M).
        L - The lower factor L. This array has to be allocated by plasma_Allocate_L.
        L_offset - The index of the first element in the array L.
        IPIV - The pivot indices; for 1 <= i <= min(M,N), row i of the matrix was interchanged with row IPIV(i). This array has to be allocated by plasma_Allocate_IPIV.
        IPIV_offset - The index of the first element in the array IPIV.
        Returns:
        Success or error code.
      • plasma_DGETRS

        public static int plasma_DGETRS(int M,                int NRHS,                int N,                double[] A,                int A_offset,                int LDA,                double[] L,                int L_offset,                int[] IPIV,                int IPIV_offset,                double[] B,                int B_offset,                int LDB)
        Solves a system of linear equations A * X = B with a general M-by-N matrix A and a M-by-NRHS B matrix using the LU factorization computed by plasma_DGETRF.
        Parameters:
        M - The number of rows of the matrix A. M >= 0.
        NRHS - The number of right hand sides, i.e., the number of columns of the matrix B. NRHS >= 0.
        N - The number of columns of the matrix A. M >= 0.
        A - An array of dimension LDA-by-N. The factors L and U from the factorization A = P*L*U as computed by plasma_DGETRF.
        A_offset - The index of the first element in the array A.
        LDA - The leading dimension of the array A. LDA >= max(1,M).
        L - The lower factor L computed by plasma_DGETRF.
        L_offset - The index of the first element in the array L.
        IPIV - The pivot indices from plasma_DGETRF; for 1<=i<=N, row i of the matrix was interchanged with row IPIV(i).
        IPIV_offset - The index of the first element in the array IPIV.
        B - An array of dimension LDB-by-NRHS. On entry, the right hand side matrix B. On exit, the solution matrix X.
        B_offset - The index of the first element in the array B.
        LDB - The leading dimension of the array B. LDB >= max(1,N).
        Returns:
        Success or error code.
      • plasma_DTRSMPL

        public static int plasma_DTRSMPL(int M,                 int NRHS,                 int N,                 double[] A,                 int A_offset,                 int LDA,                 double[] L,                 int L_offset,                 int[] IPIV,                 int IPIV_offset,                 double[] B,                 int B_offset,                 int LDB)
        Applies the factor L with the pivot IPIV from plasma_DGETRF to solve P*L*X = B, where L is an M-by-N matrix and X and B are N-by-NRHS matrices.
        Parameters:
        M - The number of rows of the matrix A. M >= 0.
        NRHS - The number of right hand sides, i.e., the number of columns of the matrix B. NRHS >= 0.
        N - The number of columns of the matrix A. M >= 0.
        A - An array of dimension LDA-by-N. The factors L and U from the factorization A = P*L*U as computed by plasma_DGETRF.
        A_offset - The index of the first element in the array A.
        LDA - The leading dimension of the array A. LDA >= max(1,M).
        L - The lower factor L computed by plasma_DGETRF.
        L_offset - The index of the first element in the array L.
        IPIV - The pivot indices from plasma_DGETRF; for 1<=i<=N, row i of the matrix was interchanged with row IPIV(i).
        IPIV_offset - The index of the first element in the array IPIV.
        B - An array of dimension LDB-by-NRHS. On entry, the right hand side matrix B. On exit, the solution matrix X.
        B_offset - The index of the first element in the array B.
        LDB - The leading dimension of the array B. LDB >= max(1,N).
        Returns:
        Success or error code.

SCaVis 2.0 © jWork.ORG