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

#include <type_traits>

#include <dune/common/hybridutilities.hh>
6
#include <dune/functions/functionspacebases/subentitydofs.hh>
7

8
#include <amdis/LinearAlgebra.hpp>
9
10
#include <amdis/functions/Interpolate.hpp>
#include <amdis/gridfunctions/GridFunction.hpp>
11
12

namespace AMDiS {
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
namespace Impl {

/// Traverse all interior boundary DOFs and apply the functor f.
/// The Functor must have the signature `void(int, typename Basis::LocalView, typename Basis::GridView::Intersection)`
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

36
37

template <class D, class R>
38
  template <class RB, class CB>
39
void DirichletBC<D,R>::
40
init(RB const& rowBasis, CB const& colBasis)
41
{
42
43
  using LV = typename CB::LocalView;
  using IS = typename CB::GridView::Intersection;
44
  dirichletNodes_.resize(colBasis.dimension());
45
46
  dirichletNodes_.assign(dirichletNodes_.size(), false);
  Impl::forEachInteriorBoundaryDOF(colBasis, [&](int localIndex, LV const& localView, IS const& intersection) {
47
    dirichletNodes_[localView.index(localIndex)] = dirichletNodes_[localView.index(localIndex)] || onBoundary(intersection);
48
49
50
51
52
  });
}


template <class D, class R>
53
  template <class Mat, class Sol, class Rhs, class RN, class RTP, class CN, class CTP>
54
void DirichletBC<D,R>::
55
fillBoundaryCondition(Mat& matrix, Sol& solution, Rhs& rhs, RN const& rowNode, RTP rowTreePath, CN const& colNode, CTP colTreePath)
56
{
57
  Dune::Hybrid::ifElse(std::is_same<RangeType_t<RN>, R>{},
58
    [&](auto id) {
59
60
      auto&& gf = makeGridFunction(values_, solution.basis()->gridView());
      AMDiS::interpolate(*solution.basis(), id(solution), gf, rowTreePath, tag::defaulted{}, dirichletNodes_);
61
62
    });

63
  dirichletBC(matrix, solution, rhs, dirichletNodes_);
64
65
66
}

} // end namespace AMDiS