3#ifndef DUNE_PDELAB_GRIDOPERATOR_FASTDG_ASSEMBLER_HH
4#define DUNE_PDELAB_GRIDOPERATOR_FASTDG_ASSEMBLER_HH
6#include <dune/common/typetraits.hh>
24 template<
typename GFSU,
typename GFSV,
typename CU,
typename CV,
bool nonoverlapping_mode=false>
31 using Element =
typename EntitySet::Element;
42 typedef typename GFSU::Traits::SizeType
SizeType;
47 FastDGAssembler (
const GFSU& gfsu_,
const GFSV& gfsv_,
const CU& cu_,
const CV& cv_)
87 template<
class LocalAssemblerEngine>
90 const bool fast =
true;
95 const bool needs_constraints_caching = assembler_engine.needsConstraintsCaching(cu,cv);
97 LFSUCache lfsu_cache(lfsu,cu,needs_constraints_caching);
98 LFSVCache lfsv_cache(lfsv,cv,needs_constraints_caching);
99 LFSUCache lfsun_cache(lfsun,cu,needs_constraints_caching);
100 LFSVCache lfsvn_cache(lfsvn,cv,needs_constraints_caching);
116 auto entity_set = gfsu.entitySet();
117 auto& index_set = entity_set.indexSet();
120 for (
const auto& element : elements(entity_set,assembler_engine.partition()))
123 auto ids = index_set.uniqueIndex(element);
131 lfsv.bind(element,std::integral_constant<bool,fast>{});
141 lfsu.bind(element,std::integral_constant<bool,fast>{});
145 assembler_engine.
onBindLFSUV(eg,lfsu_cache,lfsv_cache);
154 if (require_uv_skeleton || require_v_skeleton ||
155 require_uv_boundary || require_v_boundary ||
156 require_uv_processor || require_v_processor)
159 unsigned int intersection_index = 0;
160 for(
const auto& intersection : intersections(entity_set,element))
166 auto intersection_type = std::get<0>(intersection_data);
167 auto& outside_element = std::get<1>(intersection_data);
169 switch (intersection_type)
175 if (require_uv_skeleton || require_v_skeleton)
179 auto idn = index_set.uniqueIndex(outside_element);
182 bool visit_face = ids > idn
183 or require_skeleton_two_sided
184 or (assembler_engine.partition() == Partitions::interiorBorder
185 and outside_element.partitionType() != InteriorEntity);
191 lfsvn.bind(outside_element,std::integral_constant<bool,fast>{});
192 lfsvn_cache.update();
200 if(require_uv_skeleton){
203 lfsun.bind(outside_element,std::integral_constant<bool,fast>{});
204 lfsun_cache.update();
208 lfsu_cache,lfsv_cache,
209 lfsun_cache,lfsvn_cache);
219 lfsu_cache,lfsv_cache,
220 lfsun_cache,lfsvn_cache);
230 if(require_uv_boundary || require_v_boundary )
236 if(require_uv_boundary){
244 if(require_uv_processor || require_v_processor )
250 if(require_uv_processor){
258 ++intersection_index;
262 if(require_uv_post_skeleton || require_v_post_skeleton){
266 if(require_uv_post_skeleton){
291 typename std::conditional<
292 std::is_same<CU,EmptyTransformation>::value,
296 typename std::conditional<
297 std::is_same<CV,EmptyTransformation>::value,
const IG & ig
Definition constraints.hh:149
For backward compatibility – Do not use this!
Definition adaptivity.hh:28
@ processor
processor boundary intersection (neighbor() == false && boundary() == false) or outside entity not in...
@ boundary
domain boundary intersection (neighbor() == false && boundary() == true)
@ skeleton
skeleton intersection (neighbor() == true && boundary() == false)
@ periodic
periodic boundary intersection (neighbor() == true && boundary() == true)
std::tuple< IntersectionType, typename EntitySet::Element > classifyIntersection(const EntitySet &entity_set, const Intersection &is)
Classifies the type of an intersection wrt to the passed EntitySet.
Definition intersectiontype.hh:37
Wrap element.
Definition geometrywrapper.hh:16
Wrap intersection.
Definition geometrywrapper.hh:57
Definition lfsindexcache.hh:979
Create a local function space from a global function space.
Definition localfunctionspace.hh:717
The local assembler engine which handles the integration parts as provided by the global assemblers.
Definition common/assembler.hh:33
void onUnbindLFSV(const EG &eg, const LFSV_S &lfsv_s)
void assembleUVVolume(const EG &eg, const LFSU &lfsu, const LFSV &lfsv)
bool requireSkeletonTwoSided() const
void onBindLFSUV(const EG &eg, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s)
void assembleVSkeleton(const IG &ig, const LFSV_S &lfsv_s, const LFSV_N &lfsv_n)
bool requireVSkeleton() const
void assembleUVProcessor(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s)
void assembleUVSkeleton(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s, const LFSU_N &lfsu_n, const LFSV_N &lfsv_n)
bool requireUVVolumePostSkeleton() const
bool requireVBoundary() const
void loadCoefficientsLFSUInside(const LFSU_S &lfsu_s)
void assembleVProcessor(const IG &ig, const LFSV_S &lfsv_s)
void preAssembly()
Called directly before assembling.
void assembleVVolumePostSkeleton(const EG &eg, const LFSV &lfsv)
bool assembleCell(const EG &eg)
void assembleUVVolumePostSkeleton(const EG &eg, const LFSU &lfsu, const LFSV &lfsv)
void onUnbindLFSUVOutside(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s, const LFSU_N &lfsu_n, const LFSV_N &lfsv_n)
void assembleVBoundary(const IG &ig, const LFSV_S &lfsv_s)
void onBindLFSV(const EG &eg, const LFSV_S &lfsv_s)
void onBindLFSVOutside(const IG &ig, const LFSV_S &lfsv_s, const LFSV_N &lfsv_n)
void onUnbindLFSUV(const EG &eg, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s)
bool requireUVSkeleton() const
void assembleVVolume(const EG &eg, const LFSV &lfsv)
void onUnbindLFSVOutside(const IG &ig, const LFSV_S &lfsv_s, const LFSV_N &lfsv_n)
void postAssembly()
Called last thing after assembling.
bool requireUVBoundary() const
void assembleUVBoundary(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s)
void onBindLFSUVOutside(const IG &ig, const LFSU_S &lfsu_s, const LFSV_S &lfsv_s, const LFSU_N &lfsu_n, const LFSV_N &lfsv_n)
void loadCoefficientsLFSUOutside(const LFSU_N &lfsu_n)
bool requireVVolumePostSkeleton() const
The fast DG assembler for standard DUNE grid.
Definition fastdg/assembler.hh:25
FastDGAssembler(const GFSU &gfsu_, const GFSV &gfsv_, const CU &cu_, const CV &cv_)
Definition fastdg/assembler.hh:47
GFSU::Traits::SizeType SizeType
Size type as used in grid function space.
Definition fastdg/assembler.hh:42
typename GFSU::Traits::EntitySet EntitySet
Definition fastdg/assembler.hh:30
const GFSV & testGridFunctionSpace() const
Get the test grid function space.
Definition fastdg/assembler.hh:76
typename EntitySet::Element Element
Definition fastdg/assembler.hh:31
void assemble(LocalAssemblerEngine &assembler_engine) const
Definition fastdg/assembler.hh:88
const GFSU & trialGridFunctionSpace() const
Get the trial grid function space.
Definition fastdg/assembler.hh:70
typename EntitySet::Intersection Intersection
Definition fastdg/assembler.hh:32
GFSV TestGridFunctionSpace
Definition fastdg/assembler.hh:38
GFSU TrialGridFunctionSpace
Definition fastdg/assembler.hh:37
FastDGAssembler(const GFSU &gfsu_, const GFSV &gfsv_)
Definition fastdg/assembler.hh:58
static const bool isGalerkinMethod
Static check on whether this is a Galerkin method.
Definition fastdg/assembler.hh:45