[distributed] missing cut_facets interfaces

src/Distributed/Distributed.jl

@@ -6,6 +6,7 @@ using PartitionedArrays
 using FillArrays
 
 using Gridap.Arrays
+using Gridap.CellData
 using Gridap.Geometry
 using Gridap.Helpers
 using Gridap.ReferenceFEs
@@ -18,13 +19,17 @@ using GridapEmbedded.Interfaces: ActiveInOrOut
 using GridapEmbedded.Interfaces: SubFacetTriangulation
 using GridapEmbedded.Interfaces: SubCellData
 using GridapEmbedded.Interfaces: SubFacetData
+using GridapEmbedded.Interfaces: AbstractEmbeddedDiscretization
 using GridapEmbedded.AgFEM: _touch_aggregated_cells!
 using GridapEmbedded.AgFEM: AggregateCutCellsByThreshold
+using GridapEmbedded.MomentFittedQuadratures: MomentFitted
 using Gridap.Geometry: AppendedTriangulation
+using Gridap.Geometry: get_face_to_parent_face
 using GridapDistributed: DistributedDiscreteModel
 using GridapDistributed: DistributedTriangulation
 using GridapDistributed: DistributedFESpace
 using GridapDistributed: DistributedSingleFieldFESpace
+using GridapDistributed: DistributedMeasure
 using GridapDistributed: add_ghost_cells
 using GridapDistributed: generate_gids
 using GridapDistributed: generate_cell_gids
@@ -33,11 +38,14 @@ using GridapDistributed: _find_vector_type
 import GridapEmbedded.AgFEM: aggregate
 import GridapEmbedded.AgFEM: AgFEMSpace
 import GridapEmbedded.Interfaces: cut
+import GridapEmbedded.Interfaces: cut_facets
 import GridapEmbedded.Interfaces: EmbeddedBoundary
 import GridapEmbedded.Interfaces: compute_bgfacet_to_inoutcut
 import GridapEmbedded.Interfaces: compute_bgcell_to_inoutcut
 import GridapEmbedded.CSG: get_geometry
 import Gridap.Geometry: Triangulation
+import Gridap.Geometry: SkeletonTriangulation
+import Gridap.Geometry: BoundaryTriangulation
 import Gridap.Geometry: get_background_model
 import GridapDistributed: local_views
 import GridapDistributed: remove_ghost_cells
@@ -46,4 +54,6 @@ include("DistributedDiscretizations.jl")
 
 include("DistributedAgFEM.jl")
 
+include("DistributedQuadratures.jl")
+
 end # module

src/Distributed/DistributedDiscretizations.jl

@@ -1,13 +1,13 @@
 
-struct DistributedEmbeddedDiscretization{Dp,T,A,B} <: GridapType
+struct DistributedEmbeddedDiscretization{A,B} <: GridapType
   discretizations::A
   model::B
   function DistributedEmbeddedDiscretization(
-    d::AbstractArray{<:EmbeddedDiscretization{Dp,T}},
-    model::DistributedDiscreteModel) where {Dp,T}
+    d::AbstractArray{<:AbstractEmbeddedDiscretization},
+    model::DistributedDiscreteModel)
     A = typeof(d)
     B = typeof(model)
-    new{Dp,T,A,B}(d,model)
+    new{A,B}(d,model)
   end
 end
 
@@ -31,8 +31,29 @@ function cut(cutter::Cutter,bgmodel::DistributedDiscreteModel,args...)
     cutgeo = cut(cutter,ownmodel,args...)
     change_bgmodel(cutgeo,bgmodel,own_to_local(gids))
   end
-  ls_to_bgcell_to_inoutcut = map(c->c.ls_to_bgcell_to_inoutcut,cuts)
-  _consistent!(ls_to_bgcell_to_inoutcut,gids)
+  consistent_bgcell_to_inoutcut!(cuts,gids)
+  DistributedEmbeddedDiscretization(cuts,bgmodel)
+end
+
+function cut_facets(bgmodel::DistributedDiscreteModel,args...)
+  cut_facets(LevelSetCutter(),bgmodel,args...)
+end
+
+function cut_facets(cutter::Cutter,bgmodel::DistributedDiscreteModel,args...)
+  D = map(num_dims,local_views(bgmodel)) |> PartitionedArrays.getany
+  cell_gids = get_cell_gids(bgmodel)
+  facet_gids = get_face_gids(bgmodel,D-1)
+  cuts = map(
+    local_views(bgmodel),
+    local_views(cell_gids),
+    local_views(facet_gids)) do bgmodel,cell_gids,facet_gids
+    ownmodel = remove_ghost_cells(bgmodel,cell_gids)
+    facet_to_pfacet = get_face_to_parent_face(ownmodel,D-1)
+    cutfacets = cut_facets(cutter,ownmodel,args...)
+    cutfacets = change_bgmodel(cutfacets,bgmodel,facet_to_pfacet)
+    remove_ghost_subfacets(cutfacets,facet_gids)
+  end
+  consistent_bgfacet_to_inoutcut!(cuts,facet_gids)
   DistributedEmbeddedDiscretization(cuts,bgmodel)
 end
 
@@ -54,6 +75,16 @@ function EmbeddedBoundary(cutgeo::DistributedEmbeddedDiscretization,args...)
   remove_ghost_cells(trian)
 end
 
+function SkeletonTriangulation(cutgeo::DistributedEmbeddedDiscretization,args...)
+  trian = distributed_embedded_triangulation(SkeletonTriangulation,cutgeo,args...)
+  remove_ghost_cells(trian)
+end
+
+function BoundaryTriangulation(cutgeo::DistributedEmbeddedDiscretization,args...)
+  trian = distributed_embedded_triangulation(BoundaryTriangulation,cutgeo,args...)
+  remove_ghost_cells(trian)
+end
+
 function distributed_embedded_triangulation(
   T,
   cutgeo::DistributedEmbeddedDiscretization,
@@ -139,6 +170,30 @@ function remove_ghost_cells(model::DiscreteModel,gids::AbstractLocalIndices)
   DiscreteModelPortion(model,own_to_local(gids))
 end
 
+function consistent_bgcell_to_inoutcut!(
+  cuts::AbstractArray{<:AbstractEmbeddedDiscretization},
+  gids::PRange)
+
+  ls_to_bgcell_to_inoutcut = map(get_ls_to_bgcell_to_inoutcut,cuts)
+  _consistent!(ls_to_bgcell_to_inoutcut,gids)
+end
+
+function get_ls_to_bgcell_to_inoutcut(cut::EmbeddedDiscretization)
+  cut.ls_to_bgcell_to_inoutcut
+end
+
+function consistent_bgfacet_to_inoutcut!(
+  cuts::AbstractArray{<:AbstractEmbeddedDiscretization},
+  gids::PRange)
+
+  ls_to_bgfacet_to_inoutcut = map(get_ls_to_bgfacet_to_inoutcut,cuts)
+  _consistent!(ls_to_bgfacet_to_inoutcut,gids)
+end
+
+function get_ls_to_bgfacet_to_inoutcut(cut::EmbeddedFacetDiscretization)
+  cut.ls_to_facet_to_inoutcut
+end
+
 function _consistent!(
   p_to_i_to_a::AbstractArray{<:Vector{<:Vector}},
   prange::PRange)
@@ -207,6 +262,28 @@ function change_bgmodel(
     cut.geo)
 end
 
+function change_bgmodel(
+  cut::EmbeddedFacetDiscretization,
+  newmodel::DiscreteModel,
+  facet_to_newfacet=1:num_facets(get_background_model(cut)))
+
+  nfacets = num_facets(newmodel)
+
+  ls_to_bgf_to_ioc = map(cut.ls_to_facet_to_inoutcut) do bgf_to_ioc
+    new_bgf_to_ioc = Vector{Int8}(undef,nfacets)
+    new_bgf_to_ioc[facet_to_newfacet] = bgf_to_ioc
+    new_bgf_to_ioc
+  end
+  subfacets = change_bgmodel(cut.subfacets,facet_to_newfacet)
+  EmbeddedFacetDiscretization(
+    newmodel,
+    ls_to_bgf_to_ioc,
+    subfacets,
+    cut.ls_to_subfacet_to_inout,
+    cut.oid_to_ls,
+    cut.geo)
+end
+
 function change_bgmodel(cells::SubCellData,cell_to_newcell)
   cell_to_bgcell = lazy_map(Reindex(cell_to_newcell),cells.cell_to_bgcell)
   SubCellData(
@@ -225,3 +302,20 @@ function change_bgmodel(facets::SubFacetData,cell_to_newcell)
     facets.point_to_coords,
     facets.point_to_rcoords)
 end
+
+function remove_ghost_subfacets(cut::EmbeddedFacetDiscretization,facet_gids)
+  bgfacet_mask = map(!iszero,local_to_owner(facet_gids))
+  subfacet_mask = map(Reindex(bgfacet_mask),cut.subfacets.cell_to_bgcell)
+  new_subfacets = findall(subfacet_mask)
+  subfacets = SubCellData(cut.subfacets,new_subfacets)
+  ls_to_subfacet_to_inout = map(cut.ls_to_subfacet_to_inout) do sf_to_io
+    map(Reindex(sf_to_io),new_subfacets)
+  end
+  EmbeddedFacetDiscretization(
+    cut.bgmodel,
+    cut.ls_to_facet_to_inoutcut,
+    subfacets,
+    ls_to_subfacet_to_inout,
+    cut.oid_to_ls,
+    cut.geo)
+end

src/Distributed/DistributedQuadratures.jl

@@ -0,0 +1,18 @@
+
+function CellData.Measure(
+  t::DistributedTriangulation,
+  quad::Tuple{MomentFitted,Vararg};
+  kwargs...)
+  @notimplemented
+  name, _args, _kwargs = quad
+  cut,cutfacets,_args... = _args
+  t = remove_ghost_cells(t)
+  measures = map(
+    local_views(t),
+    local_views(cut),
+    local_views(cutfacets)) do trian,cut,cutfacets
+      quad = name, (cut,cutfacets,_args...), _kwargs
+      Measure(trian,quad;kwargs...)
+  end
+  DistributedMeasure(measures)
+end

src/MomentFittedQuadratures/MomentFittedQuadratures.jl

@@ -43,15 +43,36 @@ function Quadrature(trian::Grid,
   Quadrature(trian,momentfitted,cut,geo,degree;in_or_out=in_or_out)
 end
 
+function Quadrature(trian::Grid,
+                    ::MomentFitted,
+                    cut::AbstractEmbeddedDiscretization,
+                    cutf::AbstractEmbeddedDiscretization,
+                    degree::Int;
+                    in_or_out=IN)
+  geo = get_geometry(cut)
+  Quadrature(trian,momentfitted,cut,cutf,geo,degree;in_or_out=in_or_out)
+end
+
+function Quadrature(trian::Grid,
+                    ::MomentFitted,
+                    cut::AbstractEmbeddedDiscretization,
+                    geo::CSG.Geometry,
+                    degree::Int;
+                    in_or_out=IN)
+  cutf = cut_facets(cut,geo)
+  Quadrature(trian,momentfitted,cut,cutf,geo,degree;in_or_out=in_or_out)
+end
+
 function Quadrature(active_mesh::Grid,
                     ::MomentFitted,
                     cut::AbstractEmbeddedDiscretization,
+                    cutf::AbstractEmbeddedDiscretization,
                     geo::CSG.Geometry,
                     degree::Int;
                     in_or_out=IN)
 
   acell_to_point_vals, acell_to_weight_vals = #
-    compute_lag_moments_from_leg(cut,geo,in_or_out,degree)
+    compute_lag_moments_from_leg(cut,cutf,geo,in_or_out,degree)
   acell_to_weight_vals = collect(get_array(acell_to_weight_vals))
 
   D = num_dims(active_mesh)
@@ -131,14 +152,15 @@ function legendreToMonomial(n::Int,d::Int)
 end
 
 function compute_lag_moments_from_leg(cut::AbstractEmbeddedDiscretization,
+                                      cutf::AbstractEmbeddedDiscretization,
                                       geo::CSG.Geometry,
                                       in_or_out,
                                       degree::Int)
   cut_trian = Triangulation(cut,CUT,geo)
   T = eltype(eltype(get_node_coordinates(cut_trian)))
   D = num_dims(cut_trian)
   b = MonomialBasis{D}(T,degree)
-  mon_contribs = compute_monomial_domain_contribution(cut,geo,in_or_out,b,degree)
+  mon_contribs = compute_monomial_domain_contribution(cut,cutf,geo,in_or_out,b,degree)
   mon_moments = compute_monomial_cut_cell_moments(cut_trian,mon_contribs,b)
   mon_to_leg = Fill(legendreToMonomial(degree,D),num_cells(cut_trian))
   leg_moments = lazy_map(*,mon_to_leg,mon_moments)
@@ -149,22 +171,23 @@ function compute_lag_moments_from_leg(cut::AbstractEmbeddedDiscretization,
 end
 
 function compute_monomial_domain_contribution(cut::AbstractEmbeddedDiscretization,
+                                              cutf::AbstractEmbeddedDiscretization,
                                               in_or_out,
                                               b::MonomialBasis,
                                               deg::Int)
   geo = get_geometry(cut)
-  compute_monomial_domain_contribution(cut,geo,in_or_out,b,deg)
+  compute_monomial_domain_contribution(cut,cutf,geo,in_or_out,b,deg)
 end
 
 function compute_monomial_domain_contribution(cut::AbstractEmbeddedDiscretization,
+                                              cutf::AbstractEmbeddedDiscretization,
                                               geo::CSG.Geometry,
                                               in_or_out::Integer,
                                               b::MonomialBasis,
                                               deg::Int)
 
   cut_io = CutInOrOut(in_or_out)
   dir_Γᵉ = (-1)^(in_or_out==OUT)
-  cutf = cut_facets(cut,geo)
   # Embedded facets
   Γᵉ = EmbeddedBoundary(cut,geo)
   # Interior fitted cut facets
@@ -354,4 +377,3 @@ function _get_terms_degrees(c::CartesianIndex)
 end
 
 end #module
-