Sundials particle integrator

[srichers]

It seems the Sundials integrator in Src/Extern/SUNDIALS/AMReX_SundialsIntegrator.H is specific to multifabs. For those who implemented sundials, would it be straightforward to extend the support to ParticleContainers as well?

[gardner48]

In principle, yes. I think the only thing that would be needed is an N_Vector wrapper like the one for a MultiFab here to enable SUNDIALS to operate on a ParticleContainer.

If the ODE for moving the particles was x' = v, then you could store the velocities returned by the right-hand side function in the particle positions of the output ParticleContainer. In this case the N_Vector wrapper should work for any ParticleContainer. Instead if the ODE was x' = v, v' = a, then the particles would also need to hold the velocity in their rdata. In this case, the N_Vector wrapper would need access the rdata making it less generic but still shouldn't be hard to wrap.

[srichers]

My particles contain many fields describing the internal state of the particle (around 20-30 fields), including the velocity, and we do hold a copy of the position and velocity in the rdata. We currently use the methods in AMReX_RKIntegrator.H for time integration, and I am working on implementing a simple adaptive timestepping scheme there, but it would be very useful to have access to things like multirate schemes.

[srichers]

See https://github.com/AMReX-Astro/Emu/blob/release/Source/FlavoredNeutrinoContainer.H - there is a line where we include a python-generated file, as we generate code for different scenarios requiring different numbers of fields.

[gardner48]

Thanks for the pointer to the code. From ApplyFlavoredNeutrinoRHS it looks like you need the integrator to update all the rdata fields. That should make things a little easier, the N_Vector wrapper would just need to treat the rdata the same way it treats the data underlying a MultiFab in the current wrapper. Two extra considerations I can think of for wrapping ParticleContainers are:

1. The wrappers will need to iterate over multiple containers simultaneously (e.g., the input and output vectors to N_VAbs_MultiFab. Maybe someone that's worked ParticleContainers more than I have can comment on if that works just like it does for MultiFabs.

2. When particles are redistributed, the internal SUNDIALS memory will need to be updated accordingly. Adding a method to call the Resize function for the integrator should be all that's needed and I can easily add that to the interface.

[atmyers]

If the two ParticleContainers have been defined with the same BoxArrays and DistributionMappings, you can index them with the same ParIter, as with MultiFab.