]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Remove lapack_templates.h from more public headers. 10287/head
authorDavid Wells <drwells@email.unc.edu>
Wed, 20 May 2020 19:01:22 +0000 (15:01 -0400)
committerDavid Wells <drwells@email.unc.edu>
Wed, 20 May 2020 19:10:19 +0000 (15:10 -0400)
As mentioned previously, this file can conflict with other packages' LAPACK
wrapper headers. Fortunately we only currently use three LAPACK functions in
headers - its easy enough to special-case these so that we can get rid of the
header inclusion.

Fixes #9135.

include/deal.II/lac/qr.h
include/deal.II/lac/utilities.h
source/lac/CMakeLists.txt
source/lac/qr.cc [new file with mode: 0644]
source/lac/utilities.cc [new file with mode: 0644]

index 69ceeeda321aa842eab1746c58b3bc2017046eb8..343038c482f8882c3f4a337b66fecb28c9fe56d5 100644 (file)
@@ -22,7 +22,6 @@
 #include <deal.II/base/std_cxx14/memory.h>
 
 #include <deal.II/lac/lapack_full_matrix.h>
-#include <deal.II/lac/lapack_templates.h>
 #include <deal.II/lac/utilities.h>
 
 #include <boost/signals2/signal.hpp>
@@ -433,6 +432,41 @@ private:
 // -------------------  inline and template functions ----------------
 #ifndef DOXYGEN
 
+namespace internal
+{
+  namespace QRImplementation
+  {
+    // We want to avoid including our own LAPACK wrapper header in any external
+    // headers to avoid possible conflicts with other packages that may define
+    // their own such header. At the same time we want to be able to call some
+    // LAPACK functions from the template functions below. To resolve both
+    // problems define some extra wrappers here that can be in the header:
+    template <typename Number>
+    void
+    call_trmv(const char            uplo,
+              const char            trans,
+              const char            diag,
+              const types::blas_int n,
+              const Number *        a,
+              const types::blas_int lda,
+              Number *              x,
+              const types::blas_int incx);
+
+    template <typename Number>
+    void
+    call_trtrs(const char            uplo,
+               const char            trans,
+               const char            diag,
+               const types::blas_int n,
+               const types::blas_int nrhs,
+               const Number *        a,
+               const types::blas_int lda,
+               Number *              b,
+               const types::blas_int ldb,
+               types::blas_int *     info);
+  } // namespace QRImplementation
+} // namespace internal
+
 
 
 template <typename VectorType>
@@ -482,16 +516,16 @@ BaseQR<VectorType>::solve(Vector<Number> &      x,
   const int N     = this->current_size;
   const int n_rhs = 1;
   int       info  = 0;
-  trtrs("U",
-        transpose ? "T" : "N",
-        "N",
-        &N,
-        &n_rhs,
-        &this->R(0, 0),
-        &lda,
-        &x(0),
-        &ldb,
-        &info);
+  internal::QRImplementation::call_trtrs('U',
+                                         transpose ? 'T' : 'N',
+                                         'N',
+                                         N,
+                                         n_rhs,
+                                         &this->R(0, 0),
+                                         lda,
+                                         &x(0),
+                                         ldb,
+                                         &info);
 }
 
 
@@ -579,8 +613,8 @@ ImplicitQR<VectorType>::append_column(const VectorType &column)
       const int N     = this->current_size;
       const int n_rhs = 1;
       int       info  = 0;
-      trtrs(
-        "U", "T", "N", &N, &n_rhs, &this->R(0, 0), &lda, &u(0), &ldb, &info);
+      internal::QRImplementation::call_trtrs(
+        'U', 'T', 'N', N, n_rhs, &this->R(0, 0), lda, &u(0), ldb, &info);
 
       // finally get the diagonal element:
       // rho2 = |column|^2 - |u|^2
@@ -852,7 +886,8 @@ QR<VectorType>::multiply_with_A(VectorType &y, const Vector<Number> &x) const
   const int      N    = this->current_size;
   const int      lda  = N;
   const int      incx = 1;
-  trmv("U", "N", "N", &N, &this->R(0, 0), &lda, &x1[0], &incx);
+  internal::QRImplementation::call_trmv(
+    'U', 'N', 'N', N, &this->R(0, 0), lda, &x1[0], incx);
 
   multiply_with_Q(y, x1);
 }
@@ -868,7 +903,8 @@ QR<VectorType>::multiply_with_AT(Vector<Number> &y, const VectorType &x) const
   const int N    = this->current_size;
   const int lda  = N;
   const int incx = 1;
-  trmv("U", "T", "N", &N, &this->R(0, 0), &lda, &y[0], &incx);
+  internal::QRImplementation::call_trmv(
+    'U', 'T', 'N', N, &this->R(0, 0), lda, &y[0], incx);
 }
 
 #endif // no DOXYGEN
index 24f7e8b05b7af8b134cc61613f6c33384191a3f1..34c84f103a982442be60978967978e48ca332192 100644 (file)
@@ -22,7 +22,6 @@
 #include <deal.II/base/signaling_nan.h>
 
 #include <deal.II/lac/lapack_support.h>
-#include <deal.II/lac/lapack_templates.h>
 #include <deal.II/lac/vector_memory.h>
 
 #include <array>
@@ -211,6 +210,28 @@ namespace Utilities
 
 #ifndef DOXYGEN
 
+namespace internal
+{
+  namespace UtilitiesImplementation
+  {
+    // We want to avoid including our own LAPACK wrapper header in any external
+    // headers to avoid possible conflicts with other packages that may define
+    // their own such header. At the same time we want to be able to call some
+    // LAPACK functions from the template functions below. To resolve both
+    // problems define some extra wrappers here that can be in the header:
+    template <typename Number>
+    void
+    call_stev(const char            jobz,
+              const types::blas_int n,
+              Number *              d,
+              Number *              e,
+              Number *              z,
+              const types::blas_int ldz,
+              Number *              work,
+              types::blas_int *     info);
+  } // namespace UtilitiesImplementation
+} // namespace internal
+
 namespace Utilities
 {
   namespace LinearAlgebra
@@ -380,14 +401,14 @@ namespace Utilities
       std::vector<double>   work;    // ^^
       types::blas_int       info;
       // call lapack_templates.h wrapper:
-      stev("N",
-           &n,
-           diagonal.data(),
-           subdiagonal.data(),
-           Z.data(),
-           &ldz,
-           work.data(),
-           &info);
+      internal::UtilitiesImplementation::call_stev('N',
+                                                   n,
+                                                   diagonal.data(),
+                                                   subdiagonal.data(),
+                                                   Z.data(),
+                                                   ldz,
+                                                   work.data(),
+                                                   &info);
 
       Assert(info == 0, LAPACKSupport::ExcErrorCode("dstev", info));
 
index 07e08b9851373c87826dce168b6b0fc2b2167d7c..4a90cc4f475bfde3552794f60c25302e219ab800 100644 (file)
@@ -53,9 +53,11 @@ SET(_unity_include_src
 SET(_separate_src
   full_matrix.cc
   lapack_full_matrix.cc
+  qr.cc
   sparse_matrix.cc
   sparse_matrix_inst2.cc
   tridiagonal_matrix.cc
+  utilities.cc
   )
 
 # Add CUDA wrapper files
diff --git a/source/lac/qr.cc b/source/lac/qr.cc
new file mode 100644 (file)
index 0000000..10612fd
--- /dev/null
@@ -0,0 +1,150 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2020 - 2020 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#include <deal.II/base/config.h>
+
+#include <deal.II/lac/lapack_templates.h>
+#include <deal.II/lac/qr.h>
+
+#include <complex>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace internal
+{
+  namespace QRImplementation
+  {
+    // see the corresponding note in the header
+    template <typename Number>
+    void
+    call_trmv(const char            uplo,
+              const char            trans,
+              const char            diag,
+              const types::blas_int n,
+              const Number *        a,
+              const types::blas_int lda,
+              Number *              x,
+              const types::blas_int incx)
+    {
+      trmv(&uplo, &trans, &diag, &n, a, &lda, x, &incx);
+    }
+
+    template <typename Number>
+    void
+    call_trtrs(const char            uplo,
+               const char            trans,
+               const char            diag,
+               const types::blas_int n,
+               const types::blas_int nrhs,
+               const Number *        a,
+               const types::blas_int lda,
+               Number *              b,
+               const types::blas_int ldb,
+               types::blas_int *     info)
+    {
+      trtrs(&uplo, &trans, &diag, &n, &nrhs, a, &lda, b, &ldb, info);
+    }
+
+    template void
+    call_trmv(const char,
+              const char,
+              const char,
+              const types::blas_int,
+              const float *,
+              const types::blas_int,
+              float *,
+              const types::blas_int);
+
+    template void
+    call_trmv(const char,
+              const char,
+              const char,
+              const types::blas_int,
+              const double *,
+              const types::blas_int,
+              double *,
+              const types::blas_int);
+
+    template void
+    call_trmv(const char,
+              const char,
+              const char,
+              const types::blas_int,
+              const std::complex<float> *,
+              const types::blas_int,
+              std::complex<float> *,
+              const types::blas_int);
+
+    template void
+    call_trmv(const char,
+              const char,
+              const char,
+              const types::blas_int,
+              const std::complex<double> *,
+              const types::blas_int,
+              std::complex<double> *,
+              const types::blas_int);
+
+    template void
+    call_trtrs(const char,
+               const char,
+               const char,
+               const types::blas_int,
+               const types::blas_int,
+               const float *,
+               const types::blas_int,
+               float *,
+               const types::blas_int,
+               types::blas_int *);
+
+    template void
+    call_trtrs(const char,
+               const char,
+               const char,
+               const types::blas_int,
+               const types::blas_int,
+               const double *,
+               const types::blas_int,
+               double *,
+               const types::blas_int,
+               types::blas_int *);
+
+    template void
+    call_trtrs(const char,
+               const char,
+               const char,
+               const types::blas_int,
+               const types::blas_int,
+               const std::complex<float> *,
+               const types::blas_int,
+               std::complex<float> *,
+               const types::blas_int,
+               types::blas_int *);
+
+    template void
+    call_trtrs(const char,
+               const char,
+               const char,
+               const types::blas_int,
+               const types::blas_int,
+               const std::complex<double> *,
+               const types::blas_int,
+               std::complex<double> *,
+               const types::blas_int,
+               types::blas_int *);
+  } // namespace QRImplementation
+} // namespace internal
+
+DEAL_II_NAMESPACE_CLOSE
diff --git a/source/lac/utilities.cc b/source/lac/utilities.cc
new file mode 100644 (file)
index 0000000..7831869
--- /dev/null
@@ -0,0 +1,89 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2020 - 2020 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+#include <deal.II/base/config.h>
+
+#include <deal.II/lac/lapack_templates.h>
+#include <deal.II/lac/utilities.h>
+
+#include <complex>
+
+DEAL_II_NAMESPACE_OPEN
+
+namespace internal
+{
+  namespace UtilitiesImplementation
+  {
+    // see the corresponding note in the header
+    template <typename Number>
+    void
+    call_stev(const char            jobz,
+              const types::blas_int n,
+              Number *              d,
+              Number *              e,
+              Number *              z,
+              const types::blas_int ldz,
+              Number *              work,
+              types::blas_int *     info)
+    {
+      stev(&jobz, &n, d, e, z, &ldz, work, info);
+    }
+
+
+    template void
+    call_stev(const char,
+              const types::blas_int,
+              float *,
+              float *,
+              float *,
+              const types::blas_int,
+              float *,
+              types::blas_int *);
+
+    template void
+    call_stev(const char,
+              const types::blas_int,
+              double *,
+              double *,
+              double *,
+              const types::blas_int,
+              double *,
+              types::blas_int *);
+
+    template void
+    call_stev(const char,
+              const types::blas_int,
+              std::complex<float> *,
+              std::complex<float> *,
+              std::complex<float> *,
+              const types::blas_int,
+              std::complex<float> *,
+              types::blas_int *);
+
+    template void
+    call_stev(const char,
+              const types::blas_int,
+              std::complex<double> *,
+              std::complex<double> *,
+              std::complex<double> *,
+              const types::blas_int,
+              std::complex<double> *,
+              types::blas_int *);
+  } // namespace UtilitiesImplementation
+} // namespace internal
+
+
+
+DEAL_II_NAMESPACE_CLOSE

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.