LinearForm.inc.hpp 1.98 KB
Newer Older
1
2
3
4
#pragma once

#include <utility>

5
#include <amdis/GridFunctionOperator.hpp>
6
7
8
#include <amdis/LocalOperator.hpp>
#include <amdis/operations/Assigner.hpp>
#include <amdis/typetree/Traversal.hpp>
9
#include <amdis/typetree/TreePath.hpp>
10
11
12

namespace AMDiS {

Praetorius, Simon's avatar
Praetorius, Simon committed
13
template <class GB, class T, class Traits>
14
  template <class ContextTag, class Expr, class TreePath>
Praetorius, Simon's avatar
Praetorius, Simon committed
15
void LinearForm<GB,T,Traits>::
16
17
18
19
20
addOperator(ContextTag contextTag, Expr const& expr, TreePath path)
{
  static_assert( Concepts::PreTreePath<TreePath>,
      "path must be a valid treepath, or an integer/index-constant");

21
  auto i = makeTreePath(path);
22
23

  using LocalContext = typename ContextTag::type;
24
25
  auto op = makeOperator<LocalContext>(expr, this->basis().gridView());
  operators_[i].push(contextTag, std::move(op));
26
27
28
}


Praetorius, Simon's avatar
Praetorius, Simon committed
29
template <class GB, class T, class Traits>
30
31
  template <class LocalView, class LocalOperators,
    REQUIRES_(Concepts::LocalView<LocalView>)>
Praetorius, Simon's avatar
Praetorius, Simon committed
32
void LinearForm<GB,T,Traits>::
33
assemble(LocalView const& localView, LocalOperators& localOperators)
34
35
36
37
{
  elementVector_.resize(localView.size());
  elementVector_ = 0;

38
  auto const& gv = this->basis().gridView();
Praetorius, Simon's avatar
Praetorius, Simon committed
39
  GlobalContext<TYPEOF(gv)> context{gv, localView.element(), localView.element().geometry()};
40

41
  Traversal::forEachNode(localView.treeCache(), [&](auto const& node, auto tp) {
42
43
44
45
    auto& rhsOp = localOperators[tp];
    rhsOp.bind(context.element());
    rhsOp.assemble(context, node, elementVector_);
    rhsOp.unbind();
46
47
48
49
50
51
  });

  this->scatter(localView, elementVector_, Assigner::plus_assign{});
}


Praetorius, Simon's avatar
Praetorius, Simon committed
52
53
template <class GB, class T, class Traits>
void LinearForm<GB,T,Traits>::
54
55
assemble()
{
56
  auto localView = this->basis().localView();
57
  auto localOperators = AMDiS::localOperators(operators_);
58

59
  this->init(sizeInfo(this->basis()), true);
60
  for (auto const& element : elements(this->basis().gridView(), typename Traits::PartitionSet{})) {
61
    localView.bind(element);
62
    this->assemble(localView, localOperators);
63
64
65
66
67
68
    localView.unbind();
  }
  this->finish();
}

} // end namespace AMDiS