Merge pull request #199 from gridap/add_hcurl_projection_tests

Add Hcurl-projection test + fix BUG

Author: amartinhuertas

NEWS.md

@@ -7,6 +7,12 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 
 ## [Unreleased]
 
+### Fixed
+- BUG fix in `_generate_sign_flips(...)` private helper function for Nedelec elements. Since PR[#199](https://github.com/gridap/GridapDistributed.jl/pull/199).
+
+### Added
+- An Hcurl projection test case. Since PR[#199](https://github.com/gridap/GridapDistributed.jl/pull/199).
+
 ## [0.4.12] - 2026-03-14
 
 ### Changed

src/DivAndCurlConformingFESpaces.jl

@@ -87,9 +87,8 @@ function _common_fe_space_constructor(model,cell_reffes;conformity,split_own_and
   DistributedSingleFieldFESpace(spaces,gids,trian,vector_type)
 end
 
-function _generate_sign_flips(model,cell_reffes)
-  cell_gids  = get_cell_gids(model)
-  sign_flips = map(local_views(model),partition(cell_gids),cell_reffes) do m, p, cell_reffe
+function _generate_sign_flips(model, cell_gids, cell_reffes)
+sign_flips = map(local_views(model),partition(cell_gids),cell_reffes) do m, p, cell_reffe
     D = num_cell_dims(model)
 
     gtopo = get_grid_topology(m)
@@ -145,4 +144,34 @@ function _generate_sign_flips(model,cell_reffes)
   cache = fetch_vector_ghost_values_cache(sign_flips,partition(cell_gids))
   fetch_vector_ghost_values!(sign_flips,cache) |> wait
   sign_flips
+end 
+
+
+function _generate_sign_flips(model::DistributedDiscreteModel,
+                             cell_reffes::AbstractVector{<:AbstractVector{<:GenericRefFE{<:RaviartThomas}}})
+  cell_gids  = get_cell_gids(model)
+  _generate_sign_flips(model,cell_gids,cell_reffes)
 end
+
+function _generate_sign_flips(model::DistributedDiscreteModel,
+                              cell_reffes::AbstractVector{<:AbstractVector{<:GenericRefFE{<:Nedelec}}})
+  Dc=num_cell_dims(model)
+  is_simp = false 
+  map(cell_reffes) do cell_reffes
+    cell_reffe = cell_reffes[1]
+    polytope = get_polytope(cell_reffe)
+    if (Gridap.Geometry.is_simplex(polytope))
+        is_simp=true
+    end 
+  end 
+  if (Dc==2 || is_simp)
+    map(local_views(model), cell_reffes) do model, cell_reffes
+        Gridap.FESpaces._no_sign_flip(model, cell_reffes)
+    end
+  else
+    @assert Dc==3
+    cell_gids  = get_cell_gids(model)
+    _generate_sign_flips(model,cell_gids,cell_reffes)
+  end
+end
+

test/HcurlProjectionTests.jl

@@ -0,0 +1,97 @@
+module HcurlProjectionTests
+
+using Gridap
+using PartitionedArrays
+using GridapDistributed
+using Test
+
+u_ex_2D((x,y)) = 2*VectorValue(-y,x)
+f_ex_2D(x) = u_ex_2D(x)
+u_ex_3D((x,y,z)) = 2*VectorValue(-y,x,0.) - VectorValue(0.,-z,y)
+f_ex_3D(x) = u_ex_3D(x)
+
+function get_analytical_functions(Dc)
+  if Dc==2
+    return u_ex_2D, f_ex_2D
+  else
+    @assert Dc==3
+    return u_ex_3D, f_ex_3D
+  end
+end
+
+function solve_hcurl_projection(model::GridapDistributed.DistributedDiscreteModel{Dc},order) where {Dc}
+    u_ex, f_ex=get_analytical_functions(Dc)
+
+    V = FESpace(model,
+                ReferenceFE(nedelec,order),
+                conformity=:Hcurl,
+                dirichlet_tags="boundary")
+    
+    U = TrialFESpace(V,u_ex)
+    
+    trian = Triangulation(model)
+    degree = 2*(order+1)
+    dΩ = Measure(trian,degree)
+        
+    a(u,v) = ∫( (∇×u)⋅(∇×v) + u⋅v )dΩ
+    l(v) = ∫(f_ex⋅v)dΩ
+
+    op = AffineFEOperator(a,l,U,V)
+    if (num_free_dofs(U)==0)
+      # UMFPACK cannot handle empty linear systems
+      uh = zero(U)
+    else
+      uh = solve(op)
+    end
+    uh,U
+  end 
+
+  function check_error_hcurl_projection(model::GridapDistributed.DistributedDiscreteModel{Dc},order,uh) where {Dc}
+    trian = Triangulation(model)
+    degree = 2*(order+1)
+    dΩ = Measure(trian,degree)
+
+    u_ex, f_ex = get_analytical_functions(Dc)
+    
+    eu = u_ex - uh
+
+    l2(v) = sqrt(sum(∫(v⋅v)*dΩ))
+    hcurl(v) = sqrt(sum(∫(v⋅v + (∇×v)⋅(∇×v))*dΩ))
+    
+    eu_l2 = l2(eu)
+    eu_hcurl = hcurl(eu)
+    
+    tol = 1.0e-6
+    @test eu_l2 < tol
+    @test eu_hcurl < tol
+end
+
+ function test_2d(ranks,parts,order)
+    domain = (0,1,0,1)
+    model = CartesianDiscreteModel(ranks,parts,domain,(4,4))
+    solve_hcurl_projection(model,order) |> x -> check_error_hcurl_projection(model,order,x[1])
+ end 
+
+ function test_3d(ranks,parts,order)
+    domain = (0,1,0,1,0,1)
+    model = CartesianDiscreteModel(ranks,parts,domain,(4,4,4))
+    solve_hcurl_projection(model,order) |> x -> check_error_hcurl_projection(model,order,x[1])
+ end 
+
+function main(distribute,parts)
+  ranks = distribute(LinearIndices((prod(parts),)))
+  
+  if length(parts)==2 
+      for order=0:2
+        test_2d(ranks,parts,order)
+      end
+  elseif length(parts)==3
+      for order=0:2
+        test_3d(ranks,parts,order)
+      end
+  else 
+      @assert false 
+  end 
+end
+  
+end #module

test/mpi/runtests_np4.jl

@@ -1,5 +1,5 @@
 module NP4
-# All test running on 4 procs here
+# All tests running on either 1 or 4 procs here
 
 using GridapDistributed
 using PartitionedArrays
@@ -35,6 +35,7 @@ include("../VisualizationTests.jl")
 include("../AutodiffTests.jl")
 include("../ConstantFESpacesTests.jl")
 include("../MacroDiscreteModelsTests.jl")
+include("../HcurlProjectionTests.jl")
 
 include("runtests_np4_body.jl")
 

test/mpi/runtests_np4_body.jl

@@ -48,6 +48,10 @@ function all_tests(distribute, parts)
 
     StokesHdivDGTests.main(distribute, parts)
     PArrays.toc!(t, "StokesHdivDG")
+
+    HcurlProjectionTests.main(distribute, parts)
+    PArrays.toc!(t, "HcurlProjection")
+    
   end
 
   if TESTCASE ∈ ("all", "mpi-transient")