BLAS Functions

Table of Contents

For vector/vector, matrix/vector and matrix/matrix multiplications with large dense matrices Blaze relies on the efficiency of BLAS libraries. For this purpose, Blaze implements several convenient C++ wrapper functions for several BLAS functions. The following sections give a complete overview of all available BLAS level 1, 2 and 3 functions.


BLAS Level 1


Dot Product (dotu)

The following wrapper functions provide a generic interface for the BLAS functions for the dot product of two dense vectors (sdot(), ddot(), cdotu_sub(), and zdotu_sub()):

namespace blaze {
float dotu( int n, const float* x, int incX, const float* y, int incY );
double dotu( int n, const double* x, int incX, const double* y, int incY );
complex<float> dotu( int n, const complex<float>* x, int incX,
const complex<float>* y, int incY );
complex<double> dotu( int n, const complex<double>* x, int incX,
const complex<double>* y, int incY );
template< typename VT1, bool TF1, typename VT2, bool TF2 >
ElementType_<VT1> dotu( const DenseVector<VT1,TF1>& x, const DenseVector<VT2,TF2>& y );
} // namespace blaze

Complex Conjugate Dot Product (dotc)

The following wrapper functions provide a generic interface for the BLAS functions for the complex conjugate dot product of two dense vectors (sdot(), ddot(), cdotc_sub(), and zdotc_sub()):

namespace blaze {
float dotc( int n, const float* x, int incX, const float* y, int incY );
double dotc( int n, const double* x, int incX, const double* y, int incY );
complex<float> dotc( int n, const complex<float>* x, int incX,
const complex<float>* y, int incY );
complex<double> dotc( int n, const complex<double>* x, int incX,
const complex<double>* y, int incY );
template< typename VT1, bool TF1, typename VT2, bool TF2 >
ElementType_<VT1> dotc( const DenseVector<VT1,TF1>& x, const DenseVector<VT2,TF2>& y );
} // namespace blaze

Axpy Product (axpy)

The following wrapper functions provide a generic interface for the BLAS functions for the axpy product of two dense vectors (saxpy(), daxpy(), caxpy(), and zaxpy()):

namespace blaze {
void axpy( int n, float alpha, const float* x, int incX, float* y, int incY );
void axpy( int n, double alpha, const double* x, int incX, double* y, int incY );
void axpy( int n, complex<float> alpha, const complex<float>* x,
int incX, complex<float>* y, int incY );
void axpy( int n, complex<double> alpha, const complex<double>* x,
int incX, complex<double>* y, int incY );
template< typename VT1, bool TF1, typename VT2, bool TF2, typename ST >
void axpy( const DenseVector<VT1,TF1>& x, const DenseVector<VT2,TF2>& y, ST alpha );
} // namespace blaze


BLAS Level 2


General Matrix/Vector Multiplication (gemv)

The following wrapper functions provide a generic interface for the BLAS functions for the general matrix/vector multiplication (sgemv(), dgemv(), cgemv(), and zgemv()):

namespace blaze {
void gemv( CBLAS_ORDER layout, CBLAS_TRANSPOSE transA, int m, int n, float alpha,
const float* A, int lda, const float* x, int incX,
float beta, float* y, int incY );
void gemv( CBLAS_ORDER layout, CBLAS_TRANSPOSE transA, int m, int n, double alpha,
const double* A, int lda, const double* x, int incX,
double beta, double* y, int incY );
void gemv( CBLAS_ORDER layout, CBLAS_TRANSPOSE transA, int m, int n, complex<float> alpha,
const complex<float>* A, int lda, const complex<float>* x, int incX,
complex<float> beta, complex<float>* y, int incY );
void gemv( CBLAS_ORDER layout, CBLAS_TRANSPOSE transA, int m, int n, complex<double> alpha,
const complex<double>* A, int lda, const complex<double>* x, int incX,
complex<double> beta, complex<double>* y, int incY );
template< typename VT1, typename MT1, bool SO, typename VT2, typename ST >
void gemv( DenseVector<VT1,false>& y, const DenseMatrix<MT1,SO>& A,
const DenseVector<VT2,false>& x, ST alpha, ST beta );
template< typename VT1, typename VT2, typename MT1, bool SO, typename ST >
void gemv( DenseVector<VT1,true>& y, const DenseVector<VT2,true>& x,
const DenseMatrix<MT1,SO>& A, ST alpha, ST beta );
} // namespace blaze


Triangular Matrix/Vector Multiplication (trmv)

The following wrapper functions provide a generic interface for the BLAS functions for the matrix/vector multiplication with a triangular matrix (strmv(), dtrmv(), ctrmv(), and ztrmv()):

namespace blaze {
void trmv( CBLAS_ORDER order, CBLAS_UPLO uplo, CBLAS_TRANSPOSE transA, CBLAS_DIAG diag,
int n, const float* A, int lda, float* x, int incX );
void trmv( CBLAS_ORDER order, CBLAS_UPLO uplo, CBLAS_TRANSPOSE transA, CBLAS_DIAG diag,
int n, const double* A, int lda, double* x, int incX );
void trmv( CBLAS_ORDER order, CBLAS_UPLO uplo, CBLAS_TRANSPOSE transA, CBLAS_DIAG diag,
int n, const complex<float>* A, int lda, complex<float>* x, int incX );
void trmv( CBLAS_ORDER order, CBLAS_UPLO uplo, CBLAS_TRANSPOSE transA, CBLAS_DIAG diag,
int n, const complex<double>* A, int lda, complex<double>* x, int incX );
template< typename VT, typename MT, bool SO >
void trmv( DenseVector<VT,false>& x, const DenseMatrix<MT,SO>& A, CBLAS_UPLO uplo );
template< typename VT, typename MT, bool SO >
void trmv( DenseVector<VT,true>& x, const DenseMatrix<MT,SO>& A, CBLAS_UPLO uplo );
} // namespace blaze


BLAS Level 3


General Matrix/Matrix Multiplication (gemm)

The following wrapper functions provide a generic interface for the BLAS functions for the general matrix/matrix multiplication (sgemm(), dgemm(), cgemm(), and zgemm()):

namespace blaze {
void gemm( CBLAS_ORDER order, CBLAS_TRANSPOSE transA, CBLAS_TRANSPOSE 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 );
void gemm( CBLAS_ORDER order, CBLAS_TRANSPOSE transA, CBLAS_TRANSPOSE transB,
int m, int n, int k, double alpha, const double* A, int lda,
const double* B, int ldb, double beta, float* C, int ldc );
void gemm( CBLAS_ORDER order, CBLAS_TRANSPOSE transA, CBLAS_TRANSPOSE transB,
int m, int n, int k, complex<float> alpha, const complex<float>* A, int lda,
const complex<float>* B, int ldb, complex<float> beta, float* C, int ldc );
void gemm( CBLAS_ORDER order, CBLAS_TRANSPOSE transA, CBLAS_TRANSPOSE transB,
int m, int n, int k, complex<double> alpha, const complex<double>* A, int lda,
const complex<double>* B, int ldb, complex<double> beta, float* C, int ldc );
template< typename MT1, bool SO1, typename MT2, bool SO2, typename MT3, bool SO3, typename ST >
void gemm( DenseMatrix<MT1,SO1>& C, const DenseMatrix<MT2,SO2>& A,
const DenseMatrix<MT3,SO3>& B, ST alpha, ST beta );
} // namespace blaze


Triangular Matrix/Matrix Multiplication (trmm)

The following wrapper functions provide a generic interface for the BLAS functions for the matrix/matrix multiplication with a triangular matrix (strmm(), dtrmm(), ctrmm(), and ztrmm()):

namespace blaze {
void trmm( CBLAS_ORDER order, CBLAS_SIDE side, CBLAS_UPLO uplo, CBLAS_TRANSPOSE transA,
CBLAS_DIAG diag, int m, int n, float alpha, const float* A,
int lda, float* B, int ldb );
void trmm( CBLAS_ORDER order, CBLAS_SIDE side, CBLAS_UPLO uplo, CBLAS_TRANSPOSE transA,
CBLAS_DIAG diag, int m, int n, double alpha, const double* A,
int lda, double* B, int ldb );
void trmm( CBLAS_ORDER order, CBLAS_SIDE side, CBLAS_UPLO uplo, CBLAS_TRANSPOSE transA,
CBLAS_DIAG diag, int m, int n, complex<float> alpha, const complex<float>* A,
int lda, complex<float>* B, int ldb );
void trmm( CBLAS_ORDER order, CBLAS_SIDE side, CBLAS_UPLO uplo, CBLAS_TRANSPOSE transA,
CBLAS_DIAG diag, int m, int n, complex<double> alpha, const complex<double>* A,
int lda, complex<double>* B, int ldb );
template< typename MT1, bool SO1, typename MT2, bool SO2, typename ST >
void trmm( DenseMatrix<MT1,SO1>& B, const DenseMatrix<MT2,SO2>& A,
CBLAS_SIDE side, CBLAS_UPLO uplo, ST alpha );
} // namespace blaze


Triangular System Solver (trsm)

The following wrapper functions provide a generic interface for the BLAS functions for solving a triangular system of equations (strsm(), dtrsm(), ctrsm(), and ztrsm()):

namespace blaze {
void trsm( CBLAS_ORDER order, CBLAS_SIDE side, CBLAS_UPLO uplo, CBLAS_TRANSPOSE transA,
CBLAS_DIAG diag, int m, int n, float alpha, const float* A,
int lda, float* B, int ldb );
void trsm( CBLAS_ORDER order, CBLAS_SIDE side, CBLAS_UPLO uplo, CBLAS_TRANSPOSE transA,
CBLAS_DIAG diag, int m, int n, double alpha, const double* A,
int lda, double* B, int ldb );
void trsm( CBLAS_ORDER order, CBLAS_SIDE side, CBLAS_UPLO uplo, CBLAS_TRANSPOSE transA,
CBLAS_DIAG diag, int m, int n, complex<float> alpha, const complex<float>* A,
int lda, complex<float>* B, int ldb );
void trsm( CBLAS_ORDER order, CBLAS_SIDE side, CBLAS_UPLO uplo, CBLAS_TRANSPOSE transA,
CBLAS_DIAG diag, int m, int n, complex<double> alpha, const complex<double>* A,
int lda, complex<double>* B, int ldb );
template< typename MT, bool SO, typename VT, bool TF, typename ST >
void trsm( const DenseMatrix<MT,SO>& A, DenseVector<VT,TF>& b,
CBLAS_SIDE side, CBLAS_UPLO uplo, ST alpha );
template< typename MT1, bool SO1, typename MT2, bool SO2, typename ST >
void trsm( const DenseMatrix<MT1,SO1>& A, DenseMatrix<MT2,SO2>& B,
CBLAS_SIDE side, CBLAS_UPLO uplo, ST alpha );
} // namespace blaze


Previous: Matrix Serialization     Next: LAPACK Functions