Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
#include <stddef.h></stddef.h>
Toshihiro Shimizu 890ddd
#include <stdlib.h></stdlib.h>
Toshihiro Shimizu 890ddd
#include <memory></memory>
Campbell Barton 107701
#include <cstring></cstring>
Campbell Barton 107701
Shinya Kitaoka fd19c4
#if defined(_MSC_VER) && (_MSC_VER >= 1900)
Shinya Kitaoka fd19c4
// dummy definition for linker
Shinya Kitaoka fd19c4
#include <conio.h></conio.h>
Shinya Kitaoka fd19c4
extern "C" {
Shinya Kitaoka fd19c4
void __imp__cprintf(char const *const _Format) { _cprintf(_Format); }
Shinya Kitaoka fd19c4
}
Shinya Kitaoka fd19c4
#endif
Shinya Kitaoka fd19c4
Shinya Kitaoka 120a6e
// blasint may either be common 4 bytes or extended 8 (long)...
Shinya Kitaoka 120a6e
// Replace this and REBUILD the CBLAS with extended int if needed.
Toshihiro Shimizu 890ddd
typedef int blasint;
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
extern "C" {
Toshihiro Shimizu 890ddd
#include "tlin/cblas.h"
Toshihiro Shimizu 890ddd
}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
#include "tlin/tlin_cblas_wrap.h"
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
//========================================================================================
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
/*
Toshihiro Shimizu 890ddd
  LITTLE RESUME OF BLAS NOMENCLATURE
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
  Function names behave like:
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
    cblas_[Precision][Matrix type][Operation specifier]
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
  Precision can be:
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
    S - REAL                C - COMPLEX
Shinya Kitaoka 120a6e
    D - DOUBLE PRECISION    Z - COMPLEX*16  (may not be supported by all
Shinya Kitaoka 120a6e
  machines)
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
    In current wrapper implementation only D routines will be adopted.
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
  Matrix type can be:
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
    GE - GEneral     GB - General Band
Toshihiro Shimizu 890ddd
    SY - SYmmetric   SB - Sym. Band     SP - Sym. Packed
Toshihiro Shimizu 890ddd
    HE - HErmitian   HB - Herm. Band    HP - Herm. Packed
Toshihiro Shimizu 890ddd
    TR - TRiangular  TB - Triang. Band  TP - Triang. Packed
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
    In current wrapper implementation only GE routines will be adopted.
Toshihiro Shimizu 890ddd
Shinya Kitaoka 120a6e
  Operation specifier depends on the operation. For example, mm means matrix *
Shinya Kitaoka 120a6e
  matrix, mv matrix * vector, etc..
Toshihiro Shimizu 890ddd
*/
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
//========================================================================================
Toshihiro Shimizu 890ddd
Campbell Barton 8c6c57
static void sum(int n, const double *x, double *&y) {
Shinya Kitaoka 120a6e
  /*
Shinya Kitaoka 120a6e
void cblas_daxpy(blasint n, double a, double *x, blasint incx, double *y,
Shinya Kitaoka 120a6e
blasint incy);
Toshihiro Shimizu 890ddd
Shinya Kitaoka 120a6e
NOTE:
Toshihiro Shimizu 890ddd
Shinya Kitaoka 120a6e
Returns a * x + y.    Output is stored in y.
Toshihiro Shimizu 890ddd
Shinya Kitaoka 120a6e
incx and incy are the increments in array access - ie x[incx * i] and y[incy *
Shinya Kitaoka 120a6e
j] values only are
Shinya Kitaoka 120a6e
considered (=> we'll use 1).
Shinya Kitaoka 120a6e
*/
Shinya Kitaoka 120a6e
  double *_x = const_cast<double *="">(x);</double>
Toshihiro Shimizu 890ddd
Shinya Kitaoka 120a6e
  cblas_daxpy(n, 1.0, _x, 1, y, 1);
Toshihiro Shimizu 890ddd
}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
//---------------------------------------------------------------------------------------------------
Toshihiro Shimizu 890ddd
Shinya Kitaoka 120a6e
void tlin::multiply(int rows, int cols, const double *A, const double *x,
Shinya Kitaoka 120a6e
                    double *&y) {
Shinya Kitaoka 120a6e
  /*
Shinya Kitaoka 120a6e
void cblas_dgemv(enum CBLAS_ORDER order,  enum CBLAS_TRANSPOSE trans,  blasint
Shinya Kitaoka 120a6e
m, blasint n,
Shinya Kitaoka 120a6e
           double alpha, double  *a, blasint lda,  double  *x, blasint incx,
Shinya Kitaoka 120a6e
double beta,  double  *y, blasint incy);
Toshihiro Shimizu 890ddd
Shinya Kitaoka 120a6e
NOTE:
Toshihiro Shimizu 890ddd
Shinya Kitaoka 120a6e
Returns alpha*A*x + beta*y.    Output is stored in y.
Shinya Kitaoka 120a6e
*/
Toshihiro Shimizu 890ddd
Shinya Kitaoka 120a6e
  if (!y) {
Shinya Kitaoka 120a6e
    y = (double *)malloc(rows * sizeof(double));
Shinya Kitaoka 120a6e
    memset(y, 0, rows * sizeof(double));
Shinya Kitaoka 120a6e
  }
Toshihiro Shimizu 890ddd
Shinya Kitaoka 120a6e
  double *_A = const_cast<double *="">(A);</double>
Shinya Kitaoka 120a6e
  double *_x = const_cast<double *="">(x);</double>
Toshihiro Shimizu 890ddd
Shinya Kitaoka 120a6e
  cblas_dgemv(CblasColMajor, CblasNoTrans, rows, cols, 1.0, _A, rows, _x, 1,
Shinya Kitaoka 120a6e
              1.0, y, 1);
Toshihiro Shimizu 890ddd
}
Toshihiro Shimizu 890ddd
Toshihiro Shimizu 890ddd
//---------------------------------------------------------------------------------------------------
Toshihiro Shimizu 890ddd
Shinya Kitaoka 120a6e
void tlin::multiply(const mat &A, const double *x, double *&y) {
Shinya Kitaoka 120a6e
  multiply(A.rows(), A.cols(), A.values(), x, y);
Toshihiro Shimizu 890ddd
}