Merge pull request #160 from gridap/bugfix-num-cells

Bugfix: Wrong num_cells for triangulations containing ghosts

NEWS.md

@@ -5,6 +5,12 @@ All notable changes to this project will be documented in this file.
 The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/),
 and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html).
 
+## [Unreleased]
+
+### Fixed
+
+- Fixed bug in `num_cells` in the case where a `DistributedTriangulation` contained ghost cells. Since PR[#160](https://github.com/gridap/GridapDistributed.jl/pull/160).
+
 ## [0.4.4] 2024-08-14
 
 ### Added

src/Geometry.jl

@@ -474,8 +474,21 @@ function Geometry.get_background_model(a::DistributedTriangulation)
   a.model
 end
 
-function Geometry.num_cells(a::DistributedTriangulation)
-  sum(map(trian->num_cells(trian),local_views(a)))
+function Geometry.num_cells(a::DistributedTriangulation{Df}) where Df
+  gids = get_face_gids(a.model,Df)
+  n_loc_ocells = map(local_views(a),partition(gids)) do a, gids
+    glue = get_glue(a,Val(Df))
+    @assert isa(glue,FaceToFaceGlue)
+    tcell_to_mcell = glue.tface_to_mface
+    if isa(tcell_to_mcell,IdentityVector)
+      own_length(gids)
+    else
+      mcell_to_owned = local_to_own(gids)
+      is_owned(mcell) = !iszero(mcell_to_owned[mcell])
+      sum(is_owned,tcell_to_mcell)
+    end
+  end
+  return sum(n_loc_ocells)
 end
 
 # Triangulation constructors
@@ -670,19 +683,23 @@ function add_ghost_cells(dtrian::DistributedTriangulation)
   add_ghost_cells(dmodel,dtrian)
 end
 
-function _covers_all_faces(dmodel::DistributedDiscreteModel{Dm},
-                           dtrian::DistributedTriangulation{Dt}) where {Dm,Dt}
-  covers_all_faces=map(local_views(dmodel),local_views(dtrian)) do model, trian
+function _covers_all_faces(
+  dmodel::DistributedDiscreteModel{Dm},
+  dtrian::DistributedTriangulation{Dt}
+) where {Dm,Dt}
+  covers_all_faces = map(local_views(dmodel),local_views(dtrian)) do model, trian
     glue = get_glue(trian,Val(Dt))
     @assert isa(glue,FaceToFaceGlue)
     isa(glue.tface_to_mface,IdentityVector)
   end
   reduce(&,covers_all_faces,init=true)
 end
 
-function add_ghost_cells(dmodel::DistributedDiscreteModel{Dm},
-                         dtrian::DistributedTriangulation{Dt}) where {Dm,Dt}
-  covers_all_faces=_covers_all_faces(dmodel,dtrian)
+function add_ghost_cells(
+  dmodel::DistributedDiscreteModel{Dm},
+  dtrian::DistributedTriangulation{Dt}
+) where {Dm,Dt}
+  covers_all_faces = _covers_all_faces(dmodel,dtrian)
   if (covers_all_faces)
     trians = map(local_views(dmodel)) do model
       Triangulation(ReferenceFE{Dt},model)
@@ -703,7 +720,7 @@ function add_ghost_cells(dmodel::DistributedDiscreteModel{Dm},
     cache=fetch_vector_ghost_values_cache(mcell_intrian,partition(gids))
     fetch_vector_ghost_values!(mcell_intrian,cache) |> wait
 
-    dreffes=map(local_views(dmodel)) do model
+    dreffes = map(local_views(dmodel)) do model
       ReferenceFE{Dt}
     end
     trians = map(Triangulation,dreffes,local_views(dmodel),mcell_intrian)

test/GeometryTests.jl

@@ -59,15 +59,19 @@ function main(distribute,parts)
 
   Ω = Triangulation(with_ghost,model)
   writevtk(Ω,joinpath(output,"Ω"))
+  @test num_cells(Ω) == num_cells(model)
 
   Ω = Triangulation(no_ghost,model)
   writevtk(Ω,joinpath(output,"Ω"))
+  @test num_cells(Ω) == num_cells(model)
 
   Γ = Boundary(with_ghost,model,tags="boundary")
   writevtk(Γ,joinpath(output,"Γ"))
+  nbfacets = num_cells(Γ)
 
   Γ = Boundary(no_ghost,model,tags="boundary")
   writevtk(Γ,joinpath(output,"Γ"))
+  @test num_cells(Γ) == nbfacets
 
   function is_in(coords)
     R = 1.6