Skip to content
Draft
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
6 changes: 6 additions & 0 deletions NEWS.md
Original file line number Diff line number Diff line change
Expand Up @@ -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

### Added

- New overloads for the `TrialFESpace` constructor where the data to be imposed is passed as a `DistributedCellField`. Since PR[#183](https://github.com/gridap/GridapDistributed.jl/pull/183).

## [0.4.10] - 2025-09-29

### Added
Expand Down
25 changes: 25 additions & 0 deletions src/FESpaces.jl
Original file line number Diff line number Diff line change
Expand Up @@ -414,6 +414,21 @@ function FESpaces.TrialFESpace!(f::DistributedSingleFieldFESpace,fun)
DistributedSingleFieldFESpace(spaces,f.gids,f.trian,f.vector_type,f.metadata)
end

function FESpaces.TrialFESpace(f::DistributedSingleFieldFESpace,cf::DistributedCellField)
spaces = map(local_views(f),local_views(cf)) do s, field
TrialFESpace(s,field)
end
DistributedSingleFieldFESpace(spaces,f.gids,f.trian,f.vector_type,f.metadata)
end

function FESpaces.TrialFESpace(cf::DistributedCellField,f::DistributedSingleFieldFESpace)
spaces = map(local_views(f),local_views(cf)) do s, field
TrialFESpace(s,field)
end
DistributedSingleFieldFESpace(spaces,f.gids,f.trian,f.vector_type,f.metadata)
end


Comment on lines +424 to +431
Copy link

Copilot AI Nov 2, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

[nitpick] The two TrialFESpace overloads at lines 417-422 and 424-429 have identical implementations. The second overload (where DistributedCellField is the first argument) simply calls the same underlying Gridap TrialFESpace(s,field) function. Consider consolidating these into a single implementation or adding a comment explaining why both argument orders are necessary.

Suggested change
function FESpaces.TrialFESpace(cf::DistributedCellField,f::DistributedSingleFieldFESpace)
spaces = map(local_views(f),local_views(cf)) do s, field
TrialFESpace(s,field)
end
DistributedSingleFieldFESpace(spaces,f.gids,f.trian,f.vector_type,f.metadata)
end

Copilot uses AI. Check for mistakes.
function FESpaces.HomogeneousTrialFESpace(f::DistributedSingleFieldFESpace)
spaces = map(f.spaces) do s
HomogeneousTrialFESpace(s)
Expand Down Expand Up @@ -478,6 +493,16 @@ function FESpaces.interpolate_dirichlet!(
FEFunction(f,free_values,dirichlet_values)
end

function FESpaces.interpolate_dirichlet!(
u::DistributedCellField, free_values::AbstractVector,
dirichlet_values::AbstractArray{<:AbstractVector},
f::DistributedSingleFieldFESpace)
map(local_views(u), f.spaces,local_views(free_values),dirichlet_values) do u,V,fvec,dvec
interpolate_dirichlet!(u,fvec,dvec,V)
end
FEFunction(f,free_values,dirichlet_values)
end

function FESpaces.interpolate_everywhere(u, f::DistributedSingleFieldFESpace)
free_values = zero_free_values(f)
dirichlet_values = get_dirichlet_dof_values(f)
Expand Down
131 changes: 76 additions & 55 deletions test/FESpacesTests.jl
Original file line number Diff line number Diff line change
Expand Up @@ -85,66 +85,87 @@ function main(distribute,parts,das)
model = CartesianDiscreteModel(ranks,parts,domain,cells)
Ω = Triangulation(model)
Γ = Boundary(model)
dΩ = Measure(Ω,3)

function _interpolation_tests(u,model,Ω,Γ,output)
reffe = ReferenceFE(lagrangian,Float64,1)
V = TestFESpace(model,reffe,dirichlet_tags="boundary")
U = TrialFESpace(u,V)
V2 = FESpace(Ω,reffe)
@test get_vector_type(V) <: PVector{<:Vector}
@test get_vector_type(U) <: PVector{<:Vector}
@test get_vector_type(V2) <: PVector{<:Vector}

free_values_partition = map(partition(V.gids)) do indices
ones(Float64,local_length(indices))
end

free_values = PVector(free_values_partition,partition(V.gids))
fh = FEFunction(U,free_values)
zh = zero(U)
uh = interpolate(u,U)
eh = u - uh

uh_dir = interpolate_dirichlet(u,U)
free_values = zero_free_values(U)
dirichlet_values = get_dirichlet_dof_values(U)
uh_dir2 = interpolate_dirichlet!(u,free_values,dirichlet_values,U)

uh_everywhere = interpolate_everywhere(u,U)
dirichlet_values0 = zero_dirichlet_values(U)
uh_everywhere_ = interpolate_everywhere!(u,free_values,dirichlet_values0,U)
eh2 = u - uh_everywhere
eh2_ = u - uh_everywhere_

uh_everywhere2 = interpolate_everywhere(uh_everywhere,U)
uh_everywhere2_ = interpolate_everywhere!(uh_everywhere,free_values,dirichlet_values,U)
eh3 = u - uh_everywhere2

dofs = get_fe_dof_basis(U)
cell_vals = dofs(uh)
gather_free_values!(free_values,U,cell_vals)
gather_free_and_dirichlet_values!(free_values,dirichlet_values,U,cell_vals)
uh4 = FEFunction(U,free_values,dirichlet_values)
eh4 = u - uh4

dΩ = Measure(Ω,3)
Copy link

Copilot AI Nov 2, 2025

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The variable is redefined inside _interpolation_tests at line 131, shadowing the outer defined at line 88. While this works, it's confusing and the outer at line 88 appears to be unused. Consider removing the outer definition or using a different variable name in the inner function.

Copilot uses AI. Check for mistakes.
cont = ∫( abs2(eh) )dΩ
cont2 = ∫( abs2(eh2) )dΩ
cont2_ = ∫( abs2(eh2_) )dΩ
cont3 = ∫( abs2(eh3) )dΩ
cont4 = ∫( abs2(eh4) )dΩ
@test sqrt(sum(cont)) < 1.0e-9
@test sqrt(sum(cont2)) < 1.0e-9
@test sqrt(sum(cont2_)) < 1.0e-9
@test sqrt(sum(cont3)) < 1.0e-9
@test sqrt(sum(cont4)) < 1.0e-9

writevtk(Ω,joinpath(output,"Ω"), nsubcells=10,
celldata=["err"=>cont[Ω]],
cellfields=["uh"=>uh,"zh"=>zh,"eh"=>eh])

writevtk(Γ,joinpath(output,"Γ"),cellfields=["uh"=>uh])
end

# Testing interpolation with a Julia function
u((x,y)) = x+y
_interpolation_tests(u,model,Ω,Γ,output)

# Testing interpolation with a DistributedCellField
# The DistributedCellField MUST be defined in all cells
# (i.e., owned + ghost) for the tests within
# _interpolation_tests(...) to pass. The problematic part is
# in the interpolation of the dirichlet values
Ωghosts = Triangulation(with_ghost, model)
fields = map(local_views(Ωghosts)) do Ω
CellField(u,Ω)
end
u_cf = GridapDistributed.DistributedCellField(fields,Ωghosts)
_interpolation_tests(u_cf,model,Ωghosts,Γ,output)

reffe = ReferenceFE(lagrangian,Float64,1)
V = TestFESpace(model,reffe,dirichlet_tags="boundary")
U = TrialFESpace(u,V)
V2 = FESpace(Ω,reffe)
@test get_vector_type(V) <: PVector{<:Vector}
@test get_vector_type(U) <: PVector{<:Vector}
@test get_vector_type(V2) <: PVector{<:Vector}

free_values_partition = map(partition(V.gids)) do indices
ones(Float64,local_length(indices))
end

free_values = PVector(free_values_partition,partition(V.gids))
fh = FEFunction(U,free_values)
zh = zero(U)
uh = interpolate(u,U)
eh = u - uh

uh_dir = interpolate_dirichlet(u,U)
free_values = zero_free_values(U)
dirichlet_values = get_dirichlet_dof_values(U)
uh_dir2 = interpolate_dirichlet!(u,free_values,dirichlet_values,U)

uh_everywhere = interpolate_everywhere(u,U)
dirichlet_values0 = zero_dirichlet_values(U)
uh_everywhere_ = interpolate_everywhere!(u,free_values,dirichlet_values0,U)
eh2 = u - uh_everywhere
eh2_ = u - uh_everywhere_

uh_everywhere2 = interpolate_everywhere(uh_everywhere,U)
uh_everywhere2_ = interpolate_everywhere!(uh_everywhere,free_values,dirichlet_values,U)
eh3 = u - uh_everywhere2

dofs = get_fe_dof_basis(U)
cell_vals = dofs(uh)
gather_free_values!(free_values,U,cell_vals)
gather_free_and_dirichlet_values!(free_values,dirichlet_values,U,cell_vals)
uh4 = FEFunction(U,free_values,dirichlet_values)
eh4 = u - uh4

dΩ = Measure(Ω,3)
cont = ∫( abs2(eh) )dΩ
cont2 = ∫( abs2(eh2) )dΩ
cont2_ = ∫( abs2(eh2_) )dΩ
cont3 = ∫( abs2(eh3) )dΩ
cont4 = ∫( abs2(eh4) )dΩ
@test sqrt(sum(cont)) < 1.0e-9
@test sqrt(sum(cont2)) < 1.0e-9
@test sqrt(sum(cont2_)) < 1.0e-9
@test sqrt(sum(cont3)) < 1.0e-9
@test sqrt(sum(cont4)) < 1.0e-9


writevtk(Ω,joinpath(output,"Ω"), nsubcells=10,
celldata=["err"=>cont[Ω]],
cellfields=["uh"=>uh,"zh"=>zh,"eh"=>eh])

writevtk(Γ,joinpath(output,"Γ"),cellfields=["uh"=>uh])

# Assembly
Ωass = Triangulation(das,model)
Expand Down
Loading