Go to the documentation of this file.
3 #ifndef DUNE_PDELAB_GRIDOPERATOR_FASTDG_HH
4 #define DUNE_PDELAB_GRIDOPERATOR_FASTDG_HH
6 #include <dune/common/tupleutility.hh>
29 template<
typename GFSU,
typename GFSV,
typename LOP,
30 typename MB,
typename DF,
typename RF,
typename JF,
54 GFSU::Traits::EntitySet::Partitions::partitionIterator() == InteriorBorder_Partition
58 typedef typename std::conditional<
59 GFSU::Traits::EntitySet::Partitions::partitionIterator() == InteriorBorder_Partition,
68 template <
typename MFT>
74 FastDGGridOperator(
const GFSU & gfsu_,
const CU & cu_,
const GFSV & gfsv_,
const CV & cv_, LOP & lop_,
const MB& mb_ = MB())
75 : global_assembler(gfsu_,gfsv_,cu_,cv_)
77 , local_assembler(lop_, cu_, cv_,dof_exchanger)
86 : global_assembler(gfsu_,gfsv_)
88 , local_assembler(lop_,dof_exchanger)
156 template <
typename Gr
idOperatorTuple>
162 template <
typename T>
164 elem.localAssembler().preProcessing(
index == 0);
165 elem.localAssembler().postProcessing(
index ==
size-1);
176 template<
typename Gr
idOperatorTuple>
180 Hybrid::forEach(tuple,
181 [&](
auto &el) { setup_visitor.
visit(el); });
185 template<
typename F,
typename X>
200 global_assembler.
assemble(pattern_engine);
208 global_assembler.
assemble(residual_engine);
216 global_assembler.
assemble(jacobian_engine);
223 DUNE_THROW(Dune::Exception,
"Your trying to use a linear jacobian apply for a non linear problem.");
231 DUNE_THROW(Dune::Exception,
"Your trying to use a non linear jacobian apply for a linear problem.");
236 void DUNE_DEPRECATED_MSG(
"nonlinear_jacobian_apply(x,z,r) is deprecated. Please use jacobian_apply(solution, update, result) instead!")
240 DUNE_THROW(Dune::Exception,
"Your trying to use a non linear jacobian apply for a linear problem.");
246 dof_exchanger->accumulateBorderEntries(*
this,a);
252 dof_exchanger->update(*
this);
263 std::shared_ptr<BorderDOFExchanger> dof_exchanger;
LocalJacobianApplyAssemblerEngine & localJacobianApplyAssemblerEngine(const typename Traits::Domain &update, typename Traits::Range &result)
Definition: fastdg/localassembler.hh:190
void jacobian_apply(const Domain &update, Range &result) const
Apply jacobian matrix to the vector update without explicitly assembling it.
Definition: fastdg.hh:220
void jacobian(const Domain &x, Jacobian &a) const
Assembler jacobian.
Definition: fastdg.hh:212
const Assembler & assembler() const
Definition: fastdg.hh:149
typename impl::BackendVectorSelector< GridFunctionSpace, FieldType >::Type Vector
alias of the return type of BackendVectorSelector
Definition: backend/interface.hh:106
const P & p
Definition: constraints.hh:148
void residual(const Domain &x, Range &r) const
Assemble residual.
Definition: fastdg.hh:204
static void setupGridOperators(GridOperatorTuple tuple)
Definition: fastdg.hh:177
void copy_nonconstrained_dofs(const CG &cg, const XG &xgin, XG &xgout)
Definition: constraints.hh:987
void const Domain Range &const result
Definition: fastdg.hh:238
Dune::PDELab::Backend::Vector< GFSU, DF > Domain
The type of the domain (solution).
Definition: fastdg.hh:41
Dune::PDELab::Backend::Matrix< MB, Domain, Range, JF > Jacobian
The type of the jacobian.
Definition: fastdg.hh:45
MB MatrixBackend
The matrix backend of the grid operator.
Definition: gridoperatorutilities.hh:51
void make_consistent(Jacobian &a) const
Definition: fastdg.hh:244
CU & trialGridFunctionSpaceConstraints()
Definition: fastdg.hh:107
const CV & testGridFunctionSpaceConstraints() const
Definition: fastdg.hh:120
void DUNE_DEPRECATED_MSG("nonlinear_jacobian_apply(x,z,r) is deprecated. Please use jacobian_apply(solution, update, result) instead!") nonlinear_jacobian_apply(const Domain &solution
Apply jacobian matrix to the vector update without explicitly assembling it.
CV & testGridFunctionSpaceConstraints()
Definition: fastdg.hh:116
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
const Traits::MatrixBackend & matrixBackend() const
Get the matrix backend for this grid operator.
Definition: fastdg.hh:256
Traits class for the grid operator.
Definition: gridoperatorutilities.hh:33
LocalJacobianAssemblerEngine & localJacobianAssemblerEngine(typename Traits::Jacobian &a, const typename Traits::Solution &x)
Definition: fastdg/localassembler.hh:180
const GFSU & trialGridFunctionSpace() const
Get the trial grid function space.
Definition: fastdg.hh:96
LocalResidualAssemblerEngine & localResidualAssemblerEngine(typename Traits::Residual &r, const typename Traits::Solution &x)
Definition: fastdg/localassembler.hh:170
const GFSU & trialGridFunctionSpace() const
Get the trial grid function space.
Definition: fastdg/assembler.hh:70
LocalAssembler & localAssembler() const
Definition: fastdg.hh:151
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139
const int size
Definition: fastdg.hh:170
LOP & localOperator()
get a reference to the local operator
Definition: fastdg/localassembler.hh:111
Dune::PDELab::Backend::Matrix< MB, Domain, Range, JF > Jacobian
The type of the jacobian.
Definition: gridoperatorutilities.hh:72
void update()
Definition: fastdg.hh:249
const CU & trialConstraints() const
get the constraints on the trial grid function space
Definition: assemblerutilities.hh:205
const GFSV & testGridFunctionSpace() const
Get the test grid function space.
Definition: fastdg.hh:102
std::conditional< GFSU::Traits::EntitySet::Partitions::partitionIterator()==InteriorBorder_Partition, NonOverlappingBorderDOFExchanger< FastDGGridOperator >, OverlappingBorderDOFExchanger< FastDGGridOperator > >::type BorderDOFExchanger
Definition: fastdg.hh:62
void const Domain & update
Definition: fastdg.hh:237
void interpolate(const F &f, const GFS &gfs, XG &xg)
interpolation from a given grid function
Definition: interpolate.hh:177
GFSU::Traits::SizeType globalSizeU() const
Get dimension of space u.
Definition: fastdg.hh:136
FastDGAssembler< GFSU, GFSV, CU, CV > Assembler
The global assembler type.
Definition: fastdg.hh:38
typename impl::BackendMatrixSelector< Backend, VU, VV, E >::Type Matrix
alias of the return type of BackendMatrixSelector
Definition: backend/interface.hh:127
void jacobian_apply(const Domain &solution, const Domain &update, Range &result) const
Apply jacobian matrix to the vector update without explicitly assembling it.
Definition: fastdg.hh:228
FastDGGridOperator(const GFSU &gfsu_, const GFSV &gfsv_, LOP &lop_, const MB &mb_=MB())
Constructor for empty constraints.
Definition: fastdg.hh:85
const GFSV & testGridFunctionSpace() const
Get the test grid function space.
Definition: fastdg/assembler.hh:76
void interpolate(const X &xold, F &f, X &x) const
Interpolate the constrained dofs from given function.
Definition: fastdg.hh:186
LOP & localOperator()
Definition: fastdg.hh:125
Definition: fastdg.hh:157
const LOP & localOperator() const
Definition: fastdg.hh:130
const CU & trialGridFunctionSpaceConstraints() const
Definition: fastdg.hh:111
Definition: constraintstransformation.hh:111
void assemble(LocalAssemblerEngine &assembler_engine) const
Definition: fastdg/assembler.hh:88
The local assembler for DUNE grids.
Definition: fastdg/localassembler.hh:37
Definition: borderdofexchanger.hh:577
SetupGridOperator()
Definition: fastdg.hh:159
void visit(T &elem)
Definition: fastdg.hh:163
Dune::PDELab::Backend::Vector< GFSV, RF > Range
The type of the range (residual).
Definition: fastdg.hh:43
Helper class for adding up matrix entries on border.
Definition: borderdofexchanger.hh:67
GFSV::Traits::SizeType globalSizeV() const
Get dimension of space v.
Definition: fastdg.hh:142
MB::template Pattern< Jacobian, GFSV, GFSU > Pattern
The sparsity pattern container for the jacobian matrix.
Definition: fastdg.hh:48
FastDGLocalAssembler< FastDGGridOperator, LOP, GFSU::Traits::EntitySet::Partitions::partitionIterator()==InteriorBorder_Partition > LocalAssembler
The local assembler type.
Definition: fastdg.hh:55
Assembler & assembler()
Definition: fastdg.hh:147
FastDGGridOperator(const GFSU &gfsu_, const CU &cu_, const GFSV &gfsv_, const CV &cv_, LOP &lop_, const MB &mb_=MB())
Constructor for non trivial constraints.
Definition: fastdg.hh:74
The fast DG assembler for standard DUNE grid.
Definition: fastdg/assembler.hh:25
int index
Definition: fastdg.hh:169
void fill_pattern(Pattern &p) const
Fill pattern of jacobian matrix.
Definition: fastdg.hh:196
Dune::PDELab::GridOperatorTraits< GFSU, GFSV, MB, DF, RF, JF, CU, CV, Assembler, LocalAssembler > Traits
The grid operator traits.
Definition: fastdg.hh:66
LocalPatternAssemblerEngine & localPatternAssemblerEngine(typename Traits::MatrixPattern &p)
Definition: fastdg/localassembler.hh:161
Traits::Jacobian Type
Definition: fastdg.hh:70