11
2- function AgFEMSpace (f:: SingleFieldFESpace ,bgcell_to_bgcellin:: AbstractVector ,g:: SingleFieldFESpace = f)
2+ function AgFEMSpace (
3+ f:: SingleFieldFESpace ,
4+ bgcell_to_bgcellin:: AbstractVector ,
5+ g:: SingleFieldFESpace = f,
6+ args... )
7+
38 @assert get_triangulation (f) === get_triangulation (g)
4- AgFEMSpace (f,bgcell_to_bgcellin,get_fe_basis (g),get_fe_dof_basis (g))
9+ AgFEMSpace (f,bgcell_to_bgcellin,get_fe_basis (g),get_fe_dof_basis (g),args ... )
510end
611
712# Note: cell is in fact bgcell in this function since f will usually be an ExtendedFESpace
813function AgFEMSpace (
914 f:: SingleFieldFESpace ,
1015 bgcell_to_bgcellin:: AbstractVector ,
1116 shfns_g:: CellField ,
12- dofs_g:: CellDof )
17+ dofs_g:: CellDof ,
18+ bgcell_to_gcell:: AbstractVector = 1 : length (bgcell_to_bgcellin))
1319
1420 # Triangulation made of active cells
1521 trian_a = get_triangulation (f)
@@ -21,6 +27,7 @@ function AgFEMSpace(
2127 bgcell_to_acell = glue. mface_to_tface
2228 acell_to_bgcellin = lazy_map (Reindex (bgcell_to_bgcellin),acell_to_bgcell)
2329 acell_to_acellin = collect (lazy_map (Reindex (bgcell_to_acell),acell_to_bgcellin))
30+ acell_to_gcell = lazy_map (Reindex (bgcell_to_gcell),acell_to_bgcell)
2431
2532 # Build shape funs of g by replacing local funs in cut cells by the ones at the root
2633 # This needs to be done with shape functions in the physical domain
@@ -41,7 +48,8 @@ function AgFEMSpace(
4148 acell_to_acellin,
4249 acell_to_dof_ids,
4350 acell_to_coeffs,
44- acell_to_proj)
51+ acell_to_proj,
52+ acell_to_gcell)
4553
4654 FESpaceWithLinearConstraints (aggdof_to_fdof,aggdof_to_dofs,aggdof_to_coeffs,f)
4755end
@@ -51,7 +59,8 @@ function _setup_agfem_constraints(
5159 acell_to_acellin,
5260 acell_to_dof_ids,
5361 acell_to_coeffs,
54- acell_to_proj)
62+ acell_to_proj,
63+ acell_to_gcell)
5564
5665 n_acells = length (acell_to_acellin)
5766 fdof_to_isagg = fill (true ,n_fdofs)
@@ -62,12 +71,16 @@ function _setup_agfem_constraints(
6271 acellin = acell_to_acellin[acell]
6372 iscut = acell != acellin
6473 dofs = getindex! (cache,acell_to_dof_ids,acell)
74+ gcell = acell_to_gcell[acell]
6575 for (ldof,dof) in enumerate (dofs)
6676 if dof > 0
6777 fdof = dof
68- fdof_to_isagg[fdof] = iscut && fdof_to_isagg[fdof]
69- fdof_to_acell[fdof] = acell
70- fdof_to_ldof[fdof] = ldof
78+ acell_dof = fdof_to_acell[fdof]
79+ if acell_dof == 0 || gcell > acell_to_gcell[acell_dof]
80+ fdof_to_acell[fdof] = acell
81+ fdof_to_isagg[fdof] = iscut && fdof_to_isagg[fdof]
82+ fdof_to_ldof[fdof] = ldof
83+ end
7184 end
7285 end
7386 end
0 commit comments