diff --git a/NEWS.md b/NEWS.md index 9d92b19..1061242 100644 --- a/NEWS.md +++ b/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 diff --git a/src/DivAndCurlConformingFESpaces.jl b/src/DivAndCurlConformingFESpaces.jl index 4779474..2217c07 100644 --- a/src/DivAndCurlConformingFESpaces.jl +++ b/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 + diff --git a/test/HcurlProjectionTests.jl b/test/HcurlProjectionTests.jl new file mode 100644 index 0000000..671be7c --- /dev/null +++ b/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 diff --git a/test/mpi/runtests_np4.jl b/test/mpi/runtests_np4.jl index 3a15c54..d294615 100644 --- a/test/mpi/runtests_np4.jl +++ b/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") diff --git a/test/mpi/runtests_np4_body.jl b/test/mpi/runtests_np4_body.jl index 161a294..cc3753e 100644 --- a/test/mpi/runtests_np4_body.jl +++ b/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")