Skip to content

Commit

Permalink
* Some bug-fixes for 2D(adaptive)+1D(uniform) refinement
Browse files Browse the repository at this point in the history 
* Added missing methods for uniformly refining/coarsening in the regular
  quad/octree
  • Loading branch information
@amartinhuertas
amartinhuertas committed Mar 5, 2024
1 parent e62494d commit b6321f1
Show file tree
Hide file tree
Showing 4 changed files with 96 additions and 37 deletions.
66 changes: 40 additions & 26 deletions src/AnisotropicallyAdapted3DDistributedDiscreteModels.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -37,7 +37,7 @@ function AnisotropicallyAdapted3DDistributedDiscreteModel(
ptr_pXest_connectivity,
ptr_pXest,
pXest_type,
nothing,
PXestHorizontalRefinementRuleType(),
true,
nothing)
end
Expand Down Expand Up @@ -407,11 +407,14 @@ function generate_face_labeling(pXest_type::P6estType,

polytope = HEX



update_face_to_entity_with_ghost_data!(vertex_to_entity,
cell_prange,
num_faces(polytope,0),
cell_to_faces(topology,Dc,0))


update_face_to_entity_with_ghost_data!(edget_to_entity,
cell_prange,
num_faces(polytope,1),
Expand All @@ -421,11 +424,13 @@ function generate_face_labeling(pXest_type::P6estType,
cell_prange,
num_faces(polytope,Dc-1),
cell_to_faces(topology,Dc,Dc-1))


update_face_to_entity_with_ghost_data!(cell_to_entity,
cell_prange,
num_faces(polytope,Dc),
cell_to_faces(topology,Dc,Dc))


faces_to_entity=[vertex_to_entity,edget_to_entity,facet_to_entity,cell_to_entity]

Expand Down Expand Up @@ -476,6 +481,21 @@ face_labeling =

face_labeling
end

# map(partition(cell_prange)) do indices
# println("[$(MPI.Comm_rank(MPI.COMM_WORLD))] l2g=$(local_to_global(indices))")
# println("[$(MPI.Comm_rank(MPI.COMM_WORLD))] l2o=$(local_to_own(indices))")
# println("[$(MPI.Comm_rank(MPI.COMM_WORLD))] l2p=$(local_to_owner(indices))")


# cache=PartitionedArrays.assembly_cache(indices)

# println("[$(MPI.Comm_rank(MPI.COMM_WORLD))] ns=$(cache.neighbors_snd[])")
# println("[$(MPI.Comm_rank(MPI.COMM_WORLD))] nr=$(cache.neighbors_rcv[])")
# println("[$(MPI.Comm_rank(MPI.COMM_WORLD))] li_snd=$(cache.local_indices_snd[])")
# println("[$(MPI.Comm_rank(MPI.COMM_WORLD))] li_rcv=$(cache.local_indices_rcv[])")
# end

face_labeling
end

Expand Down Expand Up @@ -507,9 +527,7 @@ function _horizontally_refine_coarsen_balance!(model::OctreeDistributedDiscreteM
map(refinement_and_coarsening_flags,num_cols) do flags, num_cols
# The length of the local flags array has to match the number of locally owned columns in the model
@assert num_cols==length(flags)
println("BBB: $(flags)")
pXest_reset_data!(pXest_type, model.ptr_pXest, Cint(sizeof(Cint)), init_fn_callback_c, pointer(flags))
println("CCC: $(flags)")
end


Expand Down Expand Up @@ -559,16 +577,23 @@ function horizontally_adapt(model::OctreeDistributedDiscreteModel{Dc,Dp},
pXest_lnodes_destroy(model.pXest_type,ptr_pXest_lnodes)
pXest_refinement_rule_type = PXestHorizontalRefinementRuleType()

_extruded_flags = _extrude_refinement_and_coarsening_flags(_refinement_and_coarsening_flags,
model.ptr_pXest,
ptr_new_pXest)
extruded_ref_coarsen_flags=map(partition(get_cell_gids(model)),refinement_and_coarsening_flags) do indices, flags
similar(flags, length(local_to_global(indices)))
end

_extrude_refinement_and_coarsening_flags!(extruded_ref_coarsen_flags,
refinement_and_coarsening_flags,
model.ptr_pXest,
ptr_new_pXest)

stride = pXest_stride_among_children(model.pXest_type,pXest_refinement_rule_type,model.ptr_pXest)

adaptivity_glue = _compute_fine_to_coarse_model_glue(model.pXest_type,
pXest_refinement_rule_type,
model.parts,
model.dmodel,
fmodel,
_extruded_flags,
extruded_ref_coarsen_flags,
stride)
adaptive_models = map(local_views(model),
local_views(fmodel),
Expand All @@ -590,10 +615,11 @@ function horizontally_adapt(model::OctreeDistributedDiscreteModel{Dc,Dp},
return ref_model, adaptivity_glue
end

function _extrude_refinement_and_coarsening_flags(
refinement_and_coarsening_flags::MPIArray{<:Vector{<:Integer}},
ptr_pXest_old,
ptr_pXest_new)
function _extrude_refinement_and_coarsening_flags!(
extruded_flags::MPIArray{<:Vector{<:Integer}},
flags::MPIArray{<:Vector{<:Integer}},
ptr_pXest_old,
ptr_pXest_new)

pXest_old = ptr_pXest_old[]
pXest_new = ptr_pXest_new[]
Expand All @@ -602,19 +628,9 @@ end
num_trees = Cint(pXest_old.columns[].connectivity[].num_trees)
@assert num_trees == Cint(pXest_new.columns[].connectivity[].num_trees)

map(refinement_and_coarsening_flags) do flags
map(flags,extruded_flags) do flags,extruded_flags
current_old_quad=1
current_cell_old=1
total=0
for itree=0:num_trees-1
tree = pXest_tree_array_index(pXest_type,pXest_old,itree)[]
num_quads = Cint(tree.quadrants.elem_count)
q = pXest_quadrant_array_index(pXest_type,tree,0)
f,l=P6EST_COLUMN_GET_RANGE(q[])
num_layers=l-f
total+=num_quads*num_layers
end
extruded_flags=similar(flags, total)
current_cell_old=1

# Go over trees
for itree=0:num_trees-1
Expand Down Expand Up @@ -654,7 +670,5 @@ end
end
end
end

extruded_flags
end
end
end
56 changes: 48 additions & 8 deletions src/PXestTypeMethods.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2409,12 +2409,52 @@ end
function pXest_stride_among_children(pXest_type::P6estType,
::PXestHorizontalRefinementRuleType,
ptr_pXest)
# Here we are assuming:
# (1) Each processor has at least one column.
# (2) The number of layers in each column is the same within and accross processors.
pXest = ptr_pXest[]
tree = pXest_tree_array_index(pXest_type, pXest, 0)[]
q = pXest_quadrant_array_index(pXest_type, tree, 0)
f,l=P6EST_COLUMN_GET_RANGE(q[])
return l-f
# Here we are assuming:
# (1) Each processor has at least one column.
# (2) The number of layers in each column is the same within and accross processors.
num_trees = ptr_pXest[].columns[].connectivity[].num_trees
for itree = 0:num_trees-1
tree = pXest_tree_array_index(pXest_type, ptr_pXest[], itree)[]
if tree.quadrants.elem_count>0
q = pXest_quadrant_array_index(pXest_type, tree, 0)
f,l=P6EST_COLUMN_GET_RANGE(q[])
return l-f
end
end
return 0
end

function pXest_uniformly_refine!(::P4estType, ptr_pXest)
# Refine callbacks
function refine_fn_2d(::Ptr{p4est_t},which_tree::p4est_topidx_t,quadrant::Ptr{p4est_quadrant_t})
return Cint(1)
end
refine_fn_c = @cfunction($refine_fn_2d,Cint,(Ptr{p4est_t}, p4est_topidx_t, Ptr{p4est_quadrant_t}))
p4est_refine(ptr_pXest, Cint(0), refine_fn_c, C_NULL)
end

function pXest_uniformly_refine!(::P8estType, ptr_pXest)
# Refine callbacks
function refine_fn_3d(::Ptr{p8est_t},which_tree::p4est_topidx_t,quadrant::Ptr{p8est_quadrant_t})
return Cint(1)
end
# C-callable refine callback 3D
refine_fn_c = @cfunction($refine_fn_3d,Cint,(Ptr{p8est_t}, p4est_topidx_t, Ptr{p8est_quadrant_t}))
p8est_refine(ptr_pXest, Cint(0), refine_fn_c, C_NULL)
end

function pXest_uniformly_coarsen!(::P4estType, ptr_pXest)
function coarsen_fn_2d(::Ptr{p4est_t},::p4est_topidx_t,::Ptr{Ptr{p4est_quadrant_t}})
return Cint(1)
end
coarsen_fn_c=@cfunction($coarsen_fn_2d,Cint,(Ptr{p4est_t}, p4est_topidx_t, Ptr{Ptr{p4est_quadrant_t}}))
p4est_coarsen(ptr_pXest, Cint(0), coarsen_fn_c, C_NULL)
end

function pXest_uniformly_coarsen!(::P8estType, ptr_pXest)
function coarsen_fn_3d(::Ptr{p8est_t},::p4est_topidx_t,::Ptr{Ptr{p8est_quadrant_t}})
return Cint(1)
end
coarsen_fn_c=@cfunction($coarsen_fn_3d,Cint,(Ptr{p8est_t}, p4est_topidx_t, Ptr{Ptr{p8est_quadrant_t}}))
p8est_coarsen(ptr_pXest, Cint(0), coarsen_fn_c, C_NULL)
end
8 changes: 5 additions & 3 deletions test/PoissonAnisotropicOctreeModelsTests.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -35,6 +35,8 @@ module PoissonAnisotropicOctreeModelsTests
UH=TrialFESpace(VH,u)
num_local_cols=GridapP4est.num_locally_owned_columns(dmodel)
ref_coarse_flags=map(ranks,num_local_cols) do rank,num_local_cols
println("[rank:$(rank)] $(num_local_cols)")

flags=zeros(Cint,num_local_cols)
flags.=nothing_flag

Expand Down Expand Up @@ -237,7 +239,7 @@ module PoissonAnisotropicOctreeModelsTests

function test_3d(ranks,order,T::Type;num_amr_steps=5)
coarse_model=CartesianDiscreteModel((0,1,0,1),(1,1))
dmodel=AnisotropicallyAdapted3DDistributedDiscreteModel(ranks, coarse_model, 1, 1)
dmodel=AnisotropicallyAdapted3DDistributedDiscreteModel(ranks, coarse_model, 2, 2)
#writevtk(dmodel.dmodel,"model")
#test_refine_and_coarsen_at_once(ranks,dmodel,order,T)
rdmodel=dmodel
Expand All @@ -261,8 +263,8 @@ module PoissonAnisotropicOctreeModelsTests
end
end
function run(distribute)
debug_logger = ConsoleLogger(stderr, Logging.Debug)
global_logger(debug_logger); # Enable the debug logger globally
# debug_logger = ConsoleLogger(stderr, Logging.Debug)
# global_logger(debug_logger); # Enable the debug logger globally
ranks = distribute(LinearIndices((MPI.Comm_size(MPI.COMM_WORLD),)))
for perm=1:4, order=1:4, scalar_or_vector in (:scalar,)
test(ranks,perm,order,_field_type(Val{3}(),scalar_or_vector))
Expand Down
3 changes: 3 additions & 0 deletions test/runtests.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -47,6 +47,9 @@ function run_tests(testdir)
elseif f in ["PoissonNonConformingOctreeModelsTests.jl"]
np = [1,2,4]
extra_args = ""
elseif f in ["PoissonAnisotropicOctreeModelsTests.jl"]
np = [1,4]
extra_args = ""
else
np = [nprocs]
extra_args = ""
Expand Down

0 comments on commit b6321f1

Please sign in to comment.