module NMatrix::LAPACKE::LAPACK

Public Class Methods

lapacke_geev(p1, p2, p3, p4, p5, p6, p7, p8, p9, p10, p11, p12) click to toggle source 
static VALUE nm_lapacke_lapacke_geev(VALUE self, VALUE order, VALUE jobvl, VALUE jobvr, VALUE n, VALUE a, VALUE lda, VALUE w, VALUE wi, VALUE vl, VALUE ldvl, VALUE vr, VALUE ldvr) {
  static int (*geev_table[nm::NUM_DTYPES])(int, char, char, int, void* a, int, void* w, void* wi, void* vl, int, void* vr, int) = {
    NULL, NULL, NULL, NULL, NULL, // no integer ops
    nm::math::lapacke::lapacke_geev<float>,
    nm::math::lapacke::lapacke_geev<double>,
    nm::math::lapacke::lapacke_geev<nm::Complex64>,
    nm::math::lapacke::lapacke_geev<nm::Complex128>,
    NULL // no Ruby objects
  };

  nm::dtype_t dtype = NM_DTYPE(a);


  if (!geev_table[dtype]) {
    rb_raise(rb_eNotImpError, "this operation not yet implemented for non-BLAS dtypes");
    return Qfalse;
  } else {
    int N = FIX2INT(n);

    char JOBVL = lapack_evd_job_sym(jobvl),
         JOBVR = lapack_evd_job_sym(jobvr);

    void* A  = NM_STORAGE_DENSE(a)->elements;
    void* W = NM_STORAGE_DENSE(w)->elements;
    void* WI = wi == Qnil ? NULL : NM_STORAGE_DENSE(wi)->elements; //For complex, wi should be nil
    void* VL = JOBVL == 'V' ? NM_STORAGE_DENSE(vl)->elements : NULL;
    void* VR = JOBVR == 'V' ? NM_STORAGE_DENSE(vr)->elements : NULL;

    // Perform the actual calculation.
    int info = geev_table[dtype](blas_order_sym(order), JOBVL, JOBVR, N, A, FIX2INT(lda), W, WI, VL, FIX2INT(ldvl), VR, FIX2INT(ldvr));

    return INT2FIX(info);
  }
}
lapacke_gesdd(p1, p2, p3, p4, p5, p6, p7, p8, p9, p10, p11) click to toggle source 
static VALUE nm_lapacke_lapacke_gesdd(VALUE self, VALUE order, VALUE jobz, VALUE m, VALUE n, VALUE a, VALUE lda, VALUE s, VALUE u, VALUE ldu, VALUE vt, VALUE ldvt) {
  static int (*gesdd_table[nm::NUM_DTYPES])(int, char, int, int, void* a, int, void* s, void* u, int, void* vt, int) = {
    NULL, NULL, NULL, NULL, NULL, // no integer ops
    nm::math::lapacke::lapacke_gesdd<float,float>,
    nm::math::lapacke::lapacke_gesdd<double,double>,
    nm::math::lapacke::lapacke_gesdd<nm::Complex64,float>,
    nm::math::lapacke::lapacke_gesdd<nm::Complex128,double>,
    NULL // no Ruby objects
  };

  nm::dtype_t dtype = NM_DTYPE(a);


  if (!gesdd_table[dtype]) {
    rb_raise(rb_eNotImpError, "this operation not yet implemented for non-BLAS dtypes");
    return Qfalse;
  } else {
    int M = FIX2INT(m),
        N = FIX2INT(n);

    char JOBZ = lapack_svd_job_sym(jobz);

    int info = gesdd_table[dtype](blas_order_sym(order),JOBZ, M, N, NM_STORAGE_DENSE(a)->elements, FIX2INT(lda),
      NM_STORAGE_DENSE(s)->elements, NM_STORAGE_DENSE(u)->elements, FIX2INT(ldu), NM_STORAGE_DENSE(vt)->elements, FIX2INT(ldvt));
    return INT2FIX(info);
  }
}
lapacke_gesvd(p1, p2, p3, p4, p5, p6, p7, p8, p9, p10, p11, p12, p13) click to toggle source

xGESVD computes the singular value decomposition (SVD) of a real M-by-N matrix A, optionally computing the left and/or right singular vectors. The SVD is written

A = U * SIGMA * transpose(V)

where SIGMA is an M-by-N matrix which is zero except for its min(m,n) diagonal elements, U is an M-by-M orthogonal matrix, and V is an N-by-N orthogonal matrix. The diagonal elements of SIGMA are the singular values of A; they are real and non-negative, and are returned in descending order. The first min(m,n) columns of U and V are the left and right singular vectors of A.

Note that the routine returns V**T, not V.

static VALUE nm_lapacke_lapacke_gesvd(VALUE self, VALUE order, VALUE jobu, VALUE jobvt, VALUE m, VALUE n, VALUE a, VALUE lda, VALUE s, VALUE u, VALUE ldu, VALUE vt, VALUE ldvt, VALUE superb) {
  static int (*gesvd_table[nm::NUM_DTYPES])(int, char, char, int, int, void* a, int, void* s, void* u, int, void* vt, int, void* superb) = {
    NULL, NULL, NULL, NULL, NULL, // no integer ops
    nm::math::lapacke::lapacke_gesvd<float,float>,
    nm::math::lapacke::lapacke_gesvd<double,double>,
    nm::math::lapacke::lapacke_gesvd<nm::Complex64,float>,
    nm::math::lapacke::lapacke_gesvd<nm::Complex128,double>,
    NULL // no Ruby objects
  };

  nm::dtype_t dtype = NM_DTYPE(a);


  if (!gesvd_table[dtype]) {
    rb_raise(rb_eNotImpError, "this operation not yet implemented for non-BLAS dtypes");
    return Qfalse;
  } else {
    int M = FIX2INT(m),
        N = FIX2INT(n);

    char JOBU = lapack_svd_job_sym(jobu),
         JOBVT = lapack_svd_job_sym(jobvt);

    int info = gesvd_table[dtype](blas_order_sym(order),JOBU, JOBVT, M, N, NM_STORAGE_DENSE(a)->elements, FIX2INT(lda),
      NM_STORAGE_DENSE(s)->elements, NM_STORAGE_DENSE(u)->elements, FIX2INT(ldu), NM_STORAGE_DENSE(vt)->elements, FIX2INT(ldvt),
      NM_STORAGE_DENSE(superb)->elements);
    return INT2FIX(info);
  }
}
lapacke_getrf(p1, p2, p3, p4, p5) click to toggle source

Call any of the lapacke_xgetrf functions as directly as possible.

The ::lapacke_getrf functions (dgetrf, sgetrf, cgetrf, and zgetrf) compute 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).

This version of getrf (the LAPACKE one) differs from the CLAPACK version. The CLAPACK has different behavior for row-major matrices (the upper matrix has unit diagonals instead of the lower and it uses column permutations instead of rows).

This is the right-looking level 3 BLAS version of the algorithm.

Arguments ↑

See: www.netlib.org/lapack/double/dgetrf.f (You don't need argument 5; this is the value returned by this function.)

This function does almost no type checking. Seriously, be really careful when you call it! There's no exception handling, so you can easily crash Ruby!

Returns an array giving the pivot indices (normally these are argument #5).

static VALUE nm_lapacke_lapacke_getrf(VALUE self, VALUE order, VALUE m, VALUE n, VALUE a, VALUE lda) {
  static int (*ttable[nm::NUM_DTYPES])(const enum CBLAS_ORDER, const int m, const int n, void* a, const int lda, int* ipiv) = {
      NULL, NULL, NULL, NULL, NULL,
      nm::math::lapacke::lapacke_getrf<float>,
      nm::math::lapacke::lapacke_getrf<double>,
      nm::math::lapacke::lapacke_getrf<nm::Complex64>,
      nm::math::lapacke::lapacke_getrf<nm::Complex128>,
      NULL
  };

  int M = FIX2INT(m),
      N = FIX2INT(n);

  // Allocate the pivot index array, which is of size MIN(M, N).
  size_t ipiv_size = std::min(M,N);
  int* ipiv = NM_ALLOCA_N(int, ipiv_size);

  if (!ttable[NM_DTYPE(a)]) {
    rb_raise(nm_eDataTypeError, "this matrix operation undefined for integer matrices");
  } else {
    ttable[NM_DTYPE(a)](blas_order_sym(order), M, N, NM_STORAGE_DENSE(a)->elements, FIX2INT(lda), ipiv);
  }

  // Result will be stored in a. We return ipiv as an array.
  VALUE ipiv_array = rb_ary_new2(ipiv_size);
  for (size_t i = 0; i < ipiv_size; ++i) {
    rb_ary_store(ipiv_array, i, INT2FIX(ipiv[i]));
  }

  return ipiv_array;
}
lapacke_getri(p1, p2, p3, p4, p5) click to toggle source

Call any of the lapacke_xgetri functions as directly as possible.

This version (the LAPACKE version) differs from the CLAPACK version in terms of the input it expects (which is the output of getrf). See getrf for details.

This function does almost no type checking. Seriously, be really careful when you call it! There's no exception handling, so you can easily crash Ruby!

Returns an array giving the pivot indices (normally these are argument #5).

static VALUE nm_lapacke_lapacke_getri(VALUE self, VALUE order, VALUE n, VALUE a, VALUE lda, VALUE ipiv) {
  static int (*ttable[nm::NUM_DTYPES])(const enum CBLAS_ORDER, const int n, void* a, const int lda, const int* ipiv) = {
      NULL, NULL, NULL, NULL, NULL,
      nm::math::lapacke::lapacke_getri<float>,
      nm::math::lapacke::lapacke_getri<double>,
      nm::math::lapacke::lapacke_getri<nm::Complex64>,
      nm::math::lapacke::lapacke_getri<nm::Complex128>,
      NULL
  };

  // Allocate the C version of the pivot index array
  int* ipiv_;
  if (TYPE(ipiv) != T_ARRAY) {
    rb_raise(rb_eArgError, "ipiv must be of type Array");
  } else {
    ipiv_ = NM_ALLOCA_N(int, RARRAY_LEN(ipiv));
    for (int index = 0; index < RARRAY_LEN(ipiv); ++index) {
      ipiv_[index] = FIX2INT( RARRAY_PTR(ipiv)[index] );
    }
  }

  if (!ttable[NM_DTYPE(a)]) {
    rb_raise(rb_eNotImpError, "this operation not yet implemented for non-BLAS dtypes");
  } else {
    ttable[NM_DTYPE(a)](blas_order_sym(order), FIX2INT(n), NM_STORAGE_DENSE(a)->elements, FIX2INT(lda), ipiv_);
  }

  return a;
}
lapacke_getrs(p1, p2, p3, p4, p5, p6, p7, p8, p9) click to toggle source

Call any of the lapacke_xgetrs functions as directly as possible.

static VALUE nm_lapacke_lapacke_getrs(VALUE self, VALUE order, VALUE trans, VALUE n, VALUE nrhs, VALUE a, VALUE lda, VALUE ipiv, VALUE b, VALUE ldb) {
  static int (*ttable[nm::NUM_DTYPES])(const enum CBLAS_ORDER Order, char Trans, const int N,
                                       const int NRHS, const void* A, const int lda, const int* ipiv, void* B,
                                       const int ldb) = {
      NULL, NULL, NULL, NULL, NULL,
      nm::math::lapacke::lapacke_getrs<float>,
      nm::math::lapacke::lapacke_getrs<double>,
      nm::math::lapacke::lapacke_getrs<nm::Complex64>,
      nm::math::lapacke::lapacke_getrs<nm::Complex128>,
      NULL
  };

  // Allocate the C version of the pivot index array
  int* ipiv_;
  if (TYPE(ipiv) != T_ARRAY) {
    rb_raise(rb_eArgError, "ipiv must be of type Array");
  } else {
    ipiv_ = NM_ALLOCA_N(int, RARRAY_LEN(ipiv));
    for (int index = 0; index < RARRAY_LEN(ipiv); ++index) {
      ipiv_[index] = FIX2INT( RARRAY_PTR(ipiv)[index] );
    }
  }

  if (!ttable[NM_DTYPE(a)]) {
    rb_raise(nm_eDataTypeError, "this matrix operation undefined for integer matrices");
  } else {
    ttable[NM_DTYPE(a)](blas_order_sym(order), lapacke_transpose_sym(trans), FIX2INT(n), FIX2INT(nrhs), NM_STORAGE_DENSE(a)->elements, FIX2INT(lda),
                        ipiv_, NM_STORAGE_DENSE(b)->elements, FIX2INT(ldb));
  }

  // b is both returned and modified directly in the argument list.
  return b;
}
lapacke_potrf(p1, p2, p3, p4, p5) click to toggle source

Call any of the LAPACKE_xpotrf functions as directly as possible.

This function does almost no type checking. Seriously, be really careful when you call it! There's no exception handling, so you can easily crash Ruby!

static VALUE nm_lapacke_lapacke_potrf(VALUE self, VALUE order, VALUE uplo, VALUE n, VALUE a, VALUE lda) {

  static int (*ttable[nm::NUM_DTYPES])(const enum CBLAS_ORDER, char, const int n, void* a, const int lda) = {
      NULL, NULL, NULL, NULL, NULL,
      nm::math::lapacke::lapacke_potrf<float>,
      nm::math::lapacke::lapacke_potrf<double>,
      nm::math::lapacke::lapacke_potrf<nm::Complex64>,
      nm::math::lapacke::lapacke_potrf<nm::Complex128>,
      NULL
  };

  if (!ttable[NM_DTYPE(a)]) {
    rb_raise(rb_eNotImpError, "this operation not yet implemented for non-BLAS dtypes");
  } else {
    ttable[NM_DTYPE(a)](blas_order_sym(order), lapacke_uplo_sym(uplo), FIX2INT(n), NM_STORAGE_DENSE(a)->elements, FIX2INT(lda));
  }

  return a;
}
lapacke_potri(p1, p2, p3, p4, p5) click to toggle source

Call any of the lapacke_xpotri functions as directly as possible.

This function does almost no type checking. Seriously, be really careful when you call it! There's no exception handling, so you can easily crash Ruby!

static VALUE nm_lapacke_lapacke_potri(VALUE self, VALUE order, VALUE uplo, VALUE n, VALUE a, VALUE lda) {

  static int (*ttable[nm::NUM_DTYPES])(const enum CBLAS_ORDER, char, const int n, void* a, const int lda) = {
      NULL, NULL, NULL, NULL, NULL,
      nm::math::lapacke::lapacke_potri<float>,
      nm::math::lapacke::lapacke_potri<double>,
      nm::math::lapacke::lapacke_potri<nm::Complex64>,
      nm::math::lapacke::lapacke_potri<nm::Complex128>,
      NULL
  };

  if (!ttable[NM_DTYPE(a)]) {
    rb_raise(rb_eNotImpError, "this operation not yet implemented for non-BLAS dtypes");
  } else {
    ttable[NM_DTYPE(a)](blas_order_sym(order), lapacke_uplo_sym(uplo), FIX2INT(n), NM_STORAGE_DENSE(a)->elements, FIX2INT(lda));
  }

  return a;
}
lapacke_potrs(p1, p2, p3, p4, p5, p6, p7, p8) click to toggle source

Call any of the LAPACKE_xpotrs functions as directly as possible.

static VALUE nm_lapacke_lapacke_potrs(VALUE self, VALUE order, VALUE uplo, VALUE n, VALUE nrhs, VALUE a, VALUE lda, VALUE b, VALUE ldb) {
  static int (*ttable[nm::NUM_DTYPES])(const enum CBLAS_ORDER Order, char Uplo, const int N,
                                       const int NRHS, const void* A, const int lda, void* B, const int ldb) = {
      NULL, NULL, NULL, NULL, NULL,
      nm::math::lapacke::lapacke_potrs<float>,
      nm::math::lapacke::lapacke_potrs<double>,
      nm::math::lapacke::lapacke_potrs<nm::Complex64>,
      nm::math::lapacke::lapacke_potrs<nm::Complex128>,
      NULL
  };


  if (!ttable[NM_DTYPE(a)]) {
    rb_raise(nm_eDataTypeError, "this matrix operation undefined for integer matrices");
  } else {

    ttable[NM_DTYPE(a)](blas_order_sym(order), lapacke_uplo_sym(uplo), FIX2INT(n), FIX2INT(nrhs), NM_STORAGE_DENSE(a)->elements, FIX2INT(lda),
                        NM_STORAGE_DENSE(b)->elements, FIX2INT(ldb));
  }

  // b is both returned and modified directly in the argument list.
  return b;
}