Added EmbeddedCollections

scripts/Embedded/MWEs/updateable_trians.jl

@@ -16,36 +16,6 @@ end
 
 φ(r) = interpolate(x->sqrt((x[1]-0.5)^2+(x[2]-0.5)^2)-r,V_φ) # Circle
 
-struct EmbeddedTriangulationCollection
-  func::Function
-  trians::Dict{Symbol,Triangulation}
-  bgmodel::DiscreteModel
-  cutgeo::EmbeddedDiscretization
-end
-
-function EmbeddedTriangulationCollection(
-  func::Function,bgmodel::DiscreteModel,φ0
-)
-  geo = DiscreteGeometry(φ0,bgmodel)
-  cutgeo = cut(bgmodel,geo)
-  trians = Dict(pairs(func(cutgeo)))
-  return EmbeddedTriangulationCollection(func,trians,bgmodel,cutgeo)
-end
-
-function update!(c::EmbeddedTriangulationCollection,φh)
-  geo = DiscreteGeometry(φh,c.bgmodel)
-  cutgeo = cut(c.bgmodel,geo)
-  trians = func(cutgeo)
-  for (key,value) in pairs(trians)
-    c.trians[key] = value
-  end
-  return c
-end
-
-function Base.getindex(c::EmbeddedTriangulationCollection,key)
-  return c.trians[key]
-end
-
 order = 1
 model = generate_model(2,10)
 
@@ -56,19 +26,22 @@ reffe = ReferenceFE(lagrangian,Float64,order)
 V_φ = TestFESpace(model,reffe)
 
 φ0 = φ(0.2)
-Ωs = EmbeddedTriangulationCollection(model,φ0) do cutgeo
+Ωs = EmbeddedCollection(model,φ0) do cutgeo
+  Ω = Triangulation(cutgeo,PHYSICAL_IN)
+  Γ = EmbeddedBoundary(cutgeo)
   (; 
-    :Ωin  => Triangulation(cutgeo,PHYSICAL_IN),
-    :Ωout => Triangulation(cutgeo,PHYSICAL_OUT),
-    :Γ    => EmbeddedBoundary(cutgeo)
+    :Ω  => Ω,
+    :Γ  => Γ,
+    :dΩ => Measure(Ω,2*order),
+    :dΓ => Measure(Γ,2*order)
   )
 end
 
-area(Ωs) = sum(∫(1.0)*Measure(Ωs[:Ωin],2))
-contour(Ωs) = sum(∫(1.0)*Measure(Ωs[:Γ],2))
+area(Ωs) = sum(∫(1.0)*Ωs[:dΩ])
+contour(Ωs) = sum(∫(1.0)*Ωs[:dΓ])
 
 for r in 0.2:0.1:0.5
-  update!(Ωs,φ(r))
+  update_collection!(Ωs,φ(r))
   println(" >> Radius: $r")
   println(" >> Area: $(area(Ωs)) (expected: $(π*r^2))")
   println(" >> Contour: $(contour(Ωs)) (expected: $(2π*r))")

src/Embedded/Embedded.jl

@@ -1,3 +1,4 @@
 
 include("DifferentiableTriangulations.jl")
 include("SubFacetBoundaryTriangulations.jl")
+include("EmbeddedCollections.jl")

src/Embedded/EmbeddedCollections.jl

@@ -0,0 +1,28 @@
+
+struct EmbeddedCollection
+  generate :: Function
+  objects  :: Dict{Symbol,Any}
+  bgmodel  :: Union{<:DiscreteModel,<:DistributedDiscreteModel}
+
+  function EmbeddedCollection(
+    generate::Function,bgmodel,φ0
+  )
+    geo = DiscreteGeometry(φ0,bgmodel)
+    cutgeo = cut(bgmodel,geo)
+    objects = Dict{Symbol,Any}(pairs(generate(cutgeo)))
+    new(generate,objects,bgmodel)
+  end
+end
+
+(c::EmbeddedCollection)(φh) = update_collection!(c,φh)
+
+function update_collection!(c::EmbeddedCollection,φh)
+  geo = DiscreteGeometry(φh,c.bgmodel)
+  cutgeo = cut(c.bgmodel,geo)
+  merge!(c.objects,pairs(c.generate(cutgeo)))
+  return c
+end
+
+function Base.getindex(c::EmbeddedCollection,key)
+  return c.objects[key]
+end

src/GridapTopOpt.jl

@@ -52,8 +52,8 @@ export evaluate_derivatives!
 
 include("Embedded/Embedded.jl")
 export DifferentiableTriangulation
-export SubFacetBoundaryTriangulation
-export SubFacetSkeletonTriangulation
+export SubFacetBoundaryTriangulation, SubFacetSkeletonTriangulation
+export EmbeddedCollection, update_collection!
 
 include("Utilities.jl")
 export SmoothErsatzMaterialInterpolation

test/seq/EmbeddedCollectionsTests.jl

@@ -0,0 +1,59 @@
+module EmbeddedCollectionsTests
+using Test
+
+using GridapTopOpt
+using Gridap
+
+using GridapEmbedded
+using GridapEmbedded.LevelSetCutters
+using Gridap.Geometry, Gridap.FESpaces, Gridap.CellData, Gridap.Adaptivity
+
+function generate_model(D,n)
+  domain = (D==2) ? (0,1,0,1) : (0,1,0,1,0,1)
+  cell_partition = (D==2) ? (n,n) : (n,n,n)
+  base_model = UnstructuredDiscreteModel((CartesianDiscreteModel(domain,cell_partition)))
+  ref_model = refine(base_model, refinement_method = "barycentric")
+  model = ref_model.model
+  return model
+end
+
+function main()
+  order = 1
+  model = generate_model(2,40)
+
+  reffe = ReferenceFE(lagrangian,Float64,order)
+  V_φ = TestFESpace(model,reffe)
+
+  φ(r) = interpolate(x->sqrt((x[1]-0.5)^2+(x[2]-0.5)^2)-r,V_φ) # Circle
+
+  φ0 = φ(0.2)
+  Ωs = EmbeddedCollection(model,φ0) do cutgeo
+    Ω = Triangulation(cutgeo,PHYSICAL_IN)
+    Γ = EmbeddedBoundary(cutgeo)
+    (; 
+      :Ω  => Ω,
+      :Γ  => Γ,
+      :dΩ => Measure(Ω,2*order),
+      :dΓ => Measure(Γ,2*order)
+    )
+  end
+
+  area(Ωs) = sum(∫(1.0)*Ωs[:dΩ])
+  contour(Ωs) = sum(∫(1.0)*Ωs[:dΓ])
+
+  for r in 0.2:0.1:0.5
+    update_collection!(Ωs,φ(r))
+    A = area(Ωs)
+    C = contour(Ωs)
+    println(" >> Radius: $r")
+    println(" >> Area: $(A) (expected: $(π*r^2))")
+    println(" >> Contour: $(C) (expected: $(2π*r))")
+    @test abs(A - π*r^2) < 1e-3
+    @test abs(C - 2π*r) < 1e-3
+    println("---------------------------------")
+  end
+end
+
+main()
+
+end

test/seq/EmbeddedTests.jl → test/seq/EmbeddedDifferentiationTests.jl

@@ -1,4 +1,4 @@
-module EmbeddedTests
+module EmbeddedDifferentiationTests
 using Test
 
 using GridapTopOpt

test/seq/runtests.jl

@@ -10,6 +10,9 @@ using Test
 @time @testset "Inverter - HPM" begin include("InverterHPMTests.jl") end
 @time @testset "PZMultiFieldRepeatingState - ALM" begin include("PZMultiFieldRepeatingStateTests.jl") end
 @time @testset "JLD2SaveLoad" begin include("JLD2SaveLoad.jl") end
-@time @testset "Embedded" begin include("EmbeddedTests.jl") end
+@time @testset "Embedded" begin 
+  include("EmbeddedDifferentiationTests.jl")
+  include("EmbeddedCollectionsTests.jl")
+end
 
 end # module