From 3cdd8624100a8465126926bb3355a51eb55b4348 Mon Sep 17 00:00:00 2001 From: kanschat Date: Fri, 15 Sep 2006 14:46:59 +0000 Subject: [PATCH] test for all LAPACK functions separately in order to cope with Petsc's LAPACK version git-svn-id: https://svn.dealii.org/trunk@13905 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/base/include/base/config.h.in | 30 ++++++ deal.II/configure | 111 +++++++++++++++++++++ deal.II/configure.in | 1 + deal.II/contrib/blastemplates/templates.pl | 29 +++++- deal.II/lac/include/lac/lapack_support.h | 12 ++- deal.II/lac/include/lac/lapack_templates.h | 86 +++++++++++++++- deal.II/lac/source/tridiagonal_matrix.cc | 10 -- 7 files changed, 261 insertions(+), 18 deletions(-) diff --git a/deal.II/base/include/base/config.h.in b/deal.II/base/include/base/config.h.in index 20b61d05b4..e1016ea246 100644 --- a/deal.II/base/include/base/config.h.in +++ b/deal.II/base/include/base/config.h.in @@ -209,6 +209,21 @@ /* 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 @@ -253,6 +268,21 @@ /* 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 header file. */ #undef HAVE_STDINT_H diff --git a/deal.II/configure b/deal.II/configure index 9cfc7e4a3c..43d6d275ce 100755 --- a/deal.II/configure +++ b/deal.II/configure @@ -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 declares $ac_func. + For example, HP-UX 11i 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 to if __STDC__ is defined, since + exists even on freestanding compilers. */ + +#ifdef __STDC__ +# include +#else +# include +#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 diff --git a/deal.II/configure.in b/deal.II/configure.in index 0b9b4838e6..ce8c038823 100644 --- a/deal.II/configure.in +++ b/deal.II/configure.in @@ -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 diff --git a/deal.II/contrib/blastemplates/templates.pl b/deal.II/contrib/blastemplates/templates.pl index a95724d87a..5de51ec6fe 100644 --- a/deal.II/contrib/blastemplates/templates.pl +++ b/deal.II/contrib/blastemplates/templates.pl @@ -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 + 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"; } } diff --git a/deal.II/lac/include/lac/lapack_support.h b/deal.II/lac/include/lac/lapack_support.h index bb9ad58ceb..967ac9e8b1 100644 --- a/deal.II/lac/include/lac/lapack_support.h +++ b/deal.II/lac/include/lac/lapack_support.h @@ -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"); } diff --git a/deal.II/lac/include/lac/lapack_templates.h b/deal.II/lac/include/lac/lapack_templates.h index 394801f0e0..5d817fb530 100644 --- a/deal.II/lac/include/lac/lapack_templates.h +++ b/deal.II/lac/include/lac/lapack_templates.h @@ -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 - - +#include 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 diff --git a/deal.II/lac/source/tridiagonal_matrix.cc b/deal.II/lac/source/tridiagonal_matrix.cc index 951537f72c..6a6b7a2ec0 100644 --- a/deal.II/lac/source/tridiagonal_matrix.cc +++ b/deal.II/lac/source/tridiagonal_matrix.cc @@ -237,16 +237,6 @@ TridiagonalMatrix::compute_eigenvalues() } -template<> -void -TridiagonalMatrix::compute_eigenvalues() -{ -// stev for float matrices should only be called if configure found -// sstev_ being included into liblapack - Assert(false, ExcNotImplemented()); -} - - template number TridiagonalMatrix::eigenvalue(const unsigned int i) const -- 2.39.5