]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
test for all LAPACK functions separately in order to cope with Petsc's LAPACK version
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 15 Sep 2006 14:46:59 +0000 (14:46 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 15 Sep 2006 14:46:59 +0000 (14:46 +0000)
git-svn-id: https://svn.dealii.org/trunk@13905 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/base/include/base/config.h.in
deal.II/configure
deal.II/configure.in
deal.II/contrib/blastemplates/templates.pl
deal.II/lac/include/lac/lapack_support.h
deal.II/lac/include/lac/lapack_templates.h
deal.II/lac/source/tridiagonal_matrix.cc

index 20b61d05b49951602deaacf167f37158629a5ec6..e1016ea246dd1e9a104097e5729d776fbb82feb3 100644 (file)
 /* Define if the compiler provides __builtin_expect */
 #undef HAVE_BUILTIN_EXPECT
 
+/* Define to 1 if you have the `dgeevx_' function. */
+#undef HAVE_DGEEVX_
+
+/* Define to 1 if you have the `dgeev_' function. */
+#undef HAVE_DGEEV_
+
+/* Define to 1 if you have the `dgemv_' function. */
+#undef HAVE_DGEMV_
+
+/* Define to 1 if you have the `dgesvd_' function. */
+#undef HAVE_DGESVD_
+
+/* Define to 1 if you have the `dstev_' function. */
+#undef HAVE_DSTEV_
+
 /* Define to 1 if you have the `gethostname' function. */
 #undef HAVE_GETHOSTNAME
 
 /* Define if you have the rand_r function */
 #undef HAVE_RAND_R
 
+/* Define to 1 if you have the `sgeevx_' function. */
+#undef HAVE_SGEEVX_
+
+/* Define to 1 if you have the `sgeev_' function. */
+#undef HAVE_SGEEV_
+
+/* Define to 1 if you have the `sgemv_' function. */
+#undef HAVE_SGEMV_
+
+/* Define to 1 if you have the `sgesvd_' function. */
+#undef HAVE_SGESVD_
+
+/* Define to 1 if you have the `sstev_' function. */
+#undef HAVE_SSTEV_
+
 /* Define to 1 if you have the <stdint.h> header file. */
 #undef HAVE_STDINT_H
 
index 9cfc7e4a3c6d6d3a5c4a3d5587c4720bcfc64a61..43d6d275ce4bdebb822e8ce035c7ee75d102f736 100755 (executable)
@@ -12065,6 +12065,117 @@ fi;
 
 
 
+
+
+for ac_func in dgemv_ sgemv_ dgeev_ sgeev_ dgeevx_ sgeevx_ dgesvd_ sgesvd_ dstev_ sstev_
+do
+as_ac_var=`echo "ac_cv_func_$ac_func" | $as_tr_sh`
+echo "$as_me:$LINENO: checking for $ac_func" >&5
+echo $ECHO_N "checking for $ac_func... $ECHO_C" >&6
+if eval "test \"\${$as_ac_var+set}\" = set"; then
+  echo $ECHO_N "(cached) $ECHO_C" >&6
+else
+  cat >conftest.$ac_ext <<_ACEOF
+/* confdefs.h.  */
+_ACEOF
+cat confdefs.h >>conftest.$ac_ext
+cat >>conftest.$ac_ext <<_ACEOF
+/* end confdefs.h.  */
+/* Define $ac_func to an innocuous variant, in case <limits.h> declares $ac_func.
+   For example, HP-UX 11i <limits.h> declares gettimeofday.  */
+#define $ac_func innocuous_$ac_func
+
+/* System header to define __stub macros and hopefully few prototypes,
+    which can conflict with char $ac_func (); below.
+    Prefer <limits.h> to <assert.h> if __STDC__ is defined, since
+    <limits.h> exists even on freestanding compilers.  */
+
+#ifdef __STDC__
+# include <limits.h>
+#else
+# include <assert.h>
+#endif
+
+#undef $ac_func
+
+/* Override any gcc2 internal prototype to avoid an error.  */
+#ifdef __cplusplus
+extern "C"
+{
+#endif
+/* We use char because int might match the return type of a gcc2
+   builtin and then its argument prototype would still apply.  */
+char $ac_func ();
+/* The GNU C library defines this for functions which it implements
+    to always fail with ENOSYS.  Some functions are actually named
+    something starting with __ and the normal name is an alias.  */
+#if defined (__stub_$ac_func) || defined (__stub___$ac_func)
+choke me
+#else
+char (*f) () = $ac_func;
+#endif
+#ifdef __cplusplus
+}
+#endif
+
+int
+main ()
+{
+return f != $ac_func;
+  ;
+  return 0;
+}
+_ACEOF
+rm -f conftest.$ac_objext conftest$ac_exeext
+if { (eval echo "$as_me:$LINENO: \"$ac_link\"") >&5
+  (eval $ac_link) 2>conftest.er1
+  ac_status=$?
+  grep -v '^ *+' conftest.er1 >conftest.err
+  rm -f conftest.er1
+  cat conftest.err >&5
+  echo "$as_me:$LINENO: \$? = $ac_status" >&5
+  (exit $ac_status); } &&
+        { ac_try='test -z "$ac_cxx_werror_flag"
+                        || test ! -s conftest.err'
+  { (eval echo "$as_me:$LINENO: \"$ac_try\"") >&5
+  (eval $ac_try) 2>&5
+  ac_status=$?
+  echo "$as_me:$LINENO: \$? = $ac_status" >&5
+  (exit $ac_status); }; } &&
+        { ac_try='test -s conftest$ac_exeext'
+  { (eval echo "$as_me:$LINENO: \"$ac_try\"") >&5
+  (eval $ac_try) 2>&5
+  ac_status=$?
+  echo "$as_me:$LINENO: \$? = $ac_status" >&5
+  (exit $ac_status); }; }; then
+  eval "$as_ac_var=yes"
+else
+  echo "$as_me: failed program was:" >&5
+sed 's/^/| /' conftest.$ac_ext >&5
+
+eval "$as_ac_var=no"
+fi
+rm -f conftest.err conftest.$ac_objext \
+      conftest$ac_exeext conftest.$ac_ext
+fi
+echo "$as_me:$LINENO: result: `eval echo '${'$as_ac_var'}'`" >&5
+echo "${ECHO_T}`eval echo '${'$as_ac_var'}'`" >&6
+if test `eval echo '${'$as_ac_var'}'` = yes; then
+  cat >>confdefs.h <<_ACEOF
+#define `echo "HAVE_$ac_func" | $as_tr_cpp` 1
+_ACEOF
+
+fi
+done
+
+
+
+
+
+
+
+
+
 echo "$as_me:$LINENO: result: " >&5
 echo "${ECHO_T}" >&6
 echo "$as_me:$LINENO: result: ------------------ checking compiler flags ------------------" >&5
index 0b9b4838e6064ba1b5c6cc3f57bdf1a6de925058..ce8c0388234dc4b4ab673ee70256cc09a7bbc84a 100644 (file)
@@ -510,6 +510,7 @@ AC_ARG_WITH(lapack,
                           argument is given, use -llapack. Default is not to
                           use -llapack.],
   DEAL_II_WITH_LAPACK($withval))
+AC_CHECK_FUNCS([dgemv_ sgemv_ dgeev_ sgeev_ dgeevx_ sgeevx_ dgesvd_ sgesvd_ dstev_ sstev_])
 
 dnl -------------------------------------------------------------
 dnl       set include paths of several libraries
index a95724d87a966809ca03c41a82832b08fead33b6..5de51ec6fef3d959c00fb743e8282f3fda800b20 100644 (file)
@@ -2,7 +2,7 @@
 #    $Id$
 #    Version: $Name$
 #
-#    Copyright (C) 2005 by the deal authors
+#    Copyright (C) 2005, 2006 by the deal authors
 #
 #    This file is subject to QPL and may not be  distributed
 #    without copyright and license information. Please refer
@@ -35,6 +35,8 @@ print << 'EOT'
 #ifndef __LAPACK_TEMPLATES_H
 #define __LAPACK_TEMPLATES_H
 
+#include <lac/lapack_support.h>
+
 extern "C"
 {
 EOT
@@ -42,16 +44,22 @@ EOT
 
 while(<>)
 {
+    # Write comment lines literally
     if (m'^\s*//')
     {
        print;
        next;
     }
+    # Lines of the form 'typename functionname (...'
+    # where functionname is of the form d..._,
+    # that is a double precision LAPACK function
     if (m'\s*(\w+)\s+d(\w+)_\s*\(')
     {
        $double = $_;
        my $type = $1;
        my $name = $2;
+       my $capname = $name;
+       $capname =~ tr/[a-z]/[A-Z]/;
        while (<>)
        {
            $double .= $_;
@@ -64,16 +72,31 @@ while(<>)
        
        $double =~ m/\(([^\)]*)/;
        my $args = $1;
+       # The arglist for the C++ function
        $args =~ s/\s+/ /g;
+       # The arglist handed down to the FORTRAN function
        $args2 = $args;
+       # Fortunately, all arguments are pointers, so we can use the *
+       # to separate data type and argument name
        $args2 =~ s/\w+\*//g;
        $args2 =~ s/const//g;
        $args2 =~ s/\s//g;
+       # The arglist of the empty C++ function
+       $args0 = $args;
+       $args0 =~ s/\*[^,]*,/\*,/g;
+       $args0 =~ s/\*[^,]*$/\*/g;
+       
+       $templates .= "\n\n#ifdef HAVE_D$capname\_";
+       $templates .= "\ninline $type\n$name ($args)\n{\n  d$name\_ ($args2);\n}\n";
+       $templates .= "#else\ninline $type\n$name ($args0)\n";
+       $templates .= "{\n  LAPACKSupport::ExcMissing(\"d$name\");\n}\n#endif\n";
        
-       $templates .= "\n\ninline $type\n$name ($args)\n{\n  d$name\_ ($args2);\n}\n";
        $args =~ s/double/float/g;
        $type =~ s/double/float/g;
-       $templates .= "\n\ninline $type\n$name ($args)\n{\n  s$name\_ ($args2);\n}\n";
+       $templates .= "\n\n#ifdef HAVE_S$capname\_";
+       $templates .= "\ninline $type\n$name ($args)\n{\n  s$name\_ ($args2);\n}\n";
+       $templates .= "#else\ninline $type\n$name ($args0)\n";
+       $templates .= "{\n  LAPACKSupport::ExcMissing(\"s$name\");\n}\n#endif\n";
     }
 }
 
index bb9ad58ceba5b026d9eb51f436e833c13bad682e..967ac9e8b1c9a787928d3e55114f4c51ae477705 100644 (file)
@@ -114,7 +114,17 @@ namespace LAPACKSupport
   DeclException1(ExcState, State,
                 << "The function cannot be called while the matrix is in state "
                 << state_name(arg1));
-  
+
+                                  /**
+                                   * This exception is thrown if a
+                                   * certain function is not
+                                   * implemented in your LAPACK
+                                   * version.
+                                   */
+  DeclException1(ExcMissing, char*,
+                << "The function "
+                << arg1
+                << " required here is missing in your LAPACK installation");
 }
 
 
index 394801f0e02dfdd5fb1e9ddceb8aa223ce574d78..5d817fb5305bb8a8797658c61b2e97f0845f352c 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 by the deal authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -17,9 +17,7 @@
 #ifndef __LAPACK_TEMPLATES_H
 #define __LAPACK_TEMPLATES_H
 
-#include <base/config.h>
-
-
+#include <lac/lapack_support.h>
 
 extern "C"
 {
@@ -100,74 +98,154 @@ void sstev_ (const char* jobz, const int* n,
 
 
 
+#ifdef HAVE_DGEMV_
 inline void
 gemv (const char* trans, const int* m, const int* n, const double* alpha, const double* A, const int* lda, const double* x, const int* incx, const double* b, double* y, const int* incy)
 {
   dgemv_ (trans,m,n,alpha,A,lda,x,incx,b,y,incy);
 }
+#else
+inline void
+gemv (const char*, const int*, const int*, const double*, const double*, const int*, const double*, const int*, const double*, double*, const int*)
+{
+  LAPACKSupport::ExcMissing("dgemv");
+}
+#endif
 
 
+#ifdef HAVE_SGEMV_
 inline void
 gemv (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)
 {
   sgemv_ (trans,m,n,alpha,A,lda,x,incx,b,y,incy);
 }
+#else
+inline void
+gemv (const char*, const int*, const int*, const double*, const double*, const int*, const double*, const int*, const double*, double*, const int*)
+{
+  LAPACKSupport::ExcMissing("sgemv");
+}
+#endif
 
 
+#ifdef HAVE_DGEEV_
 inline void
 geev (const char* jobvl, const char* jobvr, const int* n, double* A, const int* lda, double* lambda_re, double* lambda_im, double* vl, const int* ldvl, double* vr, const int* ldva, double* work, const int* lwork, int* info)
 {
   dgeev_ (jobvl,jobvr,n,A,lda,lambda_re,lambda_im,vl,ldvl,vr,ldva,work,lwork,info);
 }
+#else
+inline void
+geev (const char*, const char*, const int*, double*, const int*, double*, double*, double*, const int*, double*, const int*, double*, const int*, int*)
+{
+  LAPACKSupport::ExcMissing("dgeev");
+}
+#endif
 
 
+#ifdef HAVE_SGEEV_
 inline void
 geev (const char* jobvl, const char* jobvr, const int* n, float* A, const int* lda, float* lambda_re, float* lambda_im, float* vl, const int* ldvl, float* vr, const int* ldva, float* work, const int* lwork, int* info)
 {
   sgeev_ (jobvl,jobvr,n,A,lda,lambda_re,lambda_im,vl,ldvl,vr,ldva,work,lwork,info);
 }
+#else
+inline void
+geev (const char*, const char*, const int*, double*, const int*, double*, double*, double*, const int*, double*, const int*, double*, const int*, int*)
+{
+  LAPACKSupport::ExcMissing("sgeev");
+}
+#endif
 
 
+#ifdef HAVE_DGEEVX_
 inline void
 geevx (const char* balanc, const char* jobvl, const char* jobvr, const char* sense, const int* n, double* A, const int* lda, double* lambda_re, double* lambda_im, double* vl, const int* ldvl, double* vr, const int* ldvr, int* ilo, int* ihi, double* scale, double* abnrm, double* rconde, double* rcondv, double* work, const int* lwork, int* iwork, int* info)
 {
   dgeevx_ (balanc,jobvl,jobvr,sense,n,A,lda,lambda_re,lambda_im,vl,ldvl,vr,ldvr,ilo,ihi,scale,abnrm,rconde,rcondv,work,lwork,iwork,info);
 }
+#else
+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*)
+{
+  LAPACKSupport::ExcMissing("dgeevx");
+}
+#endif
 
 
+#ifdef HAVE_SGEEVX_
 inline void
 geevx (const char* balanc, const char* jobvl, const char* jobvr, const char* sense, const int* n, float* A, const int* lda, float* lambda_re, float* lambda_im, float* vl, const int* ldvl, float* vr, const int* ldvr, int* ilo, int* ihi, float* scale, float* abnrm, float* rconde, float* rcondv, float* work, const int* lwork, int* iwork, int* info)
 {
   sgeevx_ (balanc,jobvl,jobvr,sense,n,A,lda,lambda_re,lambda_im,vl,ldvl,vr,ldvr,ilo,ihi,scale,abnrm,rconde,rcondv,work,lwork,iwork,info);
 }
+#else
+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*)
+{
+  LAPACKSupport::ExcMissing("sgeevx");
+}
+#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)
 {
   dgesvd_ (jobu,jobvt,n,m,A,lda,s,u,ldu,vt,ldvt,work,lwork,info);
 }
+#else
+inline void
+gesvd (int*, int*, const int*, const int*, double*, const int*, double*, double*, const int*, double*, const int*, double*, const int*, int*)
+{
+  LAPACKSupport::ExcMissing("dgesvd");
+}
+#endif
 
 
+#ifdef HAVE_SGESVD_
 inline void
 gesvd (int* jobu, int* jobvt, const int* n, const int* m, float* A, const int* lda, float* s, float* u, const int* ldu, float* vt, const int* ldvt, float* work, const int* lwork, int* info)
 {
   sgesvd_ (jobu,jobvt,n,m,A,lda,s,u,ldu,vt,ldvt,work,lwork,info);
 }
+#else
+inline void
+gesvd (int*, int*, const int*, const int*, double*, const int*, double*, double*, const int*, double*, const int*, double*, const int*, int*)
+{
+  LAPACKSupport::ExcMissing("sgesvd");
+}
+#endif
 
 
+#ifdef HAVE_DSTEV_
 inline void
 stev (const char* jobz, const int* n, double* d, double* e, double* z, const int* ldz, double* work, int* info)
 {
   dstev_ (jobz,n,d,e,z,ldz,work,info);
 }
+#else
+inline void
+stev (const char*, const int*, double*, double*, double*, const int*, double*, int*)
+{
+  LAPACKSupport::ExcMissing("dstev");
+}
+#endif
 
 
+#ifdef HAVE_SSTEV_
 inline void
 stev (const char* jobz, const int* n, float* d, float* e, float* z, const int* ldz, float* work, int* info)
 {
   sstev_ (jobz,n,d,e,z,ldz,work,info);
 }
+#else
+inline void
+stev (const char*, const int*, double*, double*, double*, const int*, double*, int*)
+{
+  LAPACKSupport::ExcMissing("sstev");
+}
+#endif
 
 
 #endif
index 951537f72c31f4b540739a93bd32cbb6b9c8fff7..6a6b7a2ec024c2aa76ac13288ce531a3ccffec77 100644 (file)
@@ -237,16 +237,6 @@ TridiagonalMatrix<double>::compute_eigenvalues()
 }
 
 
-template<>
-void
-TridiagonalMatrix<float>::compute_eigenvalues()
-{
-// stev for float matrices should only be called if configure found
-// sstev_ being included into liblapack
-  Assert(false, ExcNotImplemented());
-}
-
-
 template<typename number>
 number
 TridiagonalMatrix<number>::eigenvalue(const unsigned int i) const

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.