DirichletBC.inc.hpp 1.92 KB
Newer Older
1
2
3
4
#pragma once

#include <type_traits>

5
#include <dune/functions/functionspacebases/subentitydofs.hh>
6

7
#include <amdis/functions/Interpolate.hpp>
8
#include <amdis/LinearAlgebra.hpp>
9
10

namespace AMDiS {
11
12
13
namespace Impl {

/// Traverse all interior boundary DOFs and apply the functor f.
14
15
/// The Functor must have the signature
/// `void(int, typename Basis::LocalView, typename Basis::GridView::Intersection)`
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
template <class Basis, class F>
void forEachInteriorBoundaryDOF(Basis const& basis, F&& f)
{
  auto localView = basis.localView();
  auto seDOFs = Dune::Functions::subEntityDOFs(basis);
  auto const& gridView = basis.gridView();
  for (auto const& element : elements(gridView, typename BackendTraits<Basis>::PartitionSet{})) {
    if (element.hasBoundaryIntersections()) {
      localView.bind(element);
      for(auto const& intersection: intersections(gridView, element))
        if (intersection.boundary())
          for(auto localIndex: seDOFs.bind(localView,intersection))
            f(localIndex, localView, intersection);
    }
  }
}

} // end namespace Impl

35

36
37
template <class Mat, class Sol, class Rhs, class RB, class CB, class V>
void DirichletBC<Mat, Sol, Rhs, RB, CB, V>::
38
init()
39
{
40
  using LV = typename CB::LocalView;
41
42
43
44
45
  dirichletNodes_.assign(colBasis_.dimension(), false);
  Impl::forEachInteriorBoundaryDOF(colBasis_,
    [&](int localIndex, LV const& localView, Intersection const& intersection) {
      dirichletNodes_[localView.index(localIndex)] = dirichletNodes_[localView.index(localIndex)]
                                                     || boundarySubset_(intersection);
46
47
48
49
  });
}


50
51
template <class Mat, class Sol, class Rhs, class RB, class CB, class V>
void DirichletBC<Mat, Sol, Rhs, RB, CB, V>::
52
apply(Mat& matrix, Sol& solution, Rhs& rhs)
53
{
54
  AMDiS::interpolate(colBasis_, solution, values_, tag::defaulted{}, tag::defaulted{}, dirichletNodes_);
55
  dirichletBC(matrix, solution, rhs, dirichletNodes_);
56
57
58
}

} // end namespace AMDiS