Newest CHARMM36 port for GROMACS

We are pleased to announce an updated port of the CHARMM36 force field for GROMACS, available at http://mackerell.umaryland.edu/charmm_ff.shtml#gromacs

The July 2020 release corresponds to the annual update to the CHARMM36 force field. Notable revisions include:

  1. NBFIX terms for K+ and Ca2+ ions with carboxylates
  2. New nucleic acid model compounds
  3. CGenFF version 4.4, including amide base model compounds, boronic acid and bicyclic boron “ester” model compounds and associated parameters

Please let me know if you encounter any problems using the force field. Happy simulating!

Link: MacKerell Lab
Contact person (name and email address): Justin Lemkul (jalemkul@vt.edu)

How the work has been tested/reviewed: Forces and energies have been compared between CHARMM and GROMACS for a set of isolated molecules (gas phase) including all amino acids, nucleotides, a subset of lipids and other compounds

1 Like

Please be advised of a new “bug fix” release of the CHARMM36 port. There were some missing NBFIXes related to Ca2+ ions in the July 2020 version that have now been correctly added. Please check your results carefully if you used the July 2020 version for simulations involving Ca2+ to see if the missing NBFIXes may have impacted your results.

The tarball is named/dated “feb2021” to reflect the date of its release, but forcefield.doc still refers to July 2020 as the version of the force field, because this does correspond to the release date of the updated CHARMM force field. We maintain access to the July 2020 version on the MacKerell lab website so you can compare the force field files and easily see the changes.

I apologize for this error and will be improving our conversion scripts so that no items like this get omitted in the future. Unfortunately, there’s a fair amount of manual work required to prepare CHARMM-formatted files for conversion, and sometimes copy-paste issues spring up, which I believe was the case here.

1 Like

A new version of the CHARMM36 is available from the MacKerell lab website. Please see MacKerell Lab to get the latest port. Of note are a few new things:

  1. Reorganization of .rtp and .tdb files so terminal patching is now much easier (no more accidental use of 5TER for proteins, hooray!) and all residues now have associated .hdb entries.
  2. Update of CGenFF to version 4.6
  3. Inclusion of a large number of nonstandard amino acids (please note you will need to add these to your residuetypes.dat before they will work natively, as with any nonstandard entity)
  4. Inclusion of recent LJ-PME parameters for lipids. Please note these are in a separate force field port than the “standard” CHARMM36 parameters, since there are differences in terms of topologies and associated parameters for the reparametrized lipids. Users can apply LJ-PME for these lipids in conjunction with the standard CHARMM protein, nucleic acid, etc. force field.
  5. Small bug fixes from previous versions like incorrect impropers in the NME terminus and inconsistent atom naming in the ACE terminus that caused issues in pdb2gmx
  6. It is possible to use a modified water model with TIP3P H ε = 0.1 kcal/mol as in the CHARMM36m force field paper, which strengthens the protein-water interactions and leads to less compact polypeptide structures. There are NBFIXes in place to keep water-water interactions at their standard LJ values. Apply define = -DUSE_MODIFIED_TIP3P_EPS in your .mdp files to make use of this modified water.

If there are any questions or problems, please let me know.

I will conclude by giving a massive shout-out to @awacha who designed a vastly better conversion program for porting CHARMM files over to GROMACS format. It is a significant upgrade over our in-house version and made my life a lot easier in producing the port.

3 Likes

Dear Justin,
Very happy to see the updates !

In regards to point 4 above " Users can apply LJ-PME for these lipids in conjunction with the standard CHARMM protein, nucleic acid, etc. force field." I’ve two queston (1) if I use Charrm-gui to create , say, a lipid membrane and the GUI’s charmm36m, can I assume that it is generally suitable for other non-lipids e.g. carbohydrates ( also made in the GUI) (2) if I use the latest charmm36 in gromacs how can I assure myself that the PME parameters are correct for lipids in the model - what metric do I use to compare ? Or do I cut and paste from the lipid ff and if so just how is this accomplished e.g. as an #include in the top?

Regards,
Paul

Yes, the lipids are the only element of the force field that required reparametrization to be compatible with LJ-PME. I do not know the current status of LJ-PME topologies from CHARMM-GUI but in this case you should generate a topology natively in GROMACS using the LJ-PME version of the port and compare.

Reproduce properties of a pure bilayer.

I don’t understand this question. You should not modify anything about the force field.

Thank you very much for the late night response… Quite useful

The last part of the questions restated is : if I produce a top file natively for most of the molecules in a model, but the model also has lipids, then I should simply import and use the lipid ff version for all molecules. ?

As a related aside, I’ve had what I believe is good success with gmx x2top even with sometimes having to modify atomnames2type, , but the Gromacs definition refers to a “primitive” result. Is there some reason to avoid x2top ?

Paul

There are topological differences in some lipids (e.g. charge) so you need to create the topology for the lipids with the LJ-PME version of the force field.

Please start a new thread in the user forum for general topics like this.

Hi, I’m new here on the forum, so I’m not sure this is the right place to present my problem.

I am using charmm36-jul2021.ff version and when using pdb2gmx

gmx pdb2gmx -f mypdb.pdb -o processed.gro -ignh

I get the following error:

Start terminus THR-333: NH3+
End terminus GLY-526: COO-
Opening force field file ./charmm36-jul2021.ff/aminoacids.arn
Checking for duplicate atoms....
Now there are 1536 atoms. Deleted 6 duplicates.
Generating any missing hydrogen atoms and/or adding termini.

-------------------------------------------------------
Program:     gmx pdb2gmx, version 2020.6
Source file: src/gromacs/gmxpreprocess/pgutil.cpp (line 130)

Fatal error:
Residue 194 named GLY of a molecule in the input file was mapped
to an entry in the topology database, but the atom CB used in
that entry is not found in the input file. Perhaps your atom
and/or residue naming needs to be fixed.

For more information and tips for troubleshooting, please check the GROMACS
website at http://www.gromacs.org/Documentation/Errors
-------------------------------------------------------

I would like to point out that this only happens with charmm36 and with no other force fields available on gromacs 2020.6. Also if I remove the GLY-526 and leave the CYS-525 as the last residue, everything works.

I suspect that it may be due to having GLY as the last residue of the chain, but I am not expert enough to understand if it is my problem, if it is a bug and, if it is, how to fix it.

Best regards,
Leonardo

Hi Leonardo,

the problem is that the COO- terminus entry in aminoacids.c.tdb uses atom “CB” as the third control atom for adding OT1 and OT2. Until it is fixed in the next release of the port, please edit the line "2 8 OT C CA CB " in aminoacids.c.tdb to "2 8 OT C CA N ". That should fix it.

Sorry for the inconvenience,

Andras

Thanks for the prompt reply Andras!
I ran a test with the edit you suggested and everything works.
I am waiting for the new updated version.
Thanks again for your time.

Greetings,
Leonardo

This announcement is a little overdue, but we have updated the CHARMM36 port for GROMACS to reflect the July 2022 version of the CHARMM toppar distribution. You can find the associated files at the same link above (MacKerell lab).

Hopefully we have rectified the few minor issues noted on the forum, and we have also greatly expanded the coverage of the force field to reflect essentially all of what CHARMM36 has to offer.

One quirk to make note of, as I have gotten questions on the forum and via email: In the case of an N-terminal MET residue in a polypeptide sequence, you MUST interactively choose the appropriate terminal patch. pdb2gmx tries to be smart and match by name but in this case, it causes an erroneous attempt to patch with a carbohydrate-specific patch.

We discovered a bug related to the LJ parameters of TIP4P-Ew, which were erroneously copied from TIP4P. The port has been updated, and if you were doing any simulations with C36 and TIP4P-Ew with these files, please re-run them. Apologies for the error.

1 Like

What does this mean " you MUST interactively choose the appropriate terminal patch"?
please explain how to do this. I am just trying to do the mdtutorials protein-ligand complex example.

based on all the posts you have linked I still have no idea what I am supposed to edit for fix this fatal error:
“”"
atom C1 not found in buiding block 1MET while combining tdb and rtp
“”"

You need to use the -ter option and interactively choose from the prompts as to which termini you want. For polypeptides, you typically want NH3+. The MET1 patch is for carbohydrates, but pdb2gmx thinks it’s the best match because it matches the MET residue name. In this case, it’s the program trying to be too smart (though there are reasons for it that are typically associated with GLY residues).

Hi professor @jalemkul ,

I’m so new at using gromacs. I exactly don’t understand what you mean with “-ter option”. I mean i don’t understand how to fix it. I use charmm36-jul2022 too and i get the same fatal error which is “atom C1 not found in building block 1MET while combining tdb and rtp”
I check out the aminoacids.n.tdb file and i see gly term like this:

[ GLY-NH3+ ]
; Glycine N-terminus
[ delete ]
HN
HT1
HN1
HT2
HN2
[ replace ]
N NH3 14.007000 -0.3000
CA CT2 12.011000 0.1300
[ add ]
3 4 H N CA C
HC 1.008000 0.3300 -1

[ NH3+ ]
; standard N-terminus
[ delete ]
HN
HT1
HN1
HT2
HN2
[ replace ]
N NH3 14.007000 -0.3000
CA CT1 12.011000 0.2100
HA HB1 1.008000 0.1000
[ add ]
3 4 H N CA C
HC 1.008000 0.3300 -1

Is this any problem?
best regards.
Gozde.

You need to interactively select termini with by invoking -ter in your pdb2gmx command. This is only true if you have an N-terminal methionine due to an incorrect attempt to match the residue name with an unrelated patch in the termini database.

1 Like

He meant, for example,
gmx pdb2gmx -f 3HTB_clean.pdb -o 3HTB_processed.gro -ter
then you can choose interactively what residues you want to have in each terminus.

As a newbee, it took very long time to understend his mention but probably many of somebody like me…
A Question is for usual protein, N term : NH3+ and C term : COO-?

1 Like

Yes, this is the typical terminus chemistry for polypeptides.

Thank you so sooo much! It works now! :)

sir i am a beginner.
how can i use this CHARMM36 forcfield in my simulation