Particle Dual Memory Spaces

[jeremylt]

This PR adds dual memory spaces so that we can start implementing GPU porting for the particle code.

This dual memory space setup is standard for GPU enabled computation, but there are a few things I want to call out in my implementation.

1. Split read-only and writable access to underlying data - In order to reduce the number of memory transfers between the host and device (very expensive), we need to know every time writable access is requested, as this is the only time a synchronization between host and device is needed.

2. Lazy initialization of device memory - We will not need to carry all values to the device at first, so I set up lazy allocation of device memory. The default is host-only memory, but then a device side mirror of this memory is allocated the first time device side memory is requested. This will be the exact same space if Kokkos is only build with CPU support and not GPU support, but we let Kokkos manage that logic via its DualView object.

3. Only touching ParticleContainer - It looked to me as though the ParticleStates object was deeply connected with MPI communication and the ParticleContainer holds the contiguous arrays of particle data that are used in the computationally heavy portions of the code. To that end, I only ported the ParticleContainer data to the dual memory space layout, as MPI communication is set up on the host. This can be revisited at a later time to enable GPU aware MPI, but that would be much more invasive of a change.

4. Minimal dev facing changes - Due to the split of pointer accessors to readable and writable, with and additional optional argument for memory space, developers working in a file not yet converted to GPU computation do not have to think about any new loop constructs or anything. Even once the file is converted, all that should change is writing the loop heads via Kokkos syntax and adding an extra argument to the pointer accessors to request device memory.

5. Internal mutability pattern - I am doing something a bit strange that I want to call out and need the opinion of someone with more C++ experience. The read-only accessors are const, but the Kokkos functions to sync the DualViews and the helper function I have to initialize the device side memory the first time it is requested are non-const functions. I would argue that this sort of internal mutation is not violating constness in any way the caller cares about - the underlying data is not changing, only the location (host or device). I have cast away the constness of these DualViews inside of the read-only pointer acessors and marked the relevant private object members as mutable (see here, here, and here). I modeled this based upon how I would approach this in C, but I don't know if there is a more idiomatic/correct way to do this in C++.
   Edit: This has been updated to use "pointer to impl" instead - it has the same effects on const-correctness but also allows us to keep the Kokkos headers out of the hpp file and keep them only in the cpp files as they get ported to Kokkos.

Note: Does compile successfully with #2012 (expected as these files still all are compiled for the CPU)

---

LLM Disclosure: Used Copilot plugin in VSCode to prompt for review the commits a few times. Because of the LLM output, I added commit 93baac7 which moved the checks earlier before the count is incremented. With some prodding (because I have discomfort with it), Copilot commented on the internal mutability pattern I'm using for the Views, and I'm not entirely satisfied with it myself.

I ended up replacing it with PImp as that also allowed me to move the Kokkos headers out of the .hpp file into the .cpp file, so that dependency isn't added to files that do not strictly need it.

---

Original description for posterity:

Very WIP at the moment, just making it easier to see the pieces come together

Managing data access is the first phase of the GPU porting. We need to know when write access is granted vs read-only to ensure the dual memory spaces (host and device) stay in sync. Kokkos will manage the sync for us, but we need to flag whenever we modify memory in either space, and I want to set things up so most developers do not have to track that detail, as its a potential failure point.

Once this is in place, then we can start porting computationally heavy loops over to Kokkos loops, which then can run on the device.

• Split direct access to particle state data read-only and writable
• Use access helpers internally for particle container object
• Switch states_ in ParticleContainer from array of arrays to arrays of Kokkos views

[jeremylt]

Ok, fixed some bad errors, these two commits should be clean now (of course I'll find ways they aren't going forward, I'm sure)

[jeremylt]

ack, I have misunderstood how much raw manipulation of the ParticleStates happens in the code instead of going through the owning ParticleContainer. Might need to make ParticleStates into a separate object since that's where the actual memory is held and passed around.

[jeremylt]

I haven't tried the latest commit yet, but I think it might be the best way to go, all things considered.Will need to think, write up, and discuss (and actually test and debug to make sure it works)

[jeremylt]

Ok, got the pieces of the design I want in place. Detailed rewrite of the PR description and some questions coming tomorrow.

[jeremylt]

@mayrmt perhaps also of interest to you for review? Not sure who all would be most interested

[Copilot]

🟡 Not ready to approve

The new “host” storage is currently a default-memory-space Kokkos view (can become device memory on GPU builds), which would break host-side pointer access and MPI assumptions.

Once you've addressed the issues Copilot identified, you can request another Copilot review.

This review doesn't count toward merge requirements. Sign up for the private preview to control whether Copilot approvals count.

Pull request overview

This PR introduces dual host/device memory handling for particle state storage to enable incremental GPU porting, primarily by upgrading ParticleContainer to manage host/device synchronization and by splitting read-only vs writable data access.

Changes:

• Added ParticleSpace (Host/Device) to distinguish memory spaces when requesting particle state pointers.
• Reworked ParticleContainer state storage to use Kokkos views/dual views with lazy device initialization and explicit sync/modify behavior.
• Updated particle algorithms/interactions to use the new writable state accessors where mutation occurs.

File summaries

File	Description
src/particle/tests/4C_particle_container_test.cpp	Updates tests to use new pointer access APIs (needs adjustment to keep read-only API coverage).
src/particle/src/rigidbody/4C_particle_rigidbody.cpp	Switches state writes to get_ptr_to_state_writable / cond_get_ptr_to_state_writable.
src/particle/src/interaction/4C_particle_interaction_sph_temperature.cpp	Uses writable conditional access for temperature derivative updates.
src/particle/src/interaction/4C_particle_interaction_sph_surface_tension.cpp	Switches multiple state updates to writable accessors.
src/particle/src/interaction/4C_particle_interaction_sph_surface_tension_recoilpressure_evaporation.cpp	Uses writable accessor for acceleration updates.
src/particle/src/interaction/4C_particle_interaction_sph_surface_tension_interface_viscosity.cpp	Uses writable accessor for acceleration updates.
src/particle/src/interaction/4C_particle_interaction_sph_surface_tension_barrier_force.cpp	Uses writable accessor for acceleration updates.
src/particle/src/interaction/4C_particle_interaction_sph_rigid_particle_contact.cpp	Uses writable conditional accessor for force updates.
src/particle/src/interaction/4C_particle_interaction_sph_pressure.cpp	Uses writable accessor for pressure assignment.
src/particle/src/interaction/4C_particle_interaction_sph_peridynamic.cpp	Uses writable accessors for bond/force/acceleration/damage updates.
src/particle/src/interaction/4C_particle_interaction_sph_open_boundary.cpp	Uses writable accessors for boundary state updates.
src/particle/src/interaction/4C_particle_interaction_sph_momentum.cpp	Uses writable accessors for acceleration/modified acceleration/force updates.
src/particle/src/interaction/4C_particle_interaction_sph_heatsource.cpp	Uses writable accessor for temperature derivative updates.
src/particle/src/interaction/4C_particle_interaction_sph_heatloss_evaporation.cpp	Uses writable accessor for temperature derivative updates.
src/particle/src/interaction/4C_particle_interaction_sph_density.cpp	Uses writable accessors for density/colorfield accumulation and updates.
src/particle/src/interaction/4C_particle_interaction_sph_boundary_particle.cpp	Uses writable accessors for boundary particle state initialization.
src/particle/src/interaction/4C_particle_interaction_dem.cpp	Uses writable accessors for initialization and acceleration computation.
src/particle/src/interaction/4C_particle_interaction_dem_contact.cpp	Uses writable accessors for force/moment updates in contacts.
src/particle/src/interaction/4C_particle_interaction_dem_adhesion.cpp	Uses writable accessors for force updates in adhesion.
src/particle/src/engine/4C_particle_engine.cpp	Uses writable accessors where engine mutates particle states.
src/particle/src/engine/4C_particle_engine_typedefs.hpp	Adds a typedef alias for the new particle memory space enum.
src/particle/src/engine/4C_particle_engine_enums.hpp	Introduces ParticleSpace and declaration for enum-to-name conversion (has a Doxygen grouping issue).
src/particle/src/engine/4C_particle_engine_enums.cpp	Implements enum_to_space_name.
src/particle/src/engine/4C_particle_engine_container.hpp	Core API change: read-only vs writable accessors + host/device selection + DualView management (needs a HostSpace fix and const-cast cleanup).
src/particle/src/engine/4C_particle_engine_container.cpp	Implements Kokkos-backed allocation/resizing and updated add/remove particle logic (has a signed/unsigned warning risk).
src/particle/src/algorithm/4C_particle_algorithm_timint.cpp	Uses writable accessor for position perturbation.
src/particle/src/algorithm/4C_particle_algorithm_temperature_bc.cpp	Uses writable accessor for temperature updates.
src/particle/src/algorithm/4C_particle_algorithm_initial_field.cpp	Uses writable accessor for initial field assignment.
src/particle/src/algorithm/4C_particle_algorithm_dirichlet_bc.cpp	Uses writable accessors for applying Dirichlet BCs.
src/particle/src/algorithm/4C_particle_algorithm_constraints.cpp	Simplifies optional-state handling via conditional writable accessors.

Review details

Comments suppressed due to low confidence (3)

src/particle/src/engine/4C_particle_engine_container.hpp:216

• Same as host sync: the C-style cast to call sync_device() is unnecessary with the current mutable members; removing it avoids misleading const-casting and simplifies the internal-mutability approach.

        // casting away const-ness because this only changes the representation of the data
        (*(Kokkos::DualView<double*>*)&states_dual_[state]).sync_device();
        return (const double*)&(states_dual_[state].view_device().data()[index * statedim_[state]]);

src/particle/tests/4C_particle_container_test.cpp:345

• Same issue as above: this should use the read-only conditional accessor and a const pointer to match the test intent.

      double* currvel = container_->cond_get_ptr_to_state_writable(Particle::Velocity, index);

src/particle/tests/4C_particle_container_test.cpp:348

• Same issue as above: keep this as a read-only conditional accessor so the test continues to validate the const API.

      double* currmass = container_->cond_get_ptr_to_state_writable(Particle::Mass, index);

• Files reviewed: 30/30 changed files
• Comments generated: 7
• Review effort level: Low

We're testing this review assessment. Please use 👍 or 👎 to tell us if it's correct.

[jeremylt]

We do not have the ability to actually test device usage yet (see #2012) and I'm hesitant to add a test implying that we do

[Copilot]

🟡 Not ready to approve

A few call sites request writable state access even though they only read, which undermines the intended host/device sync minimization and should be corrected before approval.

Once you've addressed the issues Copilot identified, you can request another Copilot review.

This review doesn't count toward merge requirements. Sign up for the private preview to control whether Copilot approvals count.

Review details

Comments suppressed due to low confidence (7)

src/particle/tests/4C_particle_container_test.cpp:349

• These assertions only read velocity/mass, but they currently request writable pointers. Using the read-only accessor better matches intent and avoids marking the DualView as modified (which can cause unnecessary syncs/transfers).

      double* currvel = container_->cond_get_ptr_to_state_writable(Particle::Velocity, index);
      FOUR_C_EXPECT_ITERABLE_NEAR(currvel, vel.begin(), 3, 1.0e-14);

      double* currmass = container_->cond_get_ptr_to_state_writable(Particle::Mass, index);
      EXPECT_NEAR(currmass[0], mass[0], 1e-14);

src/particle/src/engine/4C_particle_engine_container.hpp:260

• Same issue as the read-only accessor: any space value other than Host falls through to the device branch. Making this an explicit else if (space == Device) and throwing otherwise avoids silently doing the wrong thing on invalid input.

      else
      {
        if (!is_states_dual_valid_[state]) init_state_dual(state);
        states_dual_[state].sync_device();
        states_dual_[state].modify_device();

src/particle/src/engine/4C_particle_engine_container.cpp:125

• This assertion compares ParticleState (an enum) against states.size() (size_t). On common warning levels this triggers a signed/unsigned comparison warning. Casting to a common signed type here avoids noisy builds.

    if (state < states.size() and not states[state].empty() and
        static_cast<int>(states[state].size()) != statedim_[state])
      FOUR_C_THROW("can not add particle: dimensions of state '{}' do not match!",
          enum_to_state_name(state));

src/particle/src/engine/4C_particle_engine.cpp:1283

• convert_pos_to_gid takes a const double* and this code only reads the position, but it currently requests a writable pointer. That will mark the state as modified and can trigger unnecessary host/device synchronizations (and expensive transfers once GPU is enabled).

      double* currpos = container->get_ptr_to_state_writable(Position, ownedindex);

src/particle/src/interaction/4C_particle_interaction_dem.cpp:554

• This function only reads angular velocity to compute kinetic energy, but it currently requests writable access. That needlessly marks the state as modified and can force extra sync/copies between memory spaces.

    double* angvel = container->cond_get_ptr_to_state_writable(Particle::AngularVelocity, 0);

src/particle/src/engine/4C_particle_engine_container.hpp:216

• get_ptr_to_state treats any space value other than Host as the device path. If an invalid/uninitialized ParticleSpace value is ever passed, this will silently take the wrong branch and may allocate/sync device memory unexpectedly. Consider making this an explicit else if (space == Device) and throwing for unknown values.

This issue also appears on line 256 of the same file.

      else
      {
        if (!is_states_dual_valid_[state]) init_state_dual(state);
        states_dual_[state].sync_device();
        return &(states_dual_[state].view_device().data()[index * statedim_[state]]);

src/particle/src/engine/4C_particle_engine_container.hpp:541

• Doc comment has a stray * at the end of the \brief line, which breaks the formatting.

     * \brief initialize DualView for on-device computation     *

• Files reviewed: 30/30 changed files
• Comments generated: 0 new
• Review effort level: Low

We're testing this review assessment. Please use 👍 or 👎 to tell us if it's correct.

[georghammerl]

I like the nice side effect of this PR to clearly distinguish between read and write access.

[ppraegla]

Thanks @jeremylt for the detailed PR description and your effort in porting the particle code to GPUs. This will be very interesting. I haven't worked with Kokkos so far. So, it will take me some time to read into Kokkos and review the PR.

Could you open an issue describing the plan and, if possible, concrete steps for porting the particle code to GPUs? I think this would be nice to keep track of the progress. Also, this would give me (and others) a better idea of what is planned.

[ppraegla]

I also think now is a good time to introduce a performance test. With @jeremylt and @vovannikov both working on performance improvements, it would be nice to see this in a dedicated test.

I will think about a DEM and SPH performance tests.

[jeremylt]

I have been putting my plan over in #2070, though looking at it now, it could be much clearer.

[ppraegla]

I have been putting my plan over in #2070, though looking at it now, it could be much clearer.

I was only looking at the issues and didn't see that there was already a discussion. I will have a look.

[slfuchs]

@jeremylt Welcome to 4C or at least welcome to the particle framework in 4C 😄 I am currently seldom on Github, but since implementing this particle framework kind of from scratch during my PhD time, I am still interested in the particle framework.

I am curious how you proceed with porting the framework for GPU computations. As @georghammerl mentioned, I also think that to distinguish between read and write access to the particle container improves the implementation a lot. And also as @ppraegla said, I think we should definitely run performance tests here and with all other upcoming PRs. We should not lose performance for pure CPU computations. @vovannikov is doing some performance changes currently, so probably he can also provide some tests.

[mayrmt]

@jeremylt I assumed hat this PR relies on Kokkos with a serial backend from Trilinos. Do you have concrete plans for other Kokkos backends, e.g. such as in #2012?

[jeremylt]

@mayrmt this PR doesn't rely upon any specific backend. You can write code with Kokkos that only works for a specific backend, but that's generally the wrong thing to do - the Kokkos team discourages that as it breaks portability and makes the code more complicated.

This PR sets the 'device' memory space to be whatever is the preferred ('default' in their parlance) execution space for Kokkos, based upon however Kokkos was compiled. So once #2012 merges the device memory space will automatically be CUDA when Kokkos is compiled with CUDA support and the work here will set the dual memory spaces as CPU and CUDA. If OpenMP is used instead, then the device memory space is a CPU allocation and Kokkos will ensure that both memory spaces will simply point to a single CPU allocation (also true with Kokkos only built for serial CPU execution).

[jeremylt]

@slfuchs the basic plan is to

1. identify a candidate simulation (simpler better initially, but needs to support ideally millions of local particles)
2. identify the computationally heavy loops (likely in interaction files for initial cases)
3. rewrite those loops as Kokkos loops
4. profile and adjust the loop bodies to improve performance

We then repeat steps 1-4, working our way through the different files and simulations.

One thing that will help us identify good candidate loops or loops that may prove difficult is the write patterns in the loop bodies. The body of each loop will become a single GPU or OpenMP thread, and we will get the best performance if a single thread only writes to a single memory location. Multiple threads writing to the same location is sometimes unavoidable, but its generally something we want to minimize as the resulting atomic memory operations are a known bottleneck.

I think if I'm understanding the code correctly, some of the interactions are saved to interaction pair objects that later write the results back out to the particles? That's good for minimizing these atomic operations, but I'll have to set up these pair objects with similar dual memory spaces too.