Understanding minimum domain size with MPI

[arnau-jr]

I'm trying to understand the limits of domain decomposition in HOOMD. We are running a system with 2 MPI CPU ranks and so the domain decomposition that we obtain is:
Using domain decomposition: n_x = 1 n_y = 1 n_z = 2.
When running the NPT simulation the system contracts until it reaches the minimum threshold we obtain the following message:

RuntimeError: Communication error -
Simulation box too small for domain decomposition.
r_ghost_max: 3.45523
d.z/2: 3.45489

Units in nm. Now I understand why we are getting this message but I'm not sure why we get that specific size.

According to the HOOMD documentation, the minimum domain size is twice the largest cutoff + cell buffer. In our case the largest cutoff is 0.994 nm and we use the default cell buffer of 0.4, which means that the minimum size should be 0.994+0.4=1.394 and not 3.45523 as we are getting.

My questions then are the following:

1. Am I misunderstanding the meaning of r_ghost_max? If it is half the minimum domain size then why I am getting 3.45523 and not 1.394, which is the value that we should get according to the documentation?
2. Is there any way to reduce the minimum domain size besides reducing the cutoff? We tried changing to the Stencil neighbor list with a width of 0.8 and a buffer of 0.4 but we got the same error with the same r_ghost_max.
   I realize that reducing the minimum domain size below the cutoff + buffer would wreck performance but I want to know if this possibility is contemplated in HOOMD.

[joaander]

Do you use constrain.Rigid?

• r_ghost_max = max(pairwise_rcuts + buffer) + max(rigid_body_diameters).
• d.z/2 is half the domain length.

The pairwise r_cut values and the rigid body diameters depend on the parameters you set in the simulation (pair.r_cut and rigid.body), not on the choice of neighbor list data structure. The requirement that r_ghost_max < min(d)/2 is a fundamental limitation of HOOMD-blue's choice of domain decomposition.

[arnau-jr]

You are absolutely right, we do have a two kinds of rigids bodies, the largest one is a nanoparticle with diameter 2.06123 which results in 0.994+0.4+2.06123=3.45523. In hindsight, this makes a lot of sense, since you need the whole rigid body to be in the same domain for the computations to work properly.

Do you think that this should be added on to the documentation? I'm pretty sure the rigid body diameter being a factor on minimum domain size is not mentioned in the domain decomposition page. It just states that

Given an interaction range r_interaction (for MD, this is the sum of the largest pair potential r_cut and the neighbor list buffer)

Thank you!

[joaander]

Yes, this should be mentioned in the documentation (and it currently is not).