grompp crashes in computeEffectiveAtomDensity
Summary
When an atom is outside the box, index calculation in computeEffectiveAtomDensity
produces negative values, and GROMACS crashes due to heap corruption.
Specifically, the problem seems to be caused by negative coordinates (-8.229) of magnitude larger than the box size (6.2), so that our + numCells[d]
is not enough to make the value positive.
Unit test
For src/gromacs/mdlib/tests/calc_verletbuf.cpp
:
// Issue #5002
TEST(EffectiveAtomDensity, LargeValuesHandledWell)
{
const std::vector<RVec> coordinates = { { 13.132, -8.229, -2.700 } };
const matrix box = { { 6.2, 0, 0 }, { 0, 6.2, 0 }, { 0, 0, 6.2 } };
const real cutoff = 1;
const real referenceDensity = (1) / (coordinates.size() * gmx::power3<real>(6.2/6));
const real density = computeEffectiveAtomDensity(coordinates, box, cutoff, MPI_COMM_NULL);
EXPECT_FLOAT_EQ(density, referenceDensity);
}
Exact steps to reproduce
Reproducible with current 2024 and 2023 versions:
$ ~/gromacs-2023/build/clang-17/bin/gmx grompp -f run.mdp -c system.gro -p system.top -o run.tpr
:-) GROMACS - gmx grompp, 2023.5-dev-20240219-db90c39977 (-:
Executable: /home/aland/gromacs-2023/build/clang-17/bin/gmx
Data prefix: /home/aland/gromacs-2023 (source tree)
Working dir: /home/aland/issue-f8399
Command line:
gmx grompp -f run.mdp -c system.gro -p system.top -o run.tpr
NOTE 1 [file run.mdp]:
For a correct single-point energy evaluation with nsteps = 0, use
continuation = yes to avoid constraining the input coordinates.
WARNING 1 [file run.mdp]:
The Berendsen thermostat does not generate the correct kinetic energy
distribution, and should not be used for new production simulations (in
our opinion). We would recommend the V-rescale thermostat.
Setting the LD random seed to -1120502286
Generated 1 of the 1 non-bonded parameter combinations
Excluding 3 bonded neighbours molecule type 'molecule'
Analysing residue names:
There are: 1 Other residues
Analysing residues not classified as Protein/DNA/RNA/Water and splitting into groups...
Number of degrees of freedom in T-Coupling group System is 0.00
Determining Verlet buffer for a tolerance of 0.005 kJ/mol/ps at 300 K
munmap_chunk(): invalid pointer
Aborted (core dumped)
Original report: https://gromacs.bioexcel.eu/t/gmx-grompp-crashes-with-a-very-simple-system/8399/1
Edited by Andrey Alekseenko