]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
add support for BLAS::gemm
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 17 Dec 2007 19:55:40 +0000 (19:55 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 17 Dec 2007 19:55:40 +0000 (19:55 +0000)
git-svn-id: https://svn.dealii.org/trunk@15604 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/contrib/blastemplates/lapack_templates.h.in
deal.II/contrib/blastemplates/templates.pl
deal.II/lac/include/lac/lapack_templates.h

index 89c8578f2b008806743debbd0e006359761b8dc8..7800f15b1c6129c34d22ff5d4072ed8960e68fac 100644 (file)
@@ -5,6 +5,13 @@ void dgemv_ (const char* trans, const int* m, const int* n,
             const double* x, const int* incx,
             const double* b, double* y, const int* incy);
 
+// Matrix matrix product
+void dgemm_ (const char* transa, const char* transb,
+       const int* m, const int* n, const int* k,
+       const double* alpha, const double* A, const int* lda,
+       const double* B, const int* ldb,
+       const double* beta, double* C, const int* ldc);
+
 // Compute LU factorization
 void dgetrf_ (const int* m, const int* n, double* A,
              const int* lda, int* ipiv, int* info);
index cc147b9460ce3958f4ff5d4adaab02bb64734113..855076e9a3a13ff9f9618549d8ee3781099de9bc 100644 (file)
@@ -2,7 +2,7 @@
 #    $Id$
 #    Version: $Name$
 #
-#    Copyright (C) 2005, 2006 by the deal authors
+#    Copyright (C) 2005, 2006, 2007 by the deal authors
 #
 #    This file is subject to QPL and may not be  distributed
 #    without copyright and license information. Please refer
@@ -23,7 +23,7 @@ print << 'EOT'
 //    This file was automatically generated from blas.h.in
 //    See blastemplates in the deal.II contrib directory
 //
-//    Copyright (C) 2005 by the deal authors
+//    Copyright (C) 2005, 2006, 2007 by the deal authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -37,6 +37,8 @@ print << 'EOT'
 
 #include <lac/lapack_support.h>
 
+using namespace dealii;
+
 extern "C"
 {
 EOT
index b157ed51ba5a84d989b89dbbd6253cc8c1d7e6ae..cf1463678533fa0917f4b2e0ffe7d39795407720 100644 (file)
@@ -5,7 +5,7 @@
 //    This file was automatically generated from blas.h.in
 //    See blastemplates in the deal.II contrib directory
 //
-//    Copyright (C) 2005, 2006 by the deal authors
+//    Copyright (C) 2005, 2006, 2007 by the deal authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -19,6 +19,8 @@
 
 #include <lac/lapack_support.h>
 
+using namespace dealii;
+
 extern "C"
 {
 // General Matrix
@@ -31,6 +33,17 @@ void sgemv_ (const char* trans, const int* m, const int* n,
             const float* alpha, const float* A, const int* lda,
             const float* x, const int* incx,
             const float* b, float* y, const int* incy);
+// Matrix matrix product
+void dgemm_ (const char* transa, const char* transb,
+       const int* m, const int* n, const int* k,
+       const double* alpha, const double* A, const int* lda,
+       const double* B, const int* ldb,
+       const double* beta, double* C, const int* ldc);
+void sgemm_ (const char* transa, const char* transb,
+       const int* m, const int* n, const int* k,
+       const float* alpha, const float* A, const int* lda,
+       const float* B, const int* ldb,
+       const float* beta, float* C, const int* ldc);
 // Compute LU factorization
 void dgetrf_ (const int* m, const int* n, double* A,
              const int* lda, int* ipiv, int* info);
@@ -120,7 +133,7 @@ gemv (const char* trans, const int* m, const int* n, const double* alpha, const
 inline void
 gemv (const char*, const int*, const int*, const double*, const double*, const int*, const double*, const int*, const double*, double*, const int*)
 {
-  dealii::LAPACKSupport::ExcMissing("dgemv");
+  LAPACKSupport::ExcMissing("dgemv");
 }
 #endif
 
@@ -135,7 +148,37 @@ gemv (const char* trans, const int* m, const int* n, const float* alpha, const f
 inline void
 gemv (const char*, const int*, const int*, const float*, const float*, const int*, const float*, const int*, const float*, float*, const int*)
 {
-  dealii::LAPACKSupport::ExcMissing("sgemv");
+  LAPACKSupport::ExcMissing("sgemv");
+}
+#endif
+
+
+#ifdef HAVE_DGEMM_
+inline void
+gemm (const char* transa, const char* transb, const int* m, const int* n, const int* k, const double* alpha, const double* A, const int* lda, const double* B, const int* ldb, const double* beta, double* C, const int* ldc)
+{
+  dgemm_ (transa,transb,m,n,k,alpha,A,lda,B,ldb,beta,C,ldc);
+}
+#else
+inline void
+gemm (const char*, const char*, const int*, const int*, const int*, const double*, const double*, const int*, const double*, const int*, const double*, double*, const int*)
+{
+  LAPACKSupport::ExcMissing("dgemm");
+}
+#endif
+
+
+#ifdef HAVE_SGEMM_
+inline void
+gemm (const char* transa, const char* transb, const int* m, const int* n, const int* k, const float* alpha, const float* A, const int* lda, const float* B, const int* ldb, const float* beta, float* C, const int* ldc)
+{
+  sgemm_ (transa,transb,m,n,k,alpha,A,lda,B,ldb,beta,C,ldc);
+}
+#else
+inline void
+gemm (const char*, const char*, const int*, const int*, const int*, const float*, const float*, const int*, const float*, const int*, const float*, float*, const int*)
+{
+  LAPACKSupport::ExcMissing("sgemm");
 }
 #endif
 
@@ -150,7 +193,7 @@ getrf (const int* m, const int* n, double* A, const int* lda, int* ipiv, int* in
 inline void
 getrf (const int*, const int*, double*, const int*, int*, int*)
 {
-  dealii::LAPACKSupport::ExcMissing("dgetrf");
+  LAPACKSupport::ExcMissing("dgetrf");
 }
 #endif
 
@@ -165,7 +208,7 @@ getrf (const int* m, const int* n, float* A, const int* lda, int* ipiv, int* inf
 inline void
 getrf (const int*, const int*, float*, const int*, int*, int*)
 {
-  dealii::LAPACKSupport::ExcMissing("sgetrf");
+  LAPACKSupport::ExcMissing("sgetrf");
 }
 #endif
 
@@ -180,7 +223,7 @@ getrs (const char* trans, const int* n, const int* nrhs, const double* A, const
 inline void
 getrs (const char*, const int*, const int*, const double*, const int*, const int*, double*, const int*, int*)
 {
-  dealii::LAPACKSupport::ExcMissing("dgetrs");
+  LAPACKSupport::ExcMissing("dgetrs");
 }
 #endif
 
@@ -195,7 +238,7 @@ getrs (const char* trans, const int* n, const int* nrhs, const float* A, const i
 inline void
 getrs (const char*, const int*, const int*, const float*, const int*, const int*, float*, const int*, int*)
 {
-  dealii::LAPACKSupport::ExcMissing("sgetrs");
+  LAPACKSupport::ExcMissing("sgetrs");
 }
 #endif
 
@@ -210,7 +253,7 @@ geev (const char* jobvl, const char* jobvr, const int* n, double* A, const int*
 inline void
 geev (const char*, const char*, const int*, double*, const int*, double*, double*, double*, const int*, double*, const int*, double*, const int*, int*)
 {
-  dealii::LAPACKSupport::ExcMissing("dgeev");
+  LAPACKSupport::ExcMissing("dgeev");
 }
 #endif
 
@@ -225,7 +268,7 @@ geev (const char* jobvl, const char* jobvr, const int* n, float* A, const int* l
 inline void
 geev (const char*, const char*, const int*, float*, const int*, float*, float*, float*, const int*, float*, const int*, float*, const int*, int*)
 {
-  dealii::LAPACKSupport::ExcMissing("sgeev");
+  LAPACKSupport::ExcMissing("sgeev");
 }
 #endif
 
@@ -240,7 +283,7 @@ geevx (const char* balanc, const char* jobvl, const char* jobvr, const char* sen
 inline void
 geevx (const char*, const char*, const char*, const char*, const int*, double*, const int*, double*, double*, double*, const int*, double*, const int*, int*, int*, double*, double*, double*, double*, double*, const int*, int*, int*)
 {
-  dealii::LAPACKSupport::ExcMissing("dgeevx");
+  LAPACKSupport::ExcMissing("dgeevx");
 }
 #endif
 
@@ -255,7 +298,7 @@ geevx (const char* balanc, const char* jobvl, const char* jobvr, const char* sen
 inline void
 geevx (const char*, const char*, const char*, const char*, const int*, float*, const int*, float*, float*, float*, const int*, float*, const int*, int*, int*, float*, float*, float*, float*, float*, const int*, int*, int*)
 {
-  dealii::LAPACKSupport::ExcMissing("sgeevx");
+  LAPACKSupport::ExcMissing("sgeevx");
 }
 #endif
 
@@ -270,7 +313,7 @@ gesvd (int* jobu, int* jobvt, const int* n, const int* m, double* A, const int*
 inline void
 gesvd (int*, int*, const int*, const int*, double*, const int*, double*, double*, const int*, double*, const int*, double*, const int*, int*)
 {
-  dealii::LAPACKSupport::ExcMissing("dgesvd");
+  LAPACKSupport::ExcMissing("dgesvd");
 }
 #endif
 
@@ -285,7 +328,7 @@ gesvd (int* jobu, int* jobvt, const int* n, const int* m, float* A, const int* l
 inline void
 gesvd (int*, int*, const int*, const int*, float*, const int*, float*, float*, const int*, float*, const int*, float*, const int*, int*)
 {
-  dealii::LAPACKSupport::ExcMissing("sgesvd");
+  LAPACKSupport::ExcMissing("sgesvd");
 }
 #endif
 
@@ -300,7 +343,7 @@ stev (const char* jobz, const int* n, double* d, double* e, double* z, const int
 inline void
 stev (const char*, const int*, double*, double*, double*, const int*, double*, int*)
 {
-  dealii::LAPACKSupport::ExcMissing("dstev");
+  LAPACKSupport::ExcMissing("dstev");
 }
 #endif
 
@@ -315,8 +358,9 @@ stev (const char* jobz, const int* n, float* d, float* e, float* z, const int* l
 inline void
 stev (const char*, const int*, float*, float*, float*, const int*, float*, int*)
 {
-  dealii::LAPACKSupport::ExcMissing("sstev");
+  LAPACKSupport::ExcMissing("sstev");
 }
 #endif
 
+
 #endif

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.