AMDiS  0.1
The Adaptive Multi-Dimensional Simulation Toolbox
BiLinearForm.inc.hpp
1 #pragma once
2 
3 #include <utility>
4 
5 #include <amdis/Assembler.hpp>
6 #include <amdis/LocalOperator.hpp>
7 #include <amdis/typetree/Traversal.hpp>
8 #include <amdis/utility/AssembleOperators.hpp>
9 
10 namespace AMDiS {
11 
12 template <class RB, class CB, class T, class Traits>
13  template <class ContextTag, class Expr, class RowTreePath, class ColTreePath>
15 addOperator(ContextTag contextTag, Expr const& expr,
16  RowTreePath row, ColTreePath col)
17 {
18  static_assert( Concepts::PreTreePath<RowTreePath>,
19  "row must be a valid treepath, or an integer/index-constant");
20  static_assert( Concepts::PreTreePath<ColTreePath>,
21  "col must be a valid treepath, or an integer/index-constant");
22 
23  auto i = makeTreePath(row);
24  auto node_i = child(this->rowBasis()->localView().treeCache(), i);
25  auto j = makeTreePath(col);
26  auto node_j = child(this->colBasis()->localView().treeCache(), j);
27 
28  using LocalContext = typename ContextTag::type;
29  using Tr = DefaultAssemblerTraits<LocalContext, ElementMatrix>;
30  auto op = makeLocalOperator<LocalContext>(expr, this->rowBasis()->gridView());
31  auto localAssembler = makeUniquePtr(makeAssembler<Tr>(std::move(op), node_i, node_j));
32 
33  operators_[i][j].push(contextTag, std::move(localAssembler));
34  updatePattern_ = true;
35 }
36 
37 
38 template <class RB, class CB, class T, class Traits>
39  template <class RowLocalView, class ColLocalView,
40  REQUIRES_(Concepts::LocalView<RowLocalView>),
41  REQUIRES_(Concepts::LocalView<ColLocalView>)>
43 assemble(RowLocalView const& rowLocalView, ColLocalView const& colLocalView)
44 {
45  elementMatrix_.resize(rowLocalView.size(), colLocalView.size());
46  elementMatrix_ = 0;
47 
48  auto const& gv = this->rowBasis()->gridView();
49  auto const& element = rowLocalView.element();
50  auto geometry = element.geometry();
51 
52  for_each_node(rowLocalView.treeCache(), [&](auto const& rowNode, auto rowTp) {
53  for_each_node(colLocalView.treeCache(), [&](auto const& colNode, auto colTp) {
54  auto& matOp = operators_[rowTp][colTp];
55  if (matOp) {
56  matOp.bind(element, geometry);
57  assembleOperators(gv, element, matOp, makeMatrixAssembler(rowNode, colNode, elementMatrix_));
58  matOp.unbind();
59  }
60  });
61  });
62 
63  this->scatter(rowLocalView, colLocalView, elementMatrix_);
64 }
65 
66 
67 template <class RB, class CB, class T, class Traits>
70 {
71  auto rowLocalView = this->rowBasis()->localView();
72  auto colLocalView = this->colBasis()->localView();
73 
74  this->init();
75  for (auto const& element : elements(this->rowBasis()->gridView(), typename Traits::PartitionSet{})) {
76  rowLocalView.bind(element);
77  if (this->rowBasis() == this->colBasis())
78  this->assemble(rowLocalView, rowLocalView);
79  else {
80  colLocalView.bind(element);
81  this->assemble(rowLocalView, colLocalView);
82  colLocalView.unbind();
83  }
84  rowLocalView.unbind();
85  }
86  this->finish();
87 }
88 
89 
90 } // end namespace AMDiS
void addOperator(ContextTag contextTag, Expr const &expr, Row row, Col col)
Associate a local operator with this BiLinearForm.
Contains all classes needed for solving linear and non linear equation systems.
Definition: AdaptBase.hpp:6
void assemble()
Assemble all matrix operators, TODO: incorporate boundary conditions.
Definition: BiLinearForm.inc.hpp:69
void for_each_node(Tree &&tree, PreFunc &&preFunc, LeafFunc &&leafFunc, PostFunc &&postFunc)
Traverse tree and visit each node.
Definition: Traversal.hpp:81
auto makeUniquePtr(Obj &&obj)
Create a unique_ptr by copy/move construction.
Definition: TypeTraits.hpp:107
void init(int &argc, char **&argv, std::string const &initFileName="")
Initialized the Environment for MPI.
Definition: AMDiS.hpp:29