本文整理汇总了C++中teuchos::BLAS类的典型用法代码示例。如果您正苦于以下问题:C++ BLAS类的具体用法?C++ BLAS怎么用?C++ BLAS使用的例子?那么, 这里精选的类代码示例或许可以为您提供帮助。
在下文中一共展示了BLAS类的15个代码示例,这些例子默认根据受欢迎程度排序。您可以为喜欢或者感觉有用的代码点赞,您的评价将有助于系统推荐出更棒的C++代码示例。
示例1: main
int main(int argc, char* argv[])
{
#ifdef HAVE_MPI
MPI_Init(&argc,&argv);
#endif
// Creating an instance of the BLAS class for double-precision kernels looks like:
Teuchos::BLAS<int, double> blas;
// This instance provides the access to all the BLAS kernels listed in Figure \ref{blas_kernels}:
const int n = 10;
double alpha = 2.0;
double x[ n ];
for ( int i=0; i<n; i++ ) { x[i] = i; }
blas.SCAL( n, alpha, x, 1 );
int max_idx = blas.IAMAX( n, x, 1 );
cout<< "The index of the maximum magnitude entry of x[] is the "
<< max_idx <<"-th and x[ " << max_idx-1 << " ] = "<< x[max_idx-1]
<< endl;
#ifdef HAVE_MPI
MPI_Finalize();
#endif
return 0;
}
示例2: main
int main(int argc, char* argv[])
{
std::cout << Teuchos::Teuchos_Version() << std::endl << std::endl;
// Creating an instance of the BLAS class for double-precision kernels looks like:
Teuchos::BLAS<int, double> blas;
// This instance provides the access to all the BLAS kernels like _SCAL.
const int n = 10;
double alpha = 2.0;
double x[ n ];
for ( int i=0; i<n; i++ ) { x[i] = i; }
blas.SCAL( n, alpha, x, 1 );
int max_idx = blas.IAMAX( n, x, 1 );
std::cout<< "The index of the maximum magnitude entry of x[] is the "
<< max_idx <<"-th and x[ " << max_idx-1 << " ] = "<< x[max_idx-1]
<< std::endl;
return 0;
}
示例3: timer
double
do_time_teuchos_fad_dot(unsigned int m, unsigned int ndot, unsigned int nloop)
{
Sacado::Random<double> urand(0.0, 1.0);
Teuchos::BLAS<int,FadType> blas;
std::vector<FadType> X(m), Y(m);
for (unsigned int i=0; i<m; i++) {
X[i] = FadType(ndot, urand.number());
Y[i] = FadType(ndot, urand.number());
for (unsigned int k=0; k<ndot; k++) {
X[i].fastAccessDx(k) = urand.number();
Y[i].fastAccessDx(k) = urand.number();
}
}
Teuchos::Time timer("Teuchos Fad DOT", false);
timer.start(true);
for (unsigned int j=0; j<nloop; j++) {
FadType z = blas.DOT(m, &X[0], 1, &Y[0], 1);
}
timer.stop();
return timer.totalElapsedTime() / nloop;
}
示例4:
void
Stokhos::GramSchmidtBasis<ordinal_type, value_type>::
transformCoeffs(const value_type *in, value_type *out) const
{
Teuchos::BLAS<ordinal_type, value_type> blas;
for (ordinal_type i=0; i<sz; i++)
out[i] = in[i];
blas.TRSM(Teuchos::LEFT_SIDE, Teuchos::LOWER_TRI, Teuchos::TRANS,
Teuchos::UNIT_DIAG, sz, 1, 1.0, gs_mat.values(), sz, out, sz);
}
示例5: computeIntegral
/*
Computes integrals of monomials over a given reference cell.
*/
double computeIntegral(shards::CellTopology & cellTopology, int cubDegree, int xDeg, int yDeg, int zDeg) {
DefaultCubatureFactory<double> cubFactory; // create factory
Teuchos::RCP<Cubature<double> > myCub = cubFactory.create(cellTopology, cubDegree); // create default cubature
double val = 0.0;
int cubDim = myCub->getDimension();
int numCubPoints = myCub->getNumPoints();
FieldContainer<double> point(cubDim);
FieldContainer<double> cubPoints(numCubPoints, cubDim);
FieldContainer<double> cubWeights(numCubPoints);
FieldContainer<double> functValues(numCubPoints);
myCub->getCubature(cubPoints, cubWeights);
for (int i=0; i<numCubPoints; i++) {
for (int j=0; j<cubDim; j++) {
point(j) = cubPoints(i,j);
}
functValues(i) = computeMonomial(point, xDeg, yDeg, zDeg);
}
Teuchos::BLAS<int, double> myblas;
int inc = 1;
val = myblas.DOT(numCubPoints, &functValues[0], inc, &cubWeights[0], inc);
return val;
}
示例6: GEMM
static void
GEMM (const Teuchos::ETransp transA,
const Teuchos::ETransp transB,
const Scalar& alpha,
const View<const Scalar**, LayoutLeft, DeviceType>& A,
const View<const Scalar**, LayoutLeft, DeviceType>& B,
const Scalar& beta,
const View<Scalar**, LayoutLeft, DeviceType>& C)
{
const int n = static_cast<int> (C.dimension_1 ());
const int lda = static_cast<int> (Impl::getStride2DView (A));
Teuchos::BLAS<int,Scalar> blas;
// For some BLAS implementations (e.g., MKL), GEMM when B has
// one column may be signficantly less efficient than GEMV.
if (n == 1 && transB == Teuchos::NO_TRANS) {
blas.GEMV (transA, A.dimension_0 (), A.dimension_1 (),
alpha, A.ptr_on_device (), lda,
B.ptr_on_device (), static_cast<int> (1),
beta, C.ptr_on_device (), static_cast<int> (1));
}
else {
const int m = static_cast<int> (C.dimension_0 ());
const int k = static_cast<int> (transA == Teuchos::NO_TRANS ?
A.dimension_1 () : A.dimension_0 ());
const int ldb = static_cast<int> (Impl::getStride2DView (B));
const int ldc = static_cast<int> (Impl::getStride2DView (C));
blas.GEMM (transA, transB, m, n, k, alpha,
A.ptr_on_device(), lda,
B.ptr_on_device(), ldb,
beta, C.ptr_on_device(), ldc);
}
}
示例7:
void
Stokhos::MonoProjPCEBasis<ordinal_type, value_type>::
transformCoeffs(const value_type *in, value_type *out) const
{
// Transform coefficients to normalized basis
Teuchos::BLAS<ordinal_type, value_type> blas;
blas.GEMV(Teuchos::NO_TRANS, pce_sz, this->p+1,
value_type(1.0), basis_vecs.values(), pce_sz,
in, ordinal_type(1), value_type(0.0), out, ordinal_type(1));
// Transform from normalized to original
for (ordinal_type i=0; i<pce_sz; i++)
out[i] /= pce_norms[i];
}
示例8: _Orthogonalize
void Arnoldi::_Orthogonalize(cSDV& q){
Teuchos::BLAS<int, complex<double> > blas;
complex<double> OrthogFactor;
for( int j = 0; j <= _k; ++j ){
OrthogFactor = blas.DOT(_length, q.values(), 1, _Q[j], 1);
_H[_k][j] = OrthogFactor;
complex<double>* iterQ = _Q[j];
complex<double>* iterq = q.values();
for( int i = 0; i < _length; ++i, ++iterQ, ++iterq ){
*iterq -= OrthogFactor*(*iterQ);
}
}
}
示例9: defined
KOKKOS_INLINE_FUNCTION
int
Gemm<Trans::ConjTranspose,Trans::NoTranspose,
AlgoGemm::ExternalBlas,Variant::One>
::invoke(PolicyType &policy,
MemberType &member,
const ScalarType alpha,
DenseExecViewTypeA &A,
DenseExecViewTypeB &B,
const ScalarType beta,
DenseExecViewTypeC &C) {
// static_assert( Kokkos::Impl::is_same<
// typename DenseMatrixTypeA::space_type,
// typename DenseMatrixTypeB::space_type
// >::value &&
// Kokkos::Impl::is_same<
// typename DenseMatrixTypeB::space_type,
// typename DenseMatrixTypeC::space_type
// >::value,
// "Space type of input matrices does not match" );
if (member.team_rank() == 0) {
#if \
defined( HAVE_SHYLUTACHO_TEUCHOS ) && \
defined( KOKKOS_ACTIVE_EXECUTION_MEMORY_SPACE_HOST )
typedef typename DenseExecViewTypeA::ordinal_type ordinal_type;
typedef typename DenseExecViewTypeA::value_type value_type;
Teuchos::BLAS<ordinal_type,value_type> blas;
const ordinal_type m = C.NumRows();
const ordinal_type n = C.NumCols();
const ordinal_type k = B.NumRows();
if (m > 0 && n > 0 && k > 0)
blas.GEMM(Teuchos::CONJ_TRANS, Teuchos::NO_TRANS,
m, n, k,
alpha,
A.ValuePtr(), A.BaseObject().ColStride(),
B.ValuePtr(), B.BaseObject().ColStride(),
beta,
C.ValuePtr(), C.BaseObject().ColStride());
#else
TACHO_TEST_FOR_ABORT( true, MSG_NOT_HAVE_PACKAGE("Teuchos") );
#endif
}
return 0;
}
示例10: alpha
double
do_time_teuchos_fad_gemm(unsigned int m, unsigned int n, unsigned int k,
unsigned int ndot, unsigned int nloop)
{
Sacado::Random<double> urand(0.0, 1.0);
Teuchos::BLAS<int,FadType> blas;
std::vector<FadType> A(m*k), B(k*n), C(m*n);
for (unsigned int j=0; j<k; j++) {
for (unsigned int i=0; i<m; i++) {
A[i+j*m] = FadType(ndot, urand.number());
for (unsigned int l=0; l<ndot; l++)
A[i+j*m].fastAccessDx(l) = urand.number();
}
}
for (unsigned int j=0; j<n; j++) {
for (unsigned int i=0; i<k; i++) {
B[i+j*k] = FadType(ndot, urand.number());
for (unsigned int l=0; l<ndot; l++)
B[i+j*k].fastAccessDx(l) = urand.number();
}
}
for (unsigned int j=0; j<n; j++) {
for (unsigned int i=0; i<m; i++) {
C[i+j*m] = FadType(ndot, urand.number());
for (unsigned int l=0; l<ndot; l++)
C[i+j*m].fastAccessDx(l) = urand.number();
}
}
FadType alpha(ndot, urand.number());
FadType beta(ndot, urand.number());
for (unsigned int l=0; l<ndot; l++) {
alpha.fastAccessDx(l) = urand.number();
beta.fastAccessDx(l) = urand.number();
}
Teuchos::Time timer("Teuchos Fad GEMM", false);
timer.start(true);
for (unsigned int j=0; j<nloop; j++) {
blas.GEMM(Teuchos::NO_TRANS, Teuchos::NO_TRANS, m, n, k, alpha, &A[0], m,
&B[0], k, beta, &C[0], m);
}
timer.stop();
return timer.totalElapsedTime() / nloop;
}
示例11:
KOKKOS_INLINE_FUNCTION
int
Trsm<Side::Left,Uplo::Upper,Trans::NoTranspose,
AlgoTrsm::ExternalBlas,Variant::One>
::invoke(PolicyType &policy,
const MemberType &member,
const int diagA,
const ScalarType alpha,
DenseExecViewTypeA &A,
DenseExecViewTypeB &B) {
// static_assert( Kokkos::Impl::is_same<
// typename DenseMatrixTypeA::space_type,
// Kokkos::Cuda
// >::value,
// "Cuda space is not available for calling external BLAS" );
// static_assert( Kokkos::Impl::is_same<
// typename DenseMatrixTypeA::space_type,
// typename DenseMatrixTypeB::space_type
// >::value,
// "Space type of input matrices does not match" );
//typedef typename DenseExecViewTypeA::space_type space_type;
typedef typename DenseExecViewTypeA::ordinal_type ordinal_type;
typedef typename DenseExecViewTypeA::value_type value_type;
if (member.team_rank() == 0) {
#ifdef HAVE_SHYLUTACHO_TEUCHOS
Teuchos::BLAS<ordinal_type,value_type> blas;
const ordinal_type m = A.NumRows();
const ordinal_type n = B.NumCols();
blas.TRSM(Teuchos::LEFT_SIDE, Teuchos::UPPER_TRI, Teuchos::NO_TRANS,
(diagA == Diag::Unit ? Teuchos::UNIT_DIAG : Teuchos::NON_UNIT_DIAG),
m, n,
alpha,
A.ValuePtr(), A.BaseObject().ColStride(),
B.ValuePtr(), B.BaseObject().ColStride());
#else
TACHO_TEST_FOR_ABORT( true, MSG_NOT_HAVE_PACKAGE("Teuchos") );
#endif
}
return 0;
}
示例12: if
static void
GEMM (const Teuchos::ETransp transA,
const Teuchos::ETransp transB,
const double& alpha,
const View<const double**, LayoutLeft, DeviceType>& A,
const View<const double**, LayoutLeft, DeviceType>& B,
const double& beta,
const View<double**, LayoutLeft, DeviceType>& C)
{
const int n = static_cast<int> (C.dimension_1 ());
// For some BLAS implementations (e.g., MKL), GEMM when B has
// one column may be signficantly less efficient than GEMV.
if (n == 1 && transB == Teuchos::NO_TRANS) {
char trans = 'N';
if (transA == Teuchos::TRANS) {
trans = 'T';
}
else if (transA == Teuchos::CONJ_TRANS) {
trans = 'C';
}
auto B_0 = Kokkos::subview (B, Kokkos::ALL (), 0);
auto C_0 = Kokkos::subview (C, Kokkos::ALL (), 0);
KokkosBlas::gemv (&trans, alpha, A, B_0, beta, C_0);
}
else {
const int m = static_cast<int> (C.dimension_0 ());
const int k = static_cast<int> (transA == Teuchos::NO_TRANS ? A.dimension_1 () : A.dimension_0 ());
const int lda = static_cast<int> (Impl::getStride2DView (A));
const int ldb = static_cast<int> (Impl::getStride2DView (B));
const int ldc = static_cast<int> (Impl::getStride2DView (C));
Teuchos::BLAS<int,double> blas;
blas.GEMM (transA, transB, m, n, k, alpha,
A.ptr_on_device(), lda,
B.ptr_on_device(), ldb,
beta, C.ptr_on_device(), ldc);
}
}
示例13: computeIntegral
/*
Computes integrals of monomials over a given reference cell.
*/
void computeIntegral(Teuchos::Array<double>& testIntFixDeg, shards::CellTopology & cellTopology, int cubDegree) {
DefaultCubatureFactory<double> cubFactory; // create factory
Teuchos::RCP<Cubature<double> > myCub = cubFactory.create(cellTopology, cubDegree); // create default cubature
int cubDim = myCub->getDimension();
int numCubPoints = myCub->getNumPoints();
int numPolys = (cubDegree+1)*(cubDegree+2)*(cubDegree+3)/6;
FieldContainer<double> point(cubDim);
FieldContainer<double> cubPoints(numCubPoints, cubDim);
FieldContainer<double> cubWeights(numCubPoints);
FieldContainer<double> functValues(numCubPoints, numPolys);
myCub->getCubature(cubPoints, cubWeights);
int polyCt = 0;
for (int xDeg=0; xDeg <= cubDegree; xDeg++) {
for (int yDeg=0; yDeg <= cubDegree-xDeg; yDeg++) {
for (int zDeg=0; zDeg <= cubDegree-xDeg-yDeg; zDeg++) {
for (int i=0; i<numCubPoints; i++) {
for (int j=0; j<cubDim; j++) {
point(j) = cubPoints(i,j);
}
functValues(i,polyCt) = computeMonomial(point, xDeg, yDeg, zDeg);
}
polyCt++;
}
}
}
Teuchos::BLAS<int, double> myblas;
int inc = 1;
double alpha = 1.0;
double beta = 0.0;
myblas.GEMV(Teuchos::NO_TRANS, numPolys, numCubPoints, alpha, &functValues(0,0), numPolys,
&cubWeights(0), inc, beta, &testIntFixDeg[0], inc);
}
示例14: computeIntegral
/*
Computes integrals of monomials over a given reference cell.
*/
void computeIntegral(Teuchos::Array<double>& testIntFixDeg, int cubDegree) {
CubatureGenSparse<double,3> myCub(cubDegree);
int cubDim = myCub.getDimension();
int numCubPoints = myCub.getNumPoints();
int numPolys = (cubDegree+1)*(cubDegree+2)*(cubDegree+3)/6;
FieldContainer<double> point(cubDim);
FieldContainer<double> cubPoints(numCubPoints, cubDim);
FieldContainer<double> cubWeights(numCubPoints);
FieldContainer<double> functValues(numCubPoints, numPolys);
myCub.getCubature(cubPoints, cubWeights);
int polyCt = 0;
for (int xDeg=0; xDeg <= cubDegree; xDeg++) {
for (int yDeg=0; yDeg <= cubDegree-xDeg; yDeg++) {
for (int zDeg=0; zDeg <= cubDegree-xDeg-yDeg; zDeg++) {
for (int i=0; i<numCubPoints; i++) {
for (int j=0; j<cubDim; j++) {
point(j) = cubPoints(i,j);
}
functValues(i,polyCt) = computeMonomial(point, xDeg, yDeg, zDeg);
}
polyCt++;
}
}
}
Teuchos::BLAS<int, double> myblas;
int inc = 1;
double alpha = 1.0;
double beta = 0.0;
myblas.GEMV(Teuchos::NO_TRANS, numPolys, numCubPoints, alpha, &functValues(0,0), numPolys,
&cubWeights(0), inc, beta, &testIntFixDeg[0], inc);
}
示例15: main
int main(int argc, char **argv)
{
const unsigned int n = 5;
Sacado::Fad::Vector<unsigned int, FadType> A(n*n,0),B(n,n), C(n,n);
for (unsigned int i=0; i<n; i++) {
for (unsigned int j=0; j<n; j++)
A[i+j*n] = FadType(Teuchos::ScalarTraits<double>::random());
B[i] = FadType(n, Teuchos::ScalarTraits<double>::random());
for (unsigned int j=0; j<n; j++)
B[i].fastAccessDx(j) = Teuchos::ScalarTraits<double>::random();
C[i] = 0.0;
}
double *a = A.vals();
double *b = B.vals();
double *bdx = B.dx();
std::vector<double> c(n), cdx(n*n);
Teuchos::BLAS<int,double> blas;
blas.GEMV(Teuchos::NO_TRANS, n, n, 1.0, &a[0], n, &b[0], 1, 0.0, &c[0], 1);
blas.GEMM(Teuchos::NO_TRANS, Teuchos::NO_TRANS, n, n, n, 1.0, &a[0], n, &bdx[0], n, 0.0, &cdx[0], n);
// Teuchos::BLAS<int,FadType> blas_fad;
// blas_fad.GEMV(Teuchos::NO_TRANS, n, n, 1.0, &A[0], n, &B[0], 1, 0.0, &C[0], 1);
Teuchos::BLAS<int,FadType> sacado_fad_blas(false);
sacado_fad_blas.GEMV(Teuchos::NO_TRANS, n, n, 1.0, &A[0], n, &B[0], 1, 0.0, &C[0], 1);
// Print the results
int p = 4;
int w = p+7;
std::cout.setf(std::ios::scientific);
std::cout.precision(p);
std::cout << "BLAS GEMV calculation:" << std::endl;
std::cout << "a = " << std::endl;
for (unsigned int i=0; i<n; i++) {
for (unsigned int j=0; j<n; j++)
std::cout << " " << std::setw(w) << a[i+j*n];
std::cout << std::endl;
}
std::cout << "b = " << std::endl;
for (unsigned int i=0; i<n; i++) {
std::cout << " " << std::setw(w) << b[i];
}
std::cout << std::endl;
std::cout << "bdot = " << std::endl;
for (unsigned int i=0; i<n; i++) {
for (unsigned int j=0; j<n; j++)
std::cout << " " << std::setw(w) << bdx[i+j*n];
std::cout << std::endl;
}
std::cout << "c = " << std::endl;
for (unsigned int i=0; i<n; i++) {
std::cout << " " << std::setw(w) << c[i];
}
std::cout << std::endl;
std::cout << "cdot = " << std::endl;
for (unsigned int i=0; i<n; i++) {
for (unsigned int j=0; j<n; j++)
std::cout << " " << std::setw(w) << cdx[i+j*n];
std::cout << std::endl;
}
std::cout << std::endl << std::endl;
std::cout << "FAD BLAS GEMV calculation:" << std::endl;
std::cout << "A.val() (should = a) = " << std::endl;
for (unsigned int i=0; i<n; i++) {
for (unsigned int j=0; j<n; j++)
std::cout << " " << std::setw(w) << A[i+j*n].val();
std::cout << std::endl;
}
std::cout << "B.val() (should = b) = " << std::endl;
for (unsigned int i=0; i<n; i++) {
std::cout << " " << std::setw(w) << B[i].val();
}
std::cout << std::endl;
std::cout << "B.dx() (should = bdot) = " << std::endl;
double *Bdx = B.dx();
for (unsigned int i=0; i<n; i++) {
for (unsigned int j=0; j<n; j++)
std::cout << " " << std::setw(w) << Bdx[i+j*n];
std::cout << std::endl;
}
std::cout << "C.val() (should = c) = " << std::endl;
for (unsigned int i=0; i<n; i++) {
std::cout << " " << std::setw(w) << C[i].val();
}
std::cout << std::endl;
std::cout << "C.dx() (should = cdot) = " << std::endl;
double *Cdx = C.dx();
for (unsigned int i=0; i<n; i++) {
for (unsigned int j=0; j<n; j++)
std::cout << " " << std::setw(w) << Cdx[i+j*n];
std::cout << std::endl;
}
double tol = 1.0e-14;
bool failed = false;
for (unsigned int i=0; i<n; i++) {
//.........这里部分代码省略.........