Head

GCC Code Coverage Report

Directory:	.		Exec	Total	Coverage
File:	frame/base/blas_lapack.cpp	Lines:	0	143	0.0 %
Date:	2019-01-14	Branches:	0	24	0.0 %


/**
 *  HMLP (High-Performance Machine Learning Primitives)
 *
 *  Copyright (C) 2014-2017, The University of Texas at Austin
 *
 *  This program is free software: you can redistribute it and/or modify
 *  it under the terms of the GNU General Public License as published by
 *  the Free Software Foundation, either version 3 of the License, or
 *  (at your option) any later version.
 *
 *  This program is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
 *  GNU General Public License for more details.
 *
 *  You should have received a copy of the GNU General Public License
 *  along with this program. If not, see the LICENSE file.
 *
 **/

#include <blas_lapack.hpp>

// #define DEBUG_XGEMM 1

#ifndef USE_BLAS
#warning BLAS/LAPACK routines are not compiled (-HMLP_USE_BLAS=false)
#endif

extern "C"
{
#include <external/blas_lapack_prototypes.h>
}; /** end extern "C" */






namespace hmlp
{

/**
 *  BLAS level-1 wrappers: DOT, NRM2
 */


/**
 *  @brief DDOT wrapper
 */
double xdot( int n, const double *dx, int incx, const double *dy, int incy )
{
  double ret_val;
#ifdef USE_BLAS
  ret_val = ddot_( &n, dx, &incx, dy, &incy );
#else
  printf( "xdot must enables USE_BLAS.\n" );
  exit( 1 );
#endif
  return ret_val;
}; /** end xdot() */



/**
 *  @brief SDOT wrapper
 */
float xdot( int n, const float *dx, int incx, const float *dy, int incy )
{
  float ret_val;
#ifdef USE_BLAS
  ret_val = sdot_( &n, dx, &incx, dy, &incy );
#else
  printf( "xdot must enables USE_BLAS.\n" );
  exit( 1 );
#endif
  return ret_val;
}; /** end xdot() */


/**
 *  @brief DNRM2 wrapper
 */
double xnrm2( int n, double *x, int incx )
{
  double ret_val;
#ifdef USE_BLAS
  ret_val = dnrm2_( &n, x, &incx );
#else
  printf( "xnrm2 must enables USE_BLAS.\n" );
  exit( 1 );
#endif
  return ret_val;
}; /** end xnrm2() */


/**
 *  @brief SNRM2 wrapper
 */
float  xnrm2( int n,  float *x, int incx )
{
  float ret_val;
#ifdef USE_BLAS
  ret_val = snrm2_( &n, x, &incx );
#else
  printf( "xnrm2 must enables USE_BLAS.\n" );
  exit( 1 );
#endif
  return ret_val;
}; /** end xnrm2() */











/**
 *  BLAS level-3 wrappers: GEMM, TRSM
 */


/**
 *  @brief DGEMM wrapper
 */
void xgemm
(
  const char *transA, const char *transB,
  int m, int n, int k,
  double alpha, const double *A, int lda,
                const double *B, int ldb,
  double beta,        double *C, int ldc
)
{
  double beg, xgemm_time = 0.0;
  double gflops = (double)( m ) * n * ( 2 * k ) / 1E+9;
  beg = omp_get_wtime();
#ifdef USE_BLAS
  dgemm_
  (
    transA, transB,
    &m, &n, &k,
    &alpha, A, &lda,
            B, &ldb,
    &beta,  C, &ldc
  );
#else
  //printf( "xgemm: configure HMLP_USE_BLAS=true to enable BLAS support.\n" );
  for ( int p = 0; p < k; p ++ )
  {
    for ( int j = 0; j < n; j ++ )
    {
      for ( int i = 0; i < m; i ++ )
      {
        double a, b;
        if ( *transA == 'T' ) a = A[ i * lda + p ];
        else                  a = A[ p * lda + i ];
        if ( *transB == 'T' ) b = B[ p * ldb + j ];
        else                  b = B[ j * ldb + p ];
        if ( p == 0 )
        {
          C[ j * ldc + i ] = beta * C[ j * ldc + i ] + alpha * a * b;
        }
        else
        {
          C[ j * ldc + i ] += alpha * a * b;
        }
      }
    }
  }
#endif
  xgemm_time = omp_get_wtime() - beg;
#ifdef DEBUG_XGEMM
  printf( "dgemm %s%s m %d n %d k %d, %5.2lf GFLOPS %5.2lf s\n",
      transA, transB, m, n, k, gflops / xgemm_time, xgemm_time );
#endif

#ifdef DEBUG_XGEMM
  printf( "hmlp::xgemm debug\n" );
  for ( int i = 0; i < m; i ++ )
  {
    for ( int j = 0; j < n; j ++ )
    {
      printf( "%E ", C[ j * ldc + i ] / alpha );
    }
    printf( "\n" );
  }
#endif
};

/**
 *  @brief SGEMM wrapper
 */
void xgemm
(
  const char *transA, const char *transB,
  int m, int n, int k,
  float alpha, const float *A, int lda,
               const float *B, int ldb,
  float beta,        float *C, int ldc
)
{
  double beg, xgemm_time = 0.0;
  double gflops = (double)( m ) * n * ( 2 * k ) / 1E+9;
  beg = omp_get_wtime();
#ifdef USE_BLAS
  sgemm_
  (
   transA, transB,
   &m, &n, &k,
   &alpha, A, &lda,
           B, &ldb,
   &beta,  C, &ldc
  );
#else
  //printf( "xgemm: configure HMLP_USE_BLAS=true to enable BLAS support.\n" );
  for ( int p = 0; p < k; p ++ )
  {
    for ( int j = 0; j < n; j ++ )
    {
      for ( int i = 0; i < m; i ++ )
      {
        double a, b;
        if ( *transA == 'T' ) a = A[ i * lda + p ];
        else                  a = A[ p * lda + i ];
        if ( *transB == 'T' ) b = B[ p * ldb + j ];
        else                  b = B[ j * ldb + p ];
        if ( p == 0 )
        {
          C[ j * ldc + i ] = beta * C[ j * ldc + i ] + alpha * a * b;
        }
        else
        {
          C[ j * ldc + i ] += alpha * a * b;
        }
      }
    }
  }
#endif
  xgemm_time = omp_get_wtime() - beg;
#ifdef DEBUG_XGEMM
  printf( "sgemm %s%s m %d n %d k %d, %5.2lf GFLOPS %5.2lf s\n",
      transA, transB, m, n, k, gflops / xgemm_time, xgemm_time );
#endif
#ifdef DEBUG_XGEMM
  printf( "hmlp::xgemm debug\n" );
  for ( int i = 0; i < m; i ++ )
  {
    for ( int j = 0; j < n; j ++ )
    {
      printf( "%E ", C[ j * ldc + i ] );
    }
    printf( "\n" );
  }
#endif
};


void xsyrk
(
  const char *uplo, const char *trans,
  int n, int k,
  double alpha, double *A, int lda,
  double beta,  double *C, int ldc
)
{
#ifdef USE_BLAS
  dsyrk_
  (
   uplo, trans,
   &n, &k,
   &alpha, A, &lda,
   &beta,  C, &ldc
  );
#else
  printf( "xsyrk must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xsyrk() */


void xsyrk
(
  const char *uplo, const char *trans,
  int n, int k,
  float alpha, float *A, int lda,
  float beta,  float *C, int ldc
)
{
#ifdef USE_BLAS
  ssyrk_
  (
   uplo, trans,
   &n, &k,
   &alpha, A, &lda,
   &beta,  C, &ldc
  );
#else
  printf( "xsyrk must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xsyrk() */





/**
 *  @brief DTRSM wrapper
 */
void xtrsm
(
  const char *side, const char *uplo,
  const char *transA, const char *diag,
  int m, int n,
  double alpha,
  double *A, int lda,
  double *B, int ldb
)
{
  double beg, xtrsm_time = 0.0;
  double gflops = (double)( m ) * ( m - 1 ) * n / 1E+9;
  beg = omp_get_wtime();

#ifdef USE_BLAS
  dtrsm_
  (
    side, uplo,
	  transA, diag,
	  &m, &n,
	  &alpha,
	  A, &lda,
	  B, &ldb
  );
#else
  printf( "xtrsm must enables USE_BLAS.\n" );
  exit( 1 );
#endif

  xtrsm_time = omp_get_wtime() - beg;
#ifdef DEBUG_XTRSM
  printf( "dtrsm m %d n %d, %5.2lf GFLOPS, %5.2lf s\n",
      m, n, gflops / xtrsm_time, xtrsm_time );
#endif
}; /** end xtrsm() */


/**
 *  @brief STRSM wrapper
 */
void xtrsm
(
  const char *side, const char *uplo,
  const char *transA, const char *diag,
  int m, int n,
  float alpha,
  float *A, int lda,
  float *B, int ldb
)
{
#ifdef USE_BLAS
  strsm_
  (
    side, uplo,
	  transA, diag,
	  &m, &n,
	  &alpha,
	  A, &lda,
	  B, &ldb
  );
#else
  printf( "xtrsm must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xtrsm() */


/**
 *  @brief DTRMM wrapper
 */
void xtrmm
(
  const char *side, const char *uplo,
  const char *transA, const char *diag,
  int m, int n,
  double alpha,
  double *A, int lda,
  double *B, int ldb
)
{
#ifdef USE_BLAS
  dtrmm_
  (
    side, uplo,
	  transA, diag,
	  &m, &n,
	  &alpha,
	  A, &lda,
	  B, &ldb
  );
#else
  printf( "xtrmm must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xtrmm() */


/**
 *  @brief DTRMM wrapper
 */
void xtrmm
(
  const char *side, const char *uplo,
  const char *transA, const char *diag,
  int m, int n,
  float alpha,
  float *A, int lda,
  float *B, int ldb
)
{
#ifdef USE_BLAS
  strmm_
  (
    side, uplo,
	  transA, diag,
	  &m, &n,
	  &alpha,
	  A, &lda,
	  B, &ldb
  );
#else
  printf( "xtrmm must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xtrmm() */


/**
 *  LAPACK routine wrappers: POTR(F,S), GETR(F,S), GECON, GEQRF,
 *  ORGQR, ORMQR, GEQP3, GELS
 */


/**
 *  @brief DLASWP wrapper
 */
void xlaswp( int n, double *A, int lda,
    int k1, int k2, int *ipiv, int incx )
{
#ifdef USE_BLAS
  dlaswp_( &n, A, &lda, &k1, &k2, ipiv, &incx );
#else
  printf( "xlaswp must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xlaswp() */


/**
 *  @brief SLASWP wrapper
 */
void xlaswp( int n, float *A, int lda,
    int k1, int k2, int *ipiv, int incx )
{
#ifdef USE_BLAS
  slaswp_( &n, A, &lda, &k1, &k2, ipiv, &incx );
#else
  printf( "xlaswp must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xlaswp() */


/**
 *  @brief DPOTRF wrapper
 */
void xpotrf( const char *uplo, int n, double *A, int lda )
{
#ifdef USE_BLAS
  int info;
  dpotrf_( uplo, &n, A, &lda, &info );
  if ( info ) printf( "xpotrf error code %d\n", info );
#else
  printf( "xpotrf must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xpotrf() */


/**
 *  @brief SPOTRF wrapper
 */
void xpotrf( const char *uplo, int n, float *A, int lda )
{
#ifdef USE_BLAS
  int info;
  spotrf_( uplo, &n, A, &lda, &info );
  if ( info ) printf( "xpotrf error code %d\n", info );
#else
  printf( "xpotrf must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xpotrf() */


/**
 *  @brief DPOTRS wrapper
 */
void xpotrs( const char *uplo,
  int n, int nrhs, double *A, int lda, double *B, int ldb )
{
#ifdef USE_BLAS
  int info;
  dpotrs_( uplo, &n, &nrhs, A, &lda, B, &ldb, &info );
  if ( info ) printf( "xpotrs error code %d\n", info );
#else
  printf( "xpotrs must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xpotrs() */


/**
 *  @brief SPOTRS wrapper
 */
void xpotrs( const char *uplo,
  int n, int nrhs, float *A, int lda, float *B, int ldb )
{
#ifdef USE_BLAS
  int info;
  spotrs_( uplo, &n, &nrhs, A, &lda, B, &ldb, &info );
  if ( info ) printf( "xpotrs error code %d\n", info );
#else
  printf( "xpotrs must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xpotrs() */


/**
 *  @brief DGETRF wrapper
 */
void xgetrf( int m, int n, double *A, int lda, int *ipiv )
{
#ifdef USE_BLAS
  int info;
  dgetrf_( &m, &n, A, &lda, ipiv, &info );
#else
  printf( "xgetrf must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xgetrf() */


/**
 *  @brief SGETRF wrapper
 */
void xgetrf( int m, int n, float *A, int lda, int *ipiv )
{
#ifdef USE_BLAS
  int info;
  sgetrf_( &m, &n, A, &lda, ipiv, &info );
#else
  printf( "xgetrf must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xgetrf() */


/**
 *  @brief DGETRS wrapper
 */
void xgetrs
(
  const char *trans,
  int m, int nrhs,
  double *A, int lda, int *ipiv,
  double *B, int ldb
)
{
#ifdef USE_BLAS
  int info;
  dgetrs_
  (
    trans,
    &m, &nrhs,
    A, &lda, ipiv,
    B, &ldb, &info
  );
#else
  printf( "xgetrs must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xgetrs() */


/**
 *  @brief SGETRS wrapper
 */
void xgetrs
(
  const char *trans,
  int m, int nrhs,
  float *A, int lda, int *ipiv,
  float *B, int ldb
)
{
#ifdef USE_BLAS
  int info;
  sgetrs_
  (
    trans,
    &m, &nrhs,
    A, &lda, ipiv,
    B, &ldb, &info
  );
#else
  printf( "xgetrs must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xgetrs() */


/**
 *  @brief DGECON wrapper
 */
void xgecon
(
  const char *norm,
  int n,
  double *A, int lda,
  double anorm,
  double *rcond,
  double *work, int *iwork
)
{
#ifdef USE_BLAS
  int info;
  dgecon_
  (
    norm,
    &n,
    A, &lda,
    &anorm,
    rcond,
    work, iwork, &info
  );
#else
  printf( "xgecon must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xgecon() */



/**
 *  @brief SGECON wrapper
 */
void xgecon
(
  const char *norm,
  int n,
  float *A, int lda,
  float anorm,
  float *rcond,
  float *work, int *iwork
)
{
#ifdef USE_BLAS
  int info;
  sgecon_
  (
    norm,
    &n,
    A, &lda,
    &anorm,
    rcond,
    work, iwork, &info
  );
#else
  printf( "xgecon must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xgecon() */


/**
 *  @brief DGEQRF wrapper
 */
void xgeqrf
(
  int m, int n,
  double *A, int lda,
  double *tau,
  double *work, int lwork
)
{
#ifdef USE_BLAS
  int info;
  dgeqrf_
  (
    &m, &n,
    A, &lda,
    tau,
    work, &lwork, &info
  );
  if ( info )
  {
    printf( "xgeqrf has illegal values at parameter %d\n", info );
  }
#else
  printf( "xgeqrf must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xgeqrf() */



/**
 *  @brief SGEQRF wrapper
 */
void xgeqrf
(
  int m, int n,
  float *A, int lda,
  float *tau,
  float *work, int lwork
)
{
#ifdef USE_BLAS
  int info;
  sgeqrf_
  (
    &m, &n,
    A, &lda,
    tau,
    work, &lwork, &info
  );
  if ( info )
  {
    printf( "xgeqrf has illegal values at parameter %d\n", info );
  }
#else
  printf( "xgeqrf must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xgeqrf() */


/**
 *  @brief SORGQR wrapper
 */
void xorgqr
(
  int m, int n, int k,
  double *A, int lda,
  double *tau,
  double *work, int lwork
)
{
#ifdef USE_BLAS
  int info;
  dorgqr_
  (
    &m, &n, &k,
    A, &lda,
    tau,
    work, &lwork, &info
  );
#else
  printf( "xorgqr must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xorgqr() */


/**
 *  @brief SORGQR wrapper
 */
void xorgqr
(
  int m, int n, int k,
  float *A, int lda,
  float *tau,
  float *work, int lwork
)
{
#ifdef USE_BLAS
  int info;
  sorgqr_
  (
    &m, &n, &k,
    A, &lda,
    tau,
    work, &lwork, &info
  );
#else
  printf( "xorgqr must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xorgqr() */


/**
 *  @brief DORMQR wrapper
 */
void xormqr
(
  const char *side, const char *trans,
  int m, int n, int k,
  double *A, int lda,
  double *tau,
  double *C, int ldc,
  double *work, int lwork
)
{
#ifdef USE_BLAS
  int info;
  dormqr_
  (
    side, trans,
    &m, &n, &k,
    A, &lda,
    tau,
    C, &ldc,
    work, &lwork, &info
  );
#else
  printf( "xormqr must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xormqr() */


/**
 *  @brief SORMQR wrapper
 */
void xormqr
(
  const char *side, const char *trans,
  int m, int n, int k,
  float *A, int lda,
  float *tau,
  float *C, int ldc,
  float *work, int lwork
)
{
#ifdef USE_BLAS
  int info;
  sormqr_
  (
    side, trans,
    &m, &n, &k,
    A, &lda,
    tau,
    C, &ldc,
    work, &lwork, &info
  );
#else
  printf( "xormqr must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xormqr() */


/**
 *  @brief DGEQP3 wrapper
 */
void xgeqp3
(
  int m, int n,
  double *A, int lda, int *jpvt,
  double *tau,
  double *work, int lwork
)
{
#ifdef USE_BLAS
  int info;
  dgeqp3_
  (
    &m, &n,
    A, &lda, jpvt,
    tau,
    work, &lwork, &info
  );
  if ( info )
  {
    printf( "xgeqp3 has illegal values at parameter %d\n", info );
  }
#else
  printf( "xgeqp3 must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end geqp3() */


/**
 *  @brief SGEQP3 wrapper
 */
void xgeqp3
(
  int m, int n,
  float *A, int lda, int *jpvt,
  float *tau,
  float *work, int lwork
)
{
#ifdef USE_BLAS
  int info;
  sgeqp3_
  (
    &m, &n,
    A, &lda, jpvt,
    tau,
    work, &lwork, &info
  );
  if ( info )
  {
    printf( "xgeqp3 has illegal values at parameter %d\n", info );
  }
#else
  printf( "xgeqp3 must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end geqp3() */


/**
 *  @brief DGEQP4 wrapper
 */
void xgeqp4
(
  int m, int n,
  double *A, int lda, int *jpvt,
  double *tau,
  double *work, int lwork
)
{
#ifdef USE_BLAS
  int info;
  dgeqp4
  (
    &m, &n,
    A, &lda, jpvt,
    tau,
    work, &lwork, &info
  );
  if ( info )
  {
    printf( "xgeqp4 has illegal values at parameter %d\n", info );
  }
#else
  printf( "xgeqp4 must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end geqp4() */


/**
 *  @brief SGEQP4 wrapper
 */
void xgeqp4
(
  int m, int n,
  float *A, int lda, int *jpvt,
  float *tau,
  float *work, int lwork
)
{
#ifdef USE_BLAS
  int info;
  sgeqp4
  (
    &m, &n,
    A, &lda, jpvt,
    tau,
    work, &lwork, &info
  );
  if ( info )
  {
    printf( "xgeqp4 has illegal values at parameter %d\n", info );
  }
#else
  printf( "xgeqp4 must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end geqp4() */


/**
 *  @brief DGELS wrapper
 */
void xgels
(
  const char *trans,
  int m, int n, int nrhs,
  double *A, int lda,
  double *B, int ldb,
  double *work, int lwork
)
{
#ifdef USE_BLAS
  int info;
  dgels_
  (
    trans,
    &m, &n, &nrhs,
    A, &lda,
    B, &ldb,
    work, &lwork, &info
  );
  if ( info )
  {
    printf( "xgels has illegal values at parameter %d\n", info );
  }
#else
  printf( "xgels must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end gels() */


/**
 *  @brief SGELS wrapper
 */
void xgels
(
  const char *trans,
  int m, int n, int nrhs,
  float *A, int lda,
  float *B, int ldb,
  float *work, int lwork
)
{
#ifdef USE_BLAS
  int info;
  sgels_
  (
    trans,
    &m, &n, &nrhs,
    A, &lda,
    B, &ldb,
    work, &lwork, &info
  );
  if ( info )
  {
    printf( "xgels has illegal values at parameter %d\n", info );
  }
#else
  printf( "xgels must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end gels() */



/**
 *  @brief DGESDD wrapper
 */
void xgesdd
(
  const char *jobz,
  int m, int n,
  double *A, int lda,
  double *S,
  double *U, int ldu,
  double *VT, int ldvt,
  double *work, int lwork, int *iwork
)
{
#ifdef USE_BLAS
  int info;
  dgesdd_
  (
    jobz,
    &m, &n,
    A, &lda,
    S,
    U, &ldu,
    VT, &ldvt,
    work, &lwork, iwork, &info
  );
#else
  printf( "xgesdd must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xgesdd() */

/**
 *  @brief SGESDD wrapper
 */
void xgesdd
(
  const char *jobz,
  int m, int n,
  float *A, int lda,
  float *S,
  float *U, int ldu,
  float *VT, int ldvt,
  float *work, int lwork, int *iwork
)
{
#ifdef USE_BLAS
  int info;
  sgesdd_
  (
    jobz,
    &m, &n,
    A, &lda,
    S,
    U, &ldu,
    VT, &ldvt,
    work, &lwork, iwork, &info
  );
#else
  printf( "xgesdd must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xgesdd() */



void xstev
(
  const char *jobz,
  int n,
  double *D,
  double *E,
  double *Z, int ldz,
  double *work
)
{
#ifdef USE_BLAS
  int info;
  dstev_( jobz, &n, D, E, Z, &ldz, work, &info );
#else
  printf( "xstev must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xstev() */


void xstev
(
  const char *jobz,
  int n,
  float *D,
  float *E,
  float *Z, int ldz,
  float *work
)
{
#ifdef USE_BLAS
  int info;
  sstev_( jobz, &n, D, E, Z, &ldz, work, &info );
#else
  printf( "xstev must enables USE_BLAS.\n" );
  exit( 1 );
#endif
}; /** end xstev() */


}; /** end namespace hmlp */


Generated by: GCOVR (Version 3.2)

Line	Exec	Source
1		/**
2		* HMLP (High-Performance Machine Learning Primitives)
3		*
4		* Copyright (C) 2014-2017, The University of Texas at Austin
5		*
6		* This program is free software: you can redistribute it and/or modify
7		* it under the terms of the GNU General Public License as published by
8		* the Free Software Foundation, either version 3 of the License, or
9		* (at your option) any later version.
10		*
11		* This program is distributed in the hope that it will be useful,
12		* but WITHOUT ANY WARRANTY; without even the implied warranty of
13		* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14		* GNU General Public License for more details.
15		*
16		* You should have received a copy of the GNU General Public License
17		* along with this program. If not, see the LICENSE file.
18		*
19		**/
20
21		#include <blas_lapack.hpp>
22
23		// #define DEBUG_XGEMM 1
24
25		#ifndef USE_BLAS
26		#warning BLAS/LAPACK routines are not compiled (-HMLP_USE_BLAS=false)
27		#endif
28
29		extern "C"
30		{
31		#include <external/blas_lapack_prototypes.h>
32		}; /** end extern "C" */
33
34
35
36
37
38
39		namespace hmlp
40		{
41
42		/**
43		* BLAS level-1 wrappers: DOT, NRM2
44		*/
45
46
47		/**
48		* @brief DDOT wrapper
49		*/
50		double xdot( int n, const double dx, int incx, const double dy, int incy )
51		{
52		double ret_val;
53		#ifdef USE_BLAS
54		ret_val = ddot_( &n, dx, &incx, dy, &incy );
55		#else
56		printf( "xdot must enables USE_BLAS.\n" );
57		exit( 1 );
58		#endif
59		return ret_val;
60		}; /** end xdot() */
61
62
63
64		/**
65		* @brief SDOT wrapper
66		*/
67		float xdot( int n, const float dx, int incx, const float dy, int incy )
68		{
69		float ret_val;
70		#ifdef USE_BLAS
71		ret_val = sdot_( &n, dx, &incx, dy, &incy );
72		#else
73		printf( "xdot must enables USE_BLAS.\n" );
74		exit( 1 );
75		#endif
76		return ret_val;
77		}; /** end xdot() */
78
79
80		/**
81		* @brief DNRM2 wrapper
82		*/
83		double xnrm2( int n, double *x, int incx )
84		{
85		double ret_val;
86		#ifdef USE_BLAS
87		ret_val = dnrm2_( &n, x, &incx );
88		#else
89		printf( "xnrm2 must enables USE_BLAS.\n" );
90		exit( 1 );
91		#endif
92		return ret_val;
93		}; /** end xnrm2() */
94
95
96		/**
97		* @brief SNRM2 wrapper
98		*/
99		float xnrm2( int n, float *x, int incx )
100		{
101		float ret_val;
102		#ifdef USE_BLAS
103		ret_val = snrm2_( &n, x, &incx );
104		#else
105		printf( "xnrm2 must enables USE_BLAS.\n" );
106		exit( 1 );
107		#endif
108		return ret_val;
109		}; /** end xnrm2() */
110
111
112
113
114
115
116
117
118
119
120
121		/**
122		* BLAS level-3 wrappers: GEMM, TRSM
123		*/
124
125
126		/**
127		* @brief DGEMM wrapper
128		*/
129		void xgemm
130		(
131		const char transA, const char transB,
132		int m, int n, int k,
133		double alpha, const double *A, int lda,
134		const double *B, int ldb,
135		double beta, double *C, int ldc
136		)
137		{
138		double beg, xgemm_time = 0.0;
139		double gflops = (double)( m ) * n * ( 2 * k ) / 1E+9;
140		beg = omp_get_wtime();
141		#ifdef USE_BLAS
142		dgemm_
143		(
144		transA, transB,
145		&m, &n, &k,
146		&alpha, A, &lda,
147		B, &ldb,
148		&beta, C, &ldc
149		);
150		#else
151		//printf( "xgemm: configure HMLP_USE_BLAS=true to enable BLAS support.\n" );
152		for ( int p = 0; p < k; p ++ )
153		{
154		for ( int j = 0; j < n; j ++ )
155		{
156		for ( int i = 0; i < m; i ++ )
157		{
158		double a, b;
159		if ( transA == 'T' ) a = A[ i lda + p ];
160		else a = A[ p * lda + i ];
161		if ( transB == 'T' ) b = B[ p ldb + j ];
162		else b = B[ j * ldb + p ];
163		if ( p == 0 )
164		{
165		C[ j * ldc + i ] = beta * C[ j * ldc + i ] + alpha * a * b;
166		}
167		else
168		{
169		C[ j * ldc + i ] += alpha * a * b;
170		}
171		}
172		}
173		}
174		#endif
175		xgemm_time = omp_get_wtime() - beg;
176		#ifdef DEBUG_XGEMM
177		printf( "dgemm %s%s m %d n %d k %d, %5.2lf GFLOPS %5.2lf s\n",
178		transA, transB, m, n, k, gflops / xgemm_time, xgemm_time );
179		#endif
180
181		#ifdef DEBUG_XGEMM
182		printf( "hmlp::xgemm debug\n" );
183		for ( int i = 0; i < m; i ++ )
184		{
185		for ( int j = 0; j < n; j ++ )
186		{
187		printf( "%E ", C[ j * ldc + i ] / alpha );
188		}
189		printf( "\n" );
190		}
191		#endif
192		};
193
194		/**
195		* @brief SGEMM wrapper
196		*/
197		void xgemm
198		(
199		const char transA, const char transB,
200		int m, int n, int k,
201		float alpha, const float *A, int lda,
202		const float *B, int ldb,
203		float beta, float *C, int ldc
204		)
205		{
206		double beg, xgemm_time = 0.0;
207		double gflops = (double)( m ) * n * ( 2 * k ) / 1E+9;
208		beg = omp_get_wtime();
209		#ifdef USE_BLAS
210		sgemm_
211		(
212		transA, transB,
213		&m, &n, &k,
214		&alpha, A, &lda,
215		B, &ldb,
216		&beta, C, &ldc
217		);
218		#else
219		//printf( "xgemm: configure HMLP_USE_BLAS=true to enable BLAS support.\n" );
220		for ( int p = 0; p < k; p ++ )
221		{
222		for ( int j = 0; j < n; j ++ )
223		{
224		for ( int i = 0; i < m; i ++ )
225		{
226		double a, b;
227		if ( transA == 'T' ) a = A[ i lda + p ];
228		else a = A[ p * lda + i ];
229		if ( transB == 'T' ) b = B[ p ldb + j ];
230		else b = B[ j * ldb + p ];
231		if ( p == 0 )
232		{
233		C[ j * ldc + i ] = beta * C[ j * ldc + i ] + alpha * a * b;
234		}
235		else
236		{
237		C[ j * ldc + i ] += alpha * a * b;
238		}
239		}
240		}
241		}
242		#endif
243		xgemm_time = omp_get_wtime() - beg;
244		#ifdef DEBUG_XGEMM
245		printf( "sgemm %s%s m %d n %d k %d, %5.2lf GFLOPS %5.2lf s\n",
246		transA, transB, m, n, k, gflops / xgemm_time, xgemm_time );
247		#endif
248		#ifdef DEBUG_XGEMM
249		printf( "hmlp::xgemm debug\n" );
250		for ( int i = 0; i < m; i ++ )
251		{
252		for ( int j = 0; j < n; j ++ )
253		{
254		printf( "%E ", C[ j * ldc + i ] );
255		}
256		printf( "\n" );
257		}
258		#endif
259		};
260
261
262		void xsyrk
263		(
264		const char uplo, const char trans,
265		int n, int k,
266		double alpha, double *A, int lda,
267		double beta, double *C, int ldc
268		)
269		{
270		#ifdef USE_BLAS
271		dsyrk_
272		(
273		uplo, trans,
274		&n, &k,
275		&alpha, A, &lda,
276		&beta, C, &ldc
277		);
278		#else
279		printf( "xsyrk must enables USE_BLAS.\n" );
280		exit( 1 );
281		#endif
282		}; /** end xsyrk() */
283
284
285		void xsyrk
286		(
287		const char uplo, const char trans,
288		int n, int k,
289		float alpha, float *A, int lda,
290		float beta, float *C, int ldc
291		)
292		{
293		#ifdef USE_BLAS
294		ssyrk_
295		(
296		uplo, trans,
297		&n, &k,
298		&alpha, A, &lda,
299		&beta, C, &ldc
300		);
301		#else
302		printf( "xsyrk must enables USE_BLAS.\n" );
303		exit( 1 );
304		#endif
305		}; /** end xsyrk() */
306
307
308
309
310
311		/**
312		* @brief DTRSM wrapper
313		*/
314		void xtrsm
315		(
316		const char side, const char uplo,
317		const char transA, const char diag,
318		int m, int n,
319		double alpha,
320		double *A, int lda,
321		double *B, int ldb
322		)
323		{
324		double beg, xtrsm_time = 0.0;
325		double gflops = (double)( m ) * ( m - 1 ) * n / 1E+9;
326		beg = omp_get_wtime();
327
328		#ifdef USE_BLAS
329		dtrsm_
330		(
331		side, uplo,
332		transA, diag,
333		&m, &n,
334		&alpha,
335		A, &lda,
336		B, &ldb
337		);
338		#else
339		printf( "xtrsm must enables USE_BLAS.\n" );
340		exit( 1 );
341		#endif
342
343		xtrsm_time = omp_get_wtime() - beg;
344		#ifdef DEBUG_XTRSM
345		printf( "dtrsm m %d n %d, %5.2lf GFLOPS, %5.2lf s\n",
346		m, n, gflops / xtrsm_time, xtrsm_time );
347		#endif
348		}; /** end xtrsm() */
349
350
351		/**
352		* @brief STRSM wrapper
353		*/
354		void xtrsm
355		(
356		const char side, const char uplo,
357		const char transA, const char diag,
358		int m, int n,
359		float alpha,
360		float *A, int lda,
361		float *B, int ldb
362		)
363		{
364		#ifdef USE_BLAS
365		strsm_
366		(
367		side, uplo,
368		transA, diag,
369		&m, &n,
370		&alpha,
371		A, &lda,
372		B, &ldb
373		);
374		#else
375		printf( "xtrsm must enables USE_BLAS.\n" );
376		exit( 1 );
377		#endif
378		}; /** end xtrsm() */
379
380
381		/**
382		* @brief DTRMM wrapper
383		*/
384		void xtrmm
385		(
386		const char side, const char uplo,
387		const char transA, const char diag,
388		int m, int n,
389		double alpha,
390		double *A, int lda,
391		double *B, int ldb
392		)
393		{
394		#ifdef USE_BLAS
395		dtrmm_
396		(
397		side, uplo,
398		transA, diag,
399		&m, &n,
400		&alpha,
401		A, &lda,
402		B, &ldb
403		);
404		#else
405		printf( "xtrmm must enables USE_BLAS.\n" );
406		exit( 1 );
407		#endif
408		}; /** end xtrmm() */
409
410
411		/**
412		* @brief DTRMM wrapper
413		*/
414		void xtrmm
415		(
416		const char side, const char uplo,
417		const char transA, const char diag,
418		int m, int n,
419		float alpha,
420		float *A, int lda,
421		float *B, int ldb
422		)
423		{
424		#ifdef USE_BLAS
425		strmm_
426		(
427		side, uplo,
428		transA, diag,
429		&m, &n,
430		&alpha,
431		A, &lda,
432		B, &ldb
433		);
434		#else
435		printf( "xtrmm must enables USE_BLAS.\n" );
436		exit( 1 );
437		#endif
438		}; /** end xtrmm() */
439
440
441		/**
442		* LAPACK routine wrappers: POTR(F,S), GETR(F,S), GECON, GEQRF,
443		* ORGQR, ORMQR, GEQP3, GELS
444		*/
445
446
447		/**
448		* @brief DLASWP wrapper
449		*/
450		void xlaswp( int n, double *A, int lda,
451		int k1, int k2, int *ipiv, int incx )
452		{
453		#ifdef USE_BLAS
454		dlaswp_( &n, A, &lda, &k1, &k2, ipiv, &incx );
455		#else
456		printf( "xlaswp must enables USE_BLAS.\n" );
457		exit( 1 );
458		#endif
459		}; /** end xlaswp() */
460
461
462		/**
463		* @brief SLASWP wrapper
464		*/
465		void xlaswp( int n, float *A, int lda,
466		int k1, int k2, int *ipiv, int incx )
467		{
468		#ifdef USE_BLAS
469		slaswp_( &n, A, &lda, &k1, &k2, ipiv, &incx );
470		#else
471		printf( "xlaswp must enables USE_BLAS.\n" );
472		exit( 1 );
473		#endif
474		}; /** end xlaswp() */
475
476
477		/**
478		* @brief DPOTRF wrapper
479		*/
480		void xpotrf( const char uplo, int n, double A, int lda )
481		{
482		#ifdef USE_BLAS
483		int info;
484		dpotrf_( uplo, &n, A, &lda, &info );
485		if ( info ) printf( "xpotrf error code %d\n", info );
486		#else
487		printf( "xpotrf must enables USE_BLAS.\n" );
488		exit( 1 );
489		#endif
490		}; /** end xpotrf() */
491
492
493		/**
494		* @brief SPOTRF wrapper
495		*/
496		void xpotrf( const char uplo, int n, float A, int lda )
497		{
498		#ifdef USE_BLAS
499		int info;
500		spotrf_( uplo, &n, A, &lda, &info );
501		if ( info ) printf( "xpotrf error code %d\n", info );
502		#else
503		printf( "xpotrf must enables USE_BLAS.\n" );
504		exit( 1 );
505		#endif
506		}; /** end xpotrf() */
507
508
509		/**
510		* @brief DPOTRS wrapper
511		*/
512		void xpotrs( const char *uplo,
513		int n, int nrhs, double A, int lda, double B, int ldb )
514		{
515		#ifdef USE_BLAS
516		int info;
517		dpotrs_( uplo, &n, &nrhs, A, &lda, B, &ldb, &info );
518		if ( info ) printf( "xpotrs error code %d\n", info );
519		#else
520		printf( "xpotrs must enables USE_BLAS.\n" );
521		exit( 1 );
522		#endif
523		}; /** end xpotrs() */
524
525
526		/**
527		* @brief SPOTRS wrapper
528		*/
529		void xpotrs( const char *uplo,
530		int n, int nrhs, float A, int lda, float B, int ldb )
531		{
532		#ifdef USE_BLAS
533		int info;
534		spotrs_( uplo, &n, &nrhs, A, &lda, B, &ldb, &info );
535		if ( info ) printf( "xpotrs error code %d\n", info );
536		#else
537		printf( "xpotrs must enables USE_BLAS.\n" );
538		exit( 1 );
539		#endif
540		}; /** end xpotrs() */
541
542
543		/**
544		* @brief DGETRF wrapper
545		*/
546		void xgetrf( int m, int n, double A, int lda, int ipiv )
547		{
548		#ifdef USE_BLAS
549		int info;
550		dgetrf_( &m, &n, A, &lda, ipiv, &info );
551		#else
552		printf( "xgetrf must enables USE_BLAS.\n" );
553		exit( 1 );
554		#endif
555		}; /** end xgetrf() */
556
557
558		/**
559		* @brief SGETRF wrapper
560		*/
561		void xgetrf( int m, int n, float A, int lda, int ipiv )
562		{
563		#ifdef USE_BLAS
564		int info;
565		sgetrf_( &m, &n, A, &lda, ipiv, &info );
566		#else
567		printf( "xgetrf must enables USE_BLAS.\n" );
568		exit( 1 );
569		#endif
570		}; /** end xgetrf() */
571
572
573		/**
574		* @brief DGETRS wrapper
575		*/
576		void xgetrs
577		(
578		const char *trans,
579		int m, int nrhs,
580		double A, int lda, int ipiv,
581		double *B, int ldb
582		)
583		{
584		#ifdef USE_BLAS
585		int info;
586		dgetrs_
587		(
588		trans,
589		&m, &nrhs,
590		A, &lda, ipiv,
591		B, &ldb, &info
592		);
593		#else
594		printf( "xgetrs must enables USE_BLAS.\n" );
595		exit( 1 );
596		#endif
597		}; /** end xgetrs() */
598
599
600		/**
601		* @brief SGETRS wrapper
602		*/
603		void xgetrs
604		(
605		const char *trans,
606		int m, int nrhs,
607		float A, int lda, int ipiv,
608		float *B, int ldb
609		)
610		{
611		#ifdef USE_BLAS
612		int info;
613		sgetrs_
614		(
615		trans,
616		&m, &nrhs,
617		A, &lda, ipiv,
618		B, &ldb, &info
619		);
620		#else
621		printf( "xgetrs must enables USE_BLAS.\n" );
622		exit( 1 );
623		#endif
624		}; /** end xgetrs() */
625
626
627		/**
628		* @brief DGECON wrapper
629		*/
630		void xgecon
631		(
632		const char *norm,
633		int n,
634		double *A, int lda,
635		double anorm,
636		double *rcond,
637		double work, int iwork
638		)
639		{
640		#ifdef USE_BLAS
641		int info;
642		dgecon_
643		(
644		norm,
645		&n,
646		A, &lda,
647		&anorm,
648		rcond,
649		work, iwork, &info
650		);
651		#else
652		printf( "xgecon must enables USE_BLAS.\n" );
653		exit( 1 );
654		#endif
655		}; /** end xgecon() */
656
657
658
659		/**
660		* @brief SGECON wrapper
661		*/
662		void xgecon
663		(
664		const char *norm,
665		int n,
666		float *A, int lda,
667		float anorm,
668		float *rcond,
669		float work, int iwork
670		)
671		{
672		#ifdef USE_BLAS
673		int info;
674		sgecon_
675		(
676		norm,
677		&n,
678		A, &lda,
679		&anorm,
680		rcond,
681		work, iwork, &info
682		);
683		#else
684		printf( "xgecon must enables USE_BLAS.\n" );
685		exit( 1 );
686		#endif
687		}; /** end xgecon() */
688
689
690		/**
691		* @brief DGEQRF wrapper
692		*/
693		void xgeqrf
694		(
695		int m, int n,
696		double *A, int lda,
697		double *tau,
698		double *work, int lwork
699		)
700		{
701		#ifdef USE_BLAS
702		int info;
703		dgeqrf_
704		(
705		&m, &n,
706		A, &lda,
707		tau,
708		work, &lwork, &info
709		);
710		if ( info )
711		{
712		printf( "xgeqrf has illegal values at parameter %d\n", info );
713		}
714		#else
715		printf( "xgeqrf must enables USE_BLAS.\n" );
716		exit( 1 );
717		#endif
718		}; /** end xgeqrf() */
719
720
721
722		/**
723		* @brief SGEQRF wrapper
724		*/
725		void xgeqrf
726		(
727		int m, int n,
728		float *A, int lda,
729		float *tau,
730		float *work, int lwork
731		)
732		{
733		#ifdef USE_BLAS
734		int info;
735		sgeqrf_
736		(
737		&m, &n,
738		A, &lda,
739		tau,
740		work, &lwork, &info
741		);
742		if ( info )
743		{
744		printf( "xgeqrf has illegal values at parameter %d\n", info );
745		}
746		#else
747		printf( "xgeqrf must enables USE_BLAS.\n" );
748		exit( 1 );
749		#endif
750		}; /** end xgeqrf() */
751
752
753		/**
754		* @brief SORGQR wrapper
755		*/
756		void xorgqr
757		(
758		int m, int n, int k,
759		double *A, int lda,
760		double *tau,
761		double *work, int lwork
762		)
763		{
764		#ifdef USE_BLAS
765		int info;
766		dorgqr_
767		(
768		&m, &n, &k,
769		A, &lda,
770		tau,
771		work, &lwork, &info
772		);
773		#else
774		printf( "xorgqr must enables USE_BLAS.\n" );
775		exit( 1 );
776		#endif
777		}; /** end xorgqr() */
778
779
780		/**
781		* @brief SORGQR wrapper
782		*/
783		void xorgqr
784		(
785		int m, int n, int k,
786		float *A, int lda,
787		float *tau,
788		float *work, int lwork
789		)
790		{
791		#ifdef USE_BLAS
792		int info;
793		sorgqr_
794		(
795		&m, &n, &k,
796		A, &lda,
797		tau,
798		work, &lwork, &info
799		);
800		#else
801		printf( "xorgqr must enables USE_BLAS.\n" );
802		exit( 1 );
803		#endif
804		}; /** end xorgqr() */
805
806
807		/**
808		* @brief DORMQR wrapper
809		*/
810		void xormqr
811		(
812		const char side, const char trans,
813		int m, int n, int k,
814		double *A, int lda,
815		double *tau,
816		double *C, int ldc,
817		double *work, int lwork
818		)
819		{
820		#ifdef USE_BLAS
821		int info;
822		dormqr_
823		(
824		side, trans,
825		&m, &n, &k,
826		A, &lda,
827		tau,
828		C, &ldc,
829		work, &lwork, &info
830		);
831		#else
832		printf( "xormqr must enables USE_BLAS.\n" );
833		exit( 1 );
834		#endif
835		}; /** end xormqr() */
836
837
838		/**
839		* @brief SORMQR wrapper
840		*/
841		void xormqr
842		(
843		const char side, const char trans,
844		int m, int n, int k,
845		float *A, int lda,
846		float *tau,
847		float *C, int ldc,
848		float *work, int lwork
849		)
850		{
851		#ifdef USE_BLAS
852		int info;
853		sormqr_
854		(
855		side, trans,
856		&m, &n, &k,
857		A, &lda,
858		tau,
859		C, &ldc,
860		work, &lwork, &info
861		);
862		#else
863		printf( "xormqr must enables USE_BLAS.\n" );
864		exit( 1 );
865		#endif
866		}; /** end xormqr() */
867
868
869		/**
870		* @brief DGEQP3 wrapper
871		*/
872		void xgeqp3
873		(
874		int m, int n,
875		double A, int lda, int jpvt,
876		double *tau,
877		double *work, int lwork
878		)
879		{
880		#ifdef USE_BLAS
881		int info;
882		dgeqp3_
883		(
884		&m, &n,
885		A, &lda, jpvt,
886		tau,
887		work, &lwork, &info
888		);
889		if ( info )
890		{
891		printf( "xgeqp3 has illegal values at parameter %d\n", info );
892		}
893		#else
894		printf( "xgeqp3 must enables USE_BLAS.\n" );
895		exit( 1 );
896		#endif
897		}; /** end geqp3() */
898
899
900		/**
901		* @brief SGEQP3 wrapper
902		*/
903		void xgeqp3
904		(
905		int m, int n,
906		float A, int lda, int jpvt,
907		float *tau,
908		float *work, int lwork
909		)
910		{
911		#ifdef USE_BLAS
912		int info;
913		sgeqp3_
914		(
915		&m, &n,
916		A, &lda, jpvt,
917		tau,
918		work, &lwork, &info
919		);
920		if ( info )
921		{
922		printf( "xgeqp3 has illegal values at parameter %d\n", info );
923		}
924		#else
925		printf( "xgeqp3 must enables USE_BLAS.\n" );
926		exit( 1 );
927		#endif
928		}; /** end geqp3() */
929
930
931		/**
932		* @brief DGEQP4 wrapper
933		*/
934		void xgeqp4
935		(
936		int m, int n,
937		double A, int lda, int jpvt,
938		double *tau,
939		double *work, int lwork
940		)
941		{
942		#ifdef USE_BLAS
943		int info;
944		dgeqp4
945		(
946		&m, &n,
947		A, &lda, jpvt,
948		tau,
949		work, &lwork, &info
950		);
951		if ( info )
952		{
953		printf( "xgeqp4 has illegal values at parameter %d\n", info );
954		}
955		#else
956		printf( "xgeqp4 must enables USE_BLAS.\n" );
957		exit( 1 );
958		#endif
959		}; /** end geqp4() */
960
961
962		/**
963		* @brief SGEQP4 wrapper
964		*/
965		void xgeqp4
966		(
967		int m, int n,
968		float A, int lda, int jpvt,
969		float *tau,
970		float *work, int lwork
971		)
972		{
973		#ifdef USE_BLAS
974		int info;
975		sgeqp4
976		(
977		&m, &n,
978		A, &lda, jpvt,
979		tau,
980		work, &lwork, &info
981		);
982		if ( info )
983		{
984		printf( "xgeqp4 has illegal values at parameter %d\n", info );
985		}
986		#else
987		printf( "xgeqp4 must enables USE_BLAS.\n" );
988		exit( 1 );
989		#endif
990		}; /** end geqp4() */
991
992
993		/**
994		* @brief DGELS wrapper
995		*/
996		void xgels
997		(
998		const char *trans,
999		int m, int n, int nrhs,
1000		double *A, int lda,
1001		double *B, int ldb,
1002		double *work, int lwork
1003		)
1004		{
1005		#ifdef USE_BLAS
1006		int info;
1007		dgels_
1008		(
1009		trans,
1010		&m, &n, &nrhs,
1011		A, &lda,
1012		B, &ldb,
1013		work, &lwork, &info
1014		);
1015		if ( info )
1016		{
1017		printf( "xgels has illegal values at parameter %d\n", info );
1018		}
1019		#else
1020		printf( "xgels must enables USE_BLAS.\n" );
1021		exit( 1 );
1022		#endif
1023		}; /** end gels() */
1024
1025
1026		/**
1027		* @brief SGELS wrapper
1028		*/
1029		void xgels
1030		(
1031		const char *trans,
1032		int m, int n, int nrhs,
1033		float *A, int lda,
1034		float *B, int ldb,
1035		float *work, int lwork
1036		)
1037		{
1038		#ifdef USE_BLAS
1039		int info;
1040		sgels_
1041		(
1042		trans,
1043		&m, &n, &nrhs,
1044		A, &lda,
1045		B, &ldb,
1046		work, &lwork, &info
1047		);
1048		if ( info )
1049		{
1050		printf( "xgels has illegal values at parameter %d\n", info );
1051		}
1052		#else
1053		printf( "xgels must enables USE_BLAS.\n" );
1054		exit( 1 );
1055		#endif
1056		}; /** end gels() */
1057
1058
1059
1060		/**
1061		* @brief DGESDD wrapper
1062		*/
1063		void xgesdd
1064		(
1065		const char *jobz,
1066		int m, int n,
1067		double *A, int lda,
1068		double *S,
1069		double *U, int ldu,
1070		double *VT, int ldvt,
1071		double work, int lwork, int iwork
1072		)
1073		{
1074		#ifdef USE_BLAS
1075		int info;
1076		dgesdd_
1077		(
1078		jobz,
1079		&m, &n,
1080		A, &lda,
1081		S,
1082		U, &ldu,
1083		VT, &ldvt,
1084		work, &lwork, iwork, &info
1085		);
1086		#else
1087		printf( "xgesdd must enables USE_BLAS.\n" );
1088		exit( 1 );
1089		#endif
1090		}; /** end xgesdd() */
1091
1092		/**
1093		* @brief SGESDD wrapper
1094		*/
1095		void xgesdd
1096		(
1097		const char *jobz,
1098		int m, int n,
1099		float *A, int lda,
1100		float *S,
1101		float *U, int ldu,
1102		float *VT, int ldvt,
1103		float work, int lwork, int iwork
1104		)
1105		{
1106		#ifdef USE_BLAS
1107		int info;
1108		sgesdd_
1109		(
1110		jobz,
1111		&m, &n,
1112		A, &lda,
1113		S,
1114		U, &ldu,
1115		VT, &ldvt,
1116		work, &lwork, iwork, &info
1117		);
1118		#else
1119		printf( "xgesdd must enables USE_BLAS.\n" );
1120		exit( 1 );
1121		#endif
1122		}; /** end xgesdd() */
1123
1124
1125
1126		void xstev
1127		(
1128		const char *jobz,
1129		int n,
1130		double *D,
1131		double *E,
1132		double *Z, int ldz,
1133		double *work
1134		)
1135		{
1136		#ifdef USE_BLAS
1137		int info;
1138		dstev_( jobz, &n, D, E, Z, &ldz, work, &info );
1139		#else
1140		printf( "xstev must enables USE_BLAS.\n" );
1141		exit( 1 );
1142		#endif
1143		}; /** end xstev() */
1144
1145
1146		void xstev
1147		(
1148		const char *jobz,
1149		int n,
1150		float *D,
1151		float *E,
1152		float *Z, int ldz,
1153		float *work
1154		)
1155		{
1156		#ifdef USE_BLAS
1157		int info;
1158		sstev_( jobz, &n, D, E, Z, &ldz, work, &info );
1159		#else
1160		printf( "xstev must enables USE_BLAS.\n" );
1161		exit( 1 );
1162		#endif
1163		}; /** end xstev() */
1164
1165
1166		}; /** end namespace hmlp */