1#ifndef DUNE_PDELAB_BACKEND_ISTL_GENEO_LIPTONBABUSKABASIS_HH
2#define DUNE_PDELAB_BACKEND_ISTL_GENEO_LIPTONBABUSKABASIS_HH
19 template<
class GFS,
class M,
class X,
class Y,
int dim>
20 class LiptonBabuskaBasis :
public SubdomainBasis<X>
26 LiptonBabuskaBasis(
const GFS& gfs,
const M& AF_exterior,
const M& AF_ovlp,
const double eigenvalue_threshold, X& part_unity,
27 int& nev,
int nev_arpack,
double shift = 0.001,
bool add_part_unity =
false,
int verbose = 0) {
31 nev_arpack = std::max(nev, 2);
33 DUNE_THROW(Dune::Exception,
"eigenvectors_compute is less then eigenvectors or not specified!");
35 auto AF_interior = AF_exterior;
38 ArpackGeneo::ArPackPlusPlus_Algorithms<ISTLM, X> arpack(
native(AF_exterior));
41 std::vector<double> eigenvalues(nev_arpack,0.0);
42 std::vector<X> eigenvectors(nev_arpack,X(gfs,0.0));
44 arpack.computeGenNonSymMinMagnitude(
native(AF_interior), eps, eigenvectors, eigenvalues, shift);
48 if (eigenvalue_threshold >= 0) {
49 for (
int i = 0; i < nev; i++) {
50 if (eigenvalues[i] > eigenvalue_threshold) {
56 std::cout <<
"Process " << gfs.gridView().comm().rank() <<
" picked " << cnt <<
" eigenvectors" << std::endl;
58 DUNE_THROW(Dune::Exception,
"No eigenvalue above threshold - not enough eigenvalues computed!");
63 this->local_basis.resize(cnt);
65 for (
int base_id = 0; base_id < cnt; base_id++) {
66 this->local_basis[base_id] = std::make_shared<X>(part_unity);
69 this->local_basis[base_id]->begin(),this->local_basis[base_id]->end(),
70 eigenvectors[base_id].begin(),
71 this->local_basis[base_id]->begin(),
77 for (
auto& v : this->local_basis) {
78 *v *= 1.0 / v->two_norm2();
81 if (add_part_unity && eigenvalues[0] > 1E-10) {
82 this->local_basis.insert (this->local_basis.begin(), std::make_shared<X>(part_unity));
83 this->local_basis.pop_back();
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
std::enable_if< std::is_base_of< impl::WrapperBase, T >::value, Native< T > & >::type native(T &t)
Definition: backend/interface.hh:192
typename native_type< T >::type Native
Alias of the native container type associated with T or T itself if it is not a backend wrapper.
Definition: backend/interface.hh:176