Reporting a series of errors with void parts + tentative solutions
amartinhuertas opened this issue · 1 comments
amartinhuertas commented
Hi @fverdugo,
while exploring in a parallel context several scenarios in which the portion of a given processor is void, I have found 3 different errors, which I report below as MWEs along with possible solutions. (Note there is also an attached file needed to reproduce the first error).
Can you take a look at these, and let me know what do you think?
Once we agree on the approach to solution, I will send a PR with the corresponding fixes.
Thanks!
@amartinhuertas
using Gridap
using Gridap.Geometry
using Gridap.ReferenceFEs
using Gridap.Arrays
using Gridap.Io
using Gridap.FESpaces
using Gridap.CellData
u(x)=x[1]
# Fixes Err1 (see below)
function Gridap.Io.from_dict(::Type{UnstructuredGrid},dict::Dict{Symbol,Any})
x = collect1d(dict[:node_coordinates])
T = eltype(x)
Dp = dict[:Dp]
if (length(x)!=0)
# Had to remove type annotation from this line (see corresponding line in Gridap)!!!
# There seems to be a BUG in the compiler (Julia v1.7.2)
node_coordinates = reinterpret(Point{Dp,T},x)
else
T1=Point{Dp,Float64}
node_coordinates=Vector{T1}(undef,0)
println(eltype(node_coordinates))
end
cell_node_ids = from_dict(Table{Int32,Vector{Int32},Vector{Int32}},dict[:cell_node_ids])
reffes = [from_dict(LagrangianRefFE,reffe) for reffe in dict[:reffes]]
cell_type::Vector{Int8} = dict[:cell_type]
O::Bool = dict[:orientation]
UnstructuredGrid(
node_coordinates,
cell_node_ids,
reffes,
cell_type,
O ? Oriented() : NonOriented())
end
# Fixes Err2 (see below)
function Base.map(::typeof(testitem),
a::Tuple{<:AbstractVector{<:AbstractVector{<:VectorValue}},<:AbstractVector{<:Gridap.Fields.LinearCombinationFieldVector}})
a2=testitem(a[2])
a1=Vector{eltype(eltype(a[1]))}(undef,size(a2,1))
a1.=zero(testitem(a1))
(a1,a2)
end
# Fixes Err3 (see below)
function Gridap.Geometry.is_change_possible(
strian::Gridap.Geometry.Triangulation,
ttrian::Gridap.Geometry.Triangulation)
if strian === ttrian || num_cells(strian)==num_cells(ttrian)==0
return true
end
Gridap.Helpers.@check get_background_model(strian) === get_background_model(ttrian) "Triangulations do not point to the same background discrete model!"
D = num_cell_dims(strian)
sglue = get_glue(strian,Val(D))
tglue = get_glue(ttrian,Val(D))
is_change_possible(sglue,tglue) # Fails here
end
# Fixes Err3 (see below)
function Gridap.CellData.change_domain(a::CellField,
::ReferenceDomain,
ttrian::Triangulation,
::ReferenceDomain)
msg = """\n
We cannot move the given CellField to the reference domain of the requested triangulation.
Make sure that the given triangulation is either the same as the triangulation on which the
CellField is defined, or that the latter triangulation is the background of the former.
"""
strian = get_triangulation(a)
if strian === ttrian || num_cells(strian)==num_cells(ttrian)==0
return a
end
@assert Gridap.Geometry.is_change_possible(strian,ttrian) msg
D = num_cell_dims(strian)
sglue = get_glue(strian,Val(D))
tglue = get_glue(ttrian,Val(D))
change_domain_ref_ref(a,ttrian,sglue,tglue)
end
model=Gridap.Io.from_json_file(DiscreteModel,"model_part_2.json") # Fails here (Err1)
# FE Spaces
order=1
reffe = ReferenceFE(lagrangian,Float64,order)
V = TestFESpace(model,reffe,dirichlet_tags="boundary")
U = TrialFESpace(u,V) # Fails here (Err2)
s = get_fe_dof_basis(V)
trian = get_triangulation(s)
f = CellField(u,trian,DomainStyle(s))
b = change_domain(f,s.domain_style)
s(f) # Fails here (Err2 at a lower level)
sd=get_data(s)
bd=get_data(b)
sdi=testitem(sd)
sbd=testitem(bd)
evaluate(sdi,sbd) # Fails here (Err2 at a lower level)
uh = interpolate(u,U)
Ω = view(get_triangulation(U),Int[])
dΩ = Measure(Ω,2)
dc=∫(uh)dΩ # Fails here (Err3)
trian_f=get_triangulation(uh)
trian_x=dΩ.quad.trian
Gridap.CellData.is_change_possible(trian_f,trian_x) # Fails here (Err3)
JordiManyer commented
I think this has been resolved.