Making ParticlesHandler work with dealii::Triangulation and mesh refinement/coarsening

[lium...@gmail.com]

Hi deal.II Community,

I am new to the library and I'm trying to implement a prototype Stokes flow solver with tracer particles using shared memory computing. The overall procedure is something like:

1) Solve velocities and pressure analogous to Step-22.

2) Execute coarsening and refinement

3) Interpolate velocities to particles and move particles.

What I found is after grid coarsening, the particles which originally belong to the finer grid loses its host cell and an exception is thrown from ParticleAccessor::get_surrounding_cell().

I read through step-19 on Github. From my understanding the implementation there does not involve AMR and therefore such problems do not occur. I also consulted step-68 and 70. Those two implementations use parallel::distributed::Triangulation and attach pre- and post-refinement callback functions to save and reload particles. However in these functions the dynamic_cast from Triangulation to p::d::Triangulation is invalid and returns a nullptr.

My question is am I correct to assume that the library does not yet fully support particles on dealii::Triangulation with AMR? If yes can I ask for some guidance on how to implement it myself? One naive way out seems to be to store and reload procedure as is done in the library for p::d::Triangulation, but is it too costly?

I really appreciate your time and help.

Cheers,

Qihang.

[Wolfgang Bangerth]

On 12/15/20 2:49 PM, lium...@gmail.com wrote:
>
> 1) Solve velocities and pressure analogous to Step-22.
> 2) Execute coarsening and refinement
> 3) Interpolate velocities to particles and move particles.
>
> What I found is after grid coarsening, the particles which originally belong
> to the finer grid loses its host cell and an exception is thrown from
> ParticleAccessor::get_surrounding_cell().

Qihang,
Interesting. Do you think you could put together a short testcase that
illustrates this issue?

Best
W.

--
------------------------------------------------------------------------
Wolfgang Bangerth email: bang...@colostate.edu
www: http://www.math.colostate.edu/~bangerth/

[blais...@gmail.com]

Dear Wu,

Have you tried without the dynamic cast? Can you bind the callback of a regular triangulation?

I must admit I have never tried to use particles with regular triangulation (and I think that outside of step-19, not much people have done it ).

[WU, Qihang]

Hi Wolfgang, Bruno and all community members,

My apologies for a delayed reply as I was pondering on a good way to replicate the problem without bogging you down with the unnecessary dependencies in my program. So what I come up with is to adapt the step-22 by adding a particle handler object and mimic what is done in step-19. I was able to reproduce the error for refinement. I am working on a way to reproduce the error for coarsening. The overall idea is:

- spawn one particle in one cell

- refine that cell

- An exception will be thrown in sort_particles_into_domains_and_cells():

--------------------------------------------------------

An error occurred in line <495> of file </home/wuqihang/.local/src/dealii-9.2.0/include/deal.II/grid/tria_iterator.templates.h> in function

    dealii::TriaActiveIterator<Accessor>& dealii::TriaActiveIterator<Accessor>::operator=(const dealii::TriaIterator<Accessor>&) [with Accessor = dealii::CellAccessor<2, 2>]

The violated condition was: 

    this->accessor.has_children() == false

Additional information: 

    (none)

Stacktrace:

-----------

#0  /home/wuqihang/.local/lib/libdeal_II.g.so.9.2.0: dealii::TriaActiveIterator<dealii::CellAccessor<2, 2> >::operator=(dealii::TriaIterator<dealii::CellAccessor<2, 2> > const&)

#1  /home/wuqihang/.local/lib/libdeal_II.g.so.9.2.0: dealii::Particles::ParticleHandler<2, 2>::update_cached_numbers()

#2  /home/wuqihang/.local/lib/libdeal_II.g.so.9.2.0: dealii::Particles::ParticleHandler<2, 2>::sort_particles_into_subdomains_and_cells()

#3  /home/wuqihang/Dropbox/Coding/DEAL.II/Step22ParticlesAMR/build/step-22_PtclAMR: Step22::StokesParticles<2, double>::move_particles()

#4  /home/wuqihang/Dropbox/Coding/DEAL.II/Step22ParticlesAMR/build/step-22_PtclAMR: Step22::StokesProblem<2>::second_step()

#5  /home/wuqihang/Dropbox/Coding/DEAL.II/Step22ParticlesAMR/build/step-22_PtclAMR: Step22::StokesProblem<2>::run()

#6  /home/wuqihang/Dropbox/Coding/DEAL.II/Step22ParticlesAMR/build/step-22_PtclAMR: main

--------------------------------------------------------

My limited understanding is that after the grid refinement the particle is not correctly associated with the finer child cell. I have not setup the callback function in dealii::Triangulation for this test because the dynamic_cast in the callback functions will throw an error anyways. What I will try now is to follow Bruno's advice to remove dynamic_cast in the callback functions and see what happens. I will report back if I shall make any progress.

I very much appreciate your time bearing with me and I look forward to helping out to fix the problem.

Thank you, 

Qihang.