# Error in GROMACS for adding ions to a solution

I am following Justin Lemkul's GROMACS simulation given on this website here for GROMACS.

I believe I have followed the process to the T. I would like to mention that in the "Generate Topology" section of the tutorial, he says to include

#ifdef POSRES_B
#include "posre_Protein_chain_B.itp"
#endif


However, this was already present in my posre_Protein_chain_B.itp file. However this statement was also present in my posre_Protein_chain_A.itp, posre_Protein_chain_C.itp, posre_Protein_chain_D.itp, posre_Protein_chain_E.itp files. For example, posre_Protein_chain_A.itp had

#ifdef POSRES_A
#include "posre_Protein_chain_A.itp"
#endif


in it. I have commented them out, so it now looks like:

; #ifdef POSRES_A
; #include "posre_Protein_chain_A.itp"
; #endif


Furthermore, every time I was asked to choose for N-terminus, I hit option 2: None, and everytime I was asked to choose a C-terminus, I hit option 0: COO-. This resulted in the system having a net negative charge I believe.

I don't think I have deviated anywhere else in this simulation. However, I am getting the error message in the attachment when I run

gmx grompp -f ions.mdp -c solv.gro -p topol.top -o ions.tpr

After running this command I received the following error message:

GROMACS is free software; you can redistribute it and/or modify it
under the terms of the GNU Lesser General Public License

GROMACS:      gmx grompp, version 2020.5
Executable:   /home/satyend/.local/bin/gmx
Data prefix:  /home/satyend/.local
Working dir:  /home/satyend/software/gromacs_gpu/umbrella-sampling-protofibril
Command line:
gmx grompp -f ions.mdp -c solv.gro -p topol.top -o ions.tpr

Ignoring obsolete mdp entry 'ns_type'

NOTE 1 [file ions.mdp]:
With Verlet lists the optimal nstlist is >= 10, with GPUs >= 20. Note
that with the Verlet scheme, nstlist has no effect on the accuracy of

Setting the LD random seed to 2139084382
Generated 165 of the 1596 non-bonded parameter combinations
Excluding 3 bonded neighbours molecule type 'Protein_chain_A'
Excluding 3 bonded neighbours molecule type 'Protein_chain_B'
Excluding 3 bonded neighbours molecule type 'Protein_chain_C'
Excluding 3 bonded neighbours molecule type 'Protein_chain_D'
Excluding 3 bonded neighbours molecule type 'Protein_chain_E'
Excluding 2 bonded neighbours molecule type 'SOL'

WARNING 1 [file topol.top, line 54]:
The GROMOS force fields have been parametrized with a physically
incorrect multiple-time-stepping scheme for a twin-range cut-off. When
used with a single-range cut-off (or a correct Trotter
multiple-time-stepping scheme), physical properties, such as the density,
might differ from the intended values. Since there are researchers
actively working on validating GROMOS with modern integrators we have not
yet removed the GROMOS force fields, but you should be aware of these
issues and check if molecules in your system are affected before
proceeding. Further information is available at
https://redmine.gromacs.org/issues/2884 , and a longer explanation of our
decision to remove physically incorrect algorithms can be found at
https://doi.org/10.26434/chemrxiv.11474583.v1 .

NOTE 2 [file topol.top, line 54]:
System has non-zero total charge: -10.000000
Total charge should normally be an integer. See
http://www.gromacs.org/Documentation/Floating_Point_Arithmetic
for discussion on how close it should be to an integer.

Number of degrees of freedom in T-Coupling group rest is 67749.00

NOTE 3 [file ions.mdp]:
You are using a plain Coulomb cut-off, which might produce artifacts.
You might want to consider using PME electrostatics.

There were 3 notes

There was 1 warning

-------------------------------------------------------
Program:     gmx grompp, version 2020.5
Source file: src/gromacs/gmxpreprocess/grompp.cpp (line 2353)

Fatal error:
Too many warnings (1).
If you are sure all warnings are harmless, use the -maxwarn option.

website at http://www.gromacs.org/Documentation/Errors
-------------------------------------------------------
Analysing residue names:
There are:   135    Protein residues
There are: 10727      Water residues
Analysing Protein...
This run will generate roughly 3 Mb of data


This is the link to the ions.mdp file Dr. Lemkul has provided: http://www.mdtutorials.com/gmx/umbrella/Files/ions.mdp.

What is going wrong here? I feel like I am running around in circles with this tutorial...

I appreciate any advice you have to offer!

• As you can see from the output, the warning is generated as you are using GROMOS force field, and because that force field was originally parameterized with incorrect options. If you think everything is ok, just add -maxwarn 1 at the end of the gmx grompp command you have written, and it will run. I have never used GROMOS force fields so I don't know what problems they are referring to, or whether you will get the right results. (You should visit the webpages they mention in the warning) Feb 10 at 12:39
• @ShoubhikRMaiti, thanks for your reply. I believe the command is not working, because it is not outputting ions.tpr Feb 10 at 13:41
• It's not working even if you put -maxwarn 1 ? Feb 11 at 12:19
• it is, i misspoke, sorry Feb 11 at 14:34
• @ShoubhikRMaiti It seems like your comment resolves unless I'm misinterpreting. Do you mind converting it to an answer?
– Tyberius
Jun 23 at 22:05

As I had previously written in the comments, the warning given by gromacs seems to indicate that the GROMOS force fields should not be used due to wrong parameterization:

  The GROMOS force fields have been parametrized with a physically
incorrect multiple-time-stepping scheme for a twin-range cut-off. When
used with a single-range cut-off (or a correct Trotter
multiple-time-stepping scheme), physical properties, such as the density,
might differ from the intended values. Since there are researchers
actively working on validating GROMOS with modern integrators we have not
yet removed the GROMOS force fields, but you should be aware of these
issues and check if molecules in your system are affected before
proceeding. Further information is available at
https://redmine.gromacs.org/issues/2884 , and a longer explanation of our
decision to remove physically incorrect algorithms can be found at
https://doi.org/10.26434/chemrxiv.11474583.v1


(copied from the question)

The redmine link returns a 404 error, but the chemrxiv paper is accessible. From a quick read of the paper, it seems that the GROMACS code inherited its algorithms from the GROMOS87 code. The old algorithm used a "twin range cut-off" scheme, where the short-range forces were calculated in every step, but the long range forces were calculated every N number of steps (N is a certain integer number, probably set by the software). The long-range forces were added to the short-range values every step to approximate a longer cut-off value. This algorithm was written at a time when computational power was not cheap so not calculating long range forces in every step saved a lot of cpu time. This algorithm has been found to cause artifacts especially for lipid bilayers with zwitterionic molecules. Another problem with the algorithm was that it was not time-reversible.[1]

A better algorithm is the Trotter multiple time stepping scheme which was adopted by GROMACS since version 4.5 (2013). However, the old GROMOS force fields were parameterized using the legacy GROMOS algorithm, which means that with the new Trotter algorithm, the GROMOS parameters would give different (and possibly wrong) results. There also seems to be some other issues with the constraints (SHAKE) for the GROMOS parameters.[1] This is why from version 2019.2, GROMACS (the grompp module) prints out a warning whenver the GROMOS force-fields are used, and stops execution.

The only utility of the warning is to ask the user to check whether the system they are using gives correct results with the GROMOS parameters when using the newer algorithm. So, if you get the error, you should check your system before proceeding. Once you are happy with the choice of force-field and system, then just add the argument -maxwarn 1 to the end of the gmx grompp command that you are attempting to run. The command will now run as normal, and ignore the warning. It's as simple as that.

In 2020, Diem and Oostenbrink studied the differences in solvation energy of various amino acid side chains in water and the physical properties of six pure liquids, when using the twin-range algorithm vs the Trotter algorithm. They found that the error in using the old algorithm is mostly minor and within the range of uncertainty in the experimental data used to validate the results. Another important finding was that the error is smaller when using group-based cut-off schemes than for an atom-based cut-off.[2]

References:

1. B. Hess, D. van der Spoel, M. J. Abraham, E. Lindahl, "On The Importance of Accurate Algorithms for Reliable Molecular Dynamics Simulations", https://doi.org/10.26434/chemrxiv.11474583.v1
2. M. Diem, C. Oostenbrink, "The Effect of Using a Twin-Range Cutoff Scheme for Nonbonded Interactions: Implications for Force-Field Parametrization?", J. Chem. Theory Comput., 2020, 16, 10, 5985–5990, https://pubs.acs.org/doi/10.1021/acs.jctc.0c00509