Go to the documentation of this file.
2 #ifndef DUNE_PDELAB_BACKEND_ISTL_CG_TO_DG_PROLONGATION_HH
3 #define DUNE_PDELAB_BACKEND_ISTL_CG_TO_DG_PROLONGATION_HH
5 #include <dune/common/exceptions.hh>
6 #include <dune/common/fvector.hh>
8 #include <dune/geometry/quadraturerules.hh>
10 #include <dune/localfunctions/common/interfaceswitch.hh>
12 #include <dune/typetree/pairtraversal.hh>
13 #include <dune/typetree/transformation.hh>
14 #include <dune/typetree/visitor.hh>
29 namespace CG2DGHelper {
31 template <
typename Imp>
33 typedef typename Imp::Traits::FiniteElementType::
34 Traits::LocalBasisType::Traits::RangeType
RangeType;
35 typedef typename Imp::Traits::FiniteElementType::
40 template<
typename Imp>
46 typedef typename Imp::Traits::FiniteElementType FEM;
47 typedef FiniteElementInterfaceSwitch<FEM> FESwitch;
48 typedef BasisInterfaceSwitch<typename FESwitch::Basis > BasisSwitch;
49 typedef typename BasisSwitch::DomainField DF;
50 typedef typename BasisSwitch::Range RT;
51 enum { dim = BasisSwitch::dimDomainLocal };
55 _imp(imp), _comp(comp) {}
57 void evaluate(
const Dune::FieldVector<DF,dim> & x,
58 Dune::FieldVector<DF,1> & y)
const
61 _imp.finiteElement().localBasis().evaluateFunction(x,v);
68 public TypeTree::DefaultPairVisitor,
69 public TypeTree::DynamicTraversal,
70 public TypeTree::VisitTree
79 template<
typename LFSU,
typename LFSV,
typename TreePath>
80 void leaf(
const LFSU& lfsu,
const LFSV& lfsv, TreePath treePath)
const
83 typedef typename LFSV::Traits::FiniteElementType DG_FEM;
84 typedef FiniteElementInterfaceSwitch<DG_FEM> FESwitch;
85 typedef BasisInterfaceSwitch<typename FESwitch::Basis > BasisSwitch;
86 typedef typename BasisSwitch::DomainField DF;
88 for (
unsigned int i=0; i<lfsu.size(); i++)
93 FESwitch::interpolation(lfsv.finiteElement()).
96 for (
unsigned int j=0; j<lfsv.size(); j++)
98 _mat(lfsv,j,lfsu,i) = v[j];
112 template<
typename LFSU,
typename LFSV,
typename R>
113 void computeCG2DG(
const LFSU & lfsu,
const LFSV & lfsv,
119 TypeTree::applyToTreePair(lfsu, lfsv, cg2dg);
135 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
typename M>
139 computeCG2DG(lfsu, lfsv, mat.container());
147 #endif // DUNE_PDELAB_BACKEND_ISTL_CG_TO_DG_PROLONGATION_HH
void jacobian_volume(const EG &, const LFSU &lfsu, const X &, const LFSV &lfsv, M &mat) const
Definition: cg_to_dg_prolongation.hh:136
@ doPatternVolume
Definition: cg_to_dg_prolongation.hh:123
Definition: cg_to_dg_prolongation.hh:67
ComputeCG2DGVisitor(LocalMatrix< R > &mat)
Definition: cg_to_dg_prolongation.hh:75
Imp::Traits::FiniteElementType::Traits::LocalBasisType::Traits::DomainType DomainType
Definition: cg_to_dg_prolongation.hh:36
Default class for additional methods in instationary local operators.
Definition: idefault.hh:89
CG2DGProlongation()
Definition: cg_to_dg_prolongation.hh:128
Definition: cg_to_dg_prolongation.hh:41
Default flags for all local operators.
Definition: flags.hh:18
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
void leaf(const LFSU &lfsu, const LFSV &lfsv, TreePath treePath) const
Definition: cg_to_dg_prolongation.hh:80
WrappedLocalShapeFunction(const Imp &imp, int comp)
Definition: cg_to_dg_prolongation.hh:54
Definition: cg_to_dg_prolongation.hh:107
sparsity pattern generator
Definition: pattern.hh:13
WrappedLocalShapeFunction< Imp > Traits
Definition: cg_to_dg_prolongation.hh:53
Imp::Traits::FiniteElementType::Traits::LocalBasisType::Traits::RangeType RangeType
Definition: cg_to_dg_prolongation.hh:34
Definition: cg_to_dg_prolongation.hh:32
void interpolate(const F &f, const GFS &gfs, XG &xg)
interpolation from a given grid function
Definition: interpolate.hh:204
void evaluate(const Dune::FieldVector< DF, dim > &x, Dune::FieldVector< DF, 1 > &y) const
Definition: cg_to_dg_prolongation.hh:57
@ doAlphaVolume
Definition: cg_to_dg_prolongation.hh:126