AMDiS 2.10
The Adaptive Multi-Dimensional Simulation Toolbox
DirichletBC.inc.hpp
1#pragma once
2
3#include <set>
4#include <type_traits>
5
6#include <dune/functions/functionspacebases/subentitydofs.hh>
7#include <dune/functions/functionspacebases/subspacebasis.hh>
8
9#include <amdis/LinearAlgebra.hpp>
10#include <amdis/functions/EntitySet.hpp>
11#include <amdis/interpolators/SimpleInterpolator.hpp>
12#include <amdis/typetree/TreePath.hpp>
13
14namespace AMDiS {
15
16template <class B, class Row, class Col, class V>
18init()
19{
20 if (row_ == col_) {
21 auto rowBasis = Dune::Functions::subspaceBasis(*basis_, col_);
22
23 std::set<typename B::MultiIndex> rowSet;
24
25 auto rowLocalView = rowBasis.localView();
26 auto rowDOFs = Dune::Functions::subEntityDOFs(rowBasis);
27 auto const& gridView = basis_->gridView();
28 for (auto const& element : entitySet(*basis_)) {
29 if (element.hasBoundaryIntersections()) {
30 rowLocalView.bind(element);
31 for(auto const& intersection: intersections(gridView, element)) {
32 if (intersection.boundary() && boundarySubset_(intersection)) {
33 for(auto localIndex: rowDOFs.bind(rowLocalView,intersection))
34 rowSet.insert(rowLocalView.index(localIndex));
35 }
36 }
37 }
38 }
39
40 rowIndices_.clear(); rowIndices_.reserve(rowSet.size());
41 rowIndices_.insert(rowIndices_.begin(), rowSet.begin(), rowSet.end());
42 colIndices_ = rowIndices_;
43 } else {
44 auto rowBasis = Dune::Functions::subspaceBasis(*basis_, row_);
45 auto colBasis = Dune::Functions::subspaceBasis(*basis_, col_);
46
47 std::set<typename B::MultiIndex> rowSet, colSet;
48
49 auto rowLocalView = rowBasis.localView();
50 auto colLocalView = colBasis.localView();
51 auto rowDOFs = Dune::Functions::subEntityDOFs(rowBasis);
52 auto colDOFs = Dune::Functions::subEntityDOFs(colBasis);
53 auto const& gridView = basis_->gridView();
54 for (auto const& element : entitySet(*basis_)) {
55 if (element.hasBoundaryIntersections()) {
56 rowLocalView.bind(element);
57 colLocalView.bind(element);
58 for(auto const& intersection: intersections(gridView, element)) {
59 if (intersection.boundary() && boundarySubset_(intersection)) {
60 for(auto localIndex: rowDOFs.bind(rowLocalView,intersection))
61 rowSet.insert(rowLocalView.index(localIndex));
62 for(auto localIndex: colDOFs.bind(colLocalView,intersection))
63 colSet.insert(colLocalView.index(localIndex));
64 }
65 }
66 }
67 }
68
69 rowIndices_.clear(); rowIndices_.reserve(rowSet.size());
70 rowIndices_.insert(rowIndices_.begin(), rowSet.begin(), rowSet.end());
71
72 colIndices_.clear(); colIndices_.reserve(colSet.size());
73 colIndices_.insert(colIndices_.begin(), colSet.begin(), colSet.end());
74 }
75}
76
77
78template <class B, class Row, class Col, class V>
79 template <class Mat, class Sol, class Rhs>
81apply(Mat& matrix, Sol& solution, Rhs& rhs, bool symmetric)
82{
83 std::vector<typename Sol::value_type> solutionValues;
84 solutionValues.reserve(colIndices_.size());
85 {
86 // create a copy the solution, because the valueGridFct_ might involve the solution vector.
87 Sol boundarySolution{solution};
88 valueOf(boundarySolution, col_).interpolate_noalias(valueGridFct_, tag::assign{});
89 boundarySolution.gather(colIndices_,solutionValues);
90 }
91
92 // copy boundary solution dirichlet data to solution vector
93 solution.scatter(colIndices_,solutionValues);
94 solution.finish();
95
96 if (symmetric)
97 matrix.zeroRowsColumns(rowIndices_, colIndices_, true, solution, rhs);
98 else
99 matrix.zeroRows(rowIndices_, colIndices_, true);
100
101 // copy boundary solution dirichlet data to rhs vector
102 rhs.scatter(rowIndices_,solutionValues);
103 rhs.finish();
104}
105
106} // end namespace AMDiS
void init()
Definition: DirichletBC.inc.hpp:18
void apply(Mat &matrix, Sol &solution, Rhs &rhs, bool symmetric=false)
Apply dirichlet BC to matrix and vector.
Definition: DirichletBC.inc.hpp:81
Definition: SimpleInterpolator.hpp:18