#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>
// ------------------- 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>
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);
}
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
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);
}
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
#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>
#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
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));
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
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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