]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Recognize the LAPACK functions dsygv and ssygv.
authorwillems <willems@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 7 Apr 2010 20:26:59 +0000 (20:26 +0000)
committerwillems <willems@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 7 Apr 2010 20:26:59 +0000 (20:26 +0000)
git-svn-id: https://svn.dealii.org/trunk@20963 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/configure.in
deal.II/lac/include/lac/lapack_templates.h

index cee2b7f5cb76fee8458c282e296dd85384ad7623..7ce8604fceb58f5d5fb42b912f6e4b2fda54d132 100644 (file)
@@ -567,7 +567,7 @@ fi
 AC_CHECK_FUNCS([daxpy_ saxpy_ dgemv_ sgemv_ dgeev_ sgeev_ dgeevx_ sgeevx_])
 AC_CHECK_FUNCS([dgesvd_ sgesvd_ dgemm_ sgemm_ dgetrf_ sgetrf_ dgetri_ sgetri_])
 AC_CHECK_FUNCS([dgeqrf_ sgeqrf_ dormqr_ sormqr_ dorgqr_ sorgqr_ dtrtrs_ strtrs_])
-AC_CHECK_FUNCS([dgetrs_ sgetrs_ dstev_ sstev_])
+AC_CHECK_FUNCS([dgetrs_ sgetrs_ dstev_ sstev_ dsygv_ ssygv_])
 
 dnl -------------------------------------------------------------
 dnl       set include paths of several libraries
index ea67a684c4b80037861d1645a026554b263b7e28..15ed5ecc6907d297f904f7216ff9cb739785d6a7 100644 (file)
@@ -142,6 +142,17 @@ void dsyev_ (const char *jobz, const char *uplo, const int *n,
 void ssyev_ (const char *jobz, const char *uplo, const int *n,
             float *A, const int *lda, float *w,
             float *work, const int *lwork, int *info);
+// General eigenvalues and eigenvectors of 
+// 1: A*x = lambda*B*x; 2: A*B*x = lambda*x; 3: B*A*x = lambda*x
+// A and B are symmetric and B is definite
+void dsygv_ (const int* itype, const char* jobz, const char* uplo,
+             const int* n, double* A, const int* lda, double* B,
+             const int* ldb, double* w, double* work,
+             const int* lwork, int* info);
+void ssygv_ (const int* itype, const char* jobz, const char* uplo,
+             const int* n, float* A, const int* lda, float* B,
+             const int* ldb, float* w, float* work,
+             const int* lwork, int* info);
 // Compute singular value decomposition
 void dgesvd_ (int* jobu, int* jobvt,
              const int* n, const int* m, double* A, const int* lda,
@@ -562,6 +573,36 @@ syev (const char *, const char *, const int *, float *, const int *, float *, fl
 #endif
 
 
+#ifdef HAVE_DSYGV_
+inline void
+sygv (const int* itype, const char* jobz, const char* uplo, const int* n, double* A, const int* lda, double* B, const int* ldb, double* w, double* work, const int* lwork, int* info)
+{
+  dsygv_(itype, jobz, uplo, n, A, lda, B, ldb, w, work, lwork, info);
+}
+#else
+inline void
+sygv (const int*, const char*, const char*, const int*, double*, const int*, double*, const int*, double*, double*, const int*, int*)
+{
+  Assert (false, LAPACKSupport::ExcMissing("dsygv"));
+}
+#endif
+
+
+#ifdef HAVE_SSYGV_
+inline void
+sygv (const int* itype, const char* jobz, const char* uplo, const int* n, float* A, const int* lda, float* B, const int* ldb, float* w, float* work, const int* lwork, int* info)
+{
+  ssygv_(itype, jobz, uplo, n, A, lda, B, ldb, w, work, lwork, info);
+}
+#else
+inline void
+sygv (const int*, const char*, const char*, const int*, float*, const int*, float*, const int*, float*, float*, const int*, int*)
+{
+  Assert (false, LAPACKSupport::ExcMissing("ssygv"));
+}
+#endif
+
+
 #ifdef HAVE_DGESVD_
 inline void
 gesvd (int* jobu, int* jobvt, const int* n, const int* m, double* A, const int* lda, double* s, double* u, const int* ldu, double* vt, const int* ldvt, double* work, const int* lwork, int* info)

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.