Abnormal particle contact time in DEM simulations with DMP parallelization

test.zip (37.4 KB)

Hello everyone,

I am currently working on a CFD–DEM coupled model for particle sintering at high temperature. In my sintering force model, the key variable is the contact (sintering) time between particle pairs, so it is essential to correctly track the particle–particle contact duration.

My implementation follows the way the van der Waals (VDW) force is handled in calc_force_dem.Specifically, I modified the neighbor-pair data structure PFT_NEIGHBOR(:, CC) used in the DEM neighbor search and added an extra component to store the particle-pair contact time:

PFT_NEIGHBOR(4, CC) = T_CONTACT

Here, T_CONTACT represents the accumulated contact time for each particle pair.
For debugging and visualization, I also store the contact time in DES_USR_VAR_1.

After extensive testing, I observed the following:

  1. Serial runs and DMP builds executed on a single core behave correctly:

    • The contact time (tcontact) is accumulated as expected. (see Fig. 1)

    • The resulting cohesive (sintering) force is physically reasonable.

  2. Multi-core DMP parallel runs show abnormal behavior:

    • At DMP subdomain boundaries (tested with DMP = 3 × 3 × 1), the contact time tcontact is unexpectedly reset.

    • This leads to discontinuities in the cohesive force.

    • In Fig. 2, DES_USR_VAR_1 clearly shows the abnormal reset of particle contact time across DMP boundaries.

To isolate the issue, I have tried the following:

  • Modified mpi_pack_des_mod and mpi_unpack_des_mod to include the new contact-time variable, but this did not resolve the problem.

  • Tested a simplified case where tcontact is replaced by a prescribed function (i.e., not accumulated from particle-pair contact history).
    In this case, DMP parallel runs work correctly, with no abnormal behavior.

Based on these results, I believe the issue is specifically related to the data exchange of particle-pair contact time (tcontact) across DMP subdomains, rather than the force model itself.

This issue has been troubling me for quite a long time. I would greatly appreciate any guidance from the community or the developers on how this problem should be handled properly, or whether there exists a simpler and more robust approach to record and update the particle-pair contact time in DEM simulations, especially under DMP parallelization.

I have attached my test case below.
The abnormal behavior can be observed within 2–3 minutes of runtime.
All my code modifications are clearly marked with **** yq **** comments for easy identification.

Thank you very much for your time and help.
Finally, I wish everyone a Merry Christmas :christmas_tree: