GROMACS version:2025.02
GROMACS modification: No
Hello everyone! I am simulating 2 Al2O3 surfaces with water between these surfaces. After 10 ns NPT simulation, I wanted to calculate the water density distribution between two surfaces by gmx density. However, there are some issues.
Command line:
gmx density -f npt2.xtc -o density-h2o-npt2 -b 5000 -sl 1500 -s npt2.tpr -ng 2
Reading file npt2.tpr, VERSION 2025.1 (single precision)
Select 2 groups to calculate density for:
Group 0 ( System) has 31422 elements
Group 1 ( Other) has 7800 elements
Group 2 ( Al2O3) has 7800 elements
Group 3 ( Water) has 23622 elements
Group 4 ( SOL) has 23622 elements
Group 5 ( non-Water) has 7800 elements
Select a group: 2
Selected 2: 'Al2O3'
Select a group: 3
Selected 3: 'Water'
Reading frame 0 time 5000.000
-------------------------------------------------------
Program: gmx density, version 2025.1
Source file: src/gromacs/pbcutil/mshift.cpp (line 829)
Fatal error:
There are inconsistent shifts over periodic boundaries in a molecule type
consisting of 31422 atoms. The longest distance involved in such interactions
is 8.752 nm which is above half the box length. This molecule type consists of
multiple parts, e.g. monomers, that are connected by interactions that are not
chemical bonds, e.g. restraints. Such systems can not be treated. The only
solution is increasing the box size.
I have tried “gmx trjconv -s npt2.tpr -f npt2.xtc -o npt2_mol.xtc -pbc mol” to solve the particles crossing boundary problem as most people have attempted to do.
Although I have read gonzaloqe’s post https://gromacs.bioexcel.eu/t/inconsistent-shifts-over-periodic-boundaries-in-a-molecule-type/4999, which is similar to my issue. In fact, I only considered the bond and angle (Al-O-H) in the Clayff force field, and did not use the so-called artificial planar dihedral as gonzaloqe said.
When I select a large step size in the nstxtcout option, I can use the gmx trjconv command to process it. Unfortunately, I cannot repeat this command when I reduce the nstxtcout.
Were the trajectories generated with the periodic molecule yes option. Without the gro file and xtc file it would be hard to diagnose the problem. If you are willing to share I can help.
Thank you, Arun. At first, I thought that the issue was caused by the use of the NPT ensemble, which led to the detachment of -OH groups on the Al₂O₃ surface, resulting in the problem mentioned above. However, after switching to the NVT ensemble and fixing all parts of Al₂O₃ except for the surface H atoms, the same problem still occurred. Below are the gro and xtc files generated from the NVT simulation at 298 K. It seems that I cannot upload gro and xtc files here.
Thanks for the file. I do not think NVT, NPT ensemble would matter here. it depends on how you have built the topology of the system (restraints, angles dihedrals and so on). Can you run the simulation of Al2O3 in vacuum (NVT) with no water, just alumina surfaces exposed to vacuum on both sides for a long time with no issues and use the post processing tools of gromacs. Are the potential energies stable ? is the Alumina surface periodic in x and y ? if so how are you handling it ? are the atoms in the alumina position restrained ? Are bonds defined for OH groups on the alumina in the topology file ?
Thank you again, Arun, for your prompt response. While examining the trajectory files, I noticed that the -OH groups on the Al₂O₃ surface, which has periodicity in the xy direction, were crossing the boundary. Regarding the treatment of the xy periodicity, I first created a supercell and then added a vacuum layer in the z direction before applying the periodic boundary conditions in the x, y, and z directions. During the simulation, I fixed all parts of Al₂O₃ except for the surface H atoms and defined the O-H groups on the surface according to the Clayff force field. I will run the simulation in a vacuum environment without water later.
From the error message I can only guess that the topology is not built correctly. can you please share the top file and the gro file for the Al2O3 alone (without water) when you have time.
Thanks
arun
I ran the case in a vacuum as well, and I still encountered the issue of H atoms crossing the boundary, which prevents the use of tools like gmx density. However, the potential energy has been ensured to converge.
Thank you. I am not sure what do you mean by H-atoms are crossing the boundaries. Which H atoms are you referring to?. H attached to O are covalently bonded. They can cross the periodic boundary in the x -y dimension and if you have vacuum on top and bottom and you are running an NVT simulation they cannot cross in the z dimension. I would check your topology and the way you are functionalizing your alumina surface
Thank you so much, arun. I observed that hydrogen atoms on the aluminum oxide surface may have crossed the boundary, causing displacement exceeding half of the box, thereby rendering gmx density unusable.
Hi wunwun. Let’s diagnose the simplest case with no liquid since the problem persists with no liquid, it has to do something with alumina surface. Can you also please provide gro file of al2O3 that corresponds to the itp file you sent. Thanks
When not using the -center option, I think molecules do not have to be made whole. We should remove this operation in that case. I think you can circumvent the issue by commenting out all lines containing gpbc in src/gromacs/gmxana/gmx_density.cpp and recompiling Gromacs.
