Pdb2gmx: Fatal error: Atom CA in residue ACE 1 was not found in rtp entry ACE

User discussions
1

GROMACS version: 2024.1
GROMACS modification: No

Hello community.

I am trying to reproduce a tutorial for umbrella sampling.

In it, the preparation of the topology is shown as follows:

gmx pdb2gmx -f 2BEG_model1_capped.pdb -ignh -ter -o complex.gro

What I am trying to do is run that command selecting a CHARMM36 force-field:

gmx pdb2gmx -f 2BEG_model1_capped.pdb -ff charmm36_ljpme-jul2022 -ignh -ter -o complex.pdb

However, that command throws this error:

Fatal error:
Atom CA in residue ACE 1 was not found in rtp entry ACE with 6 atoms
while sorting atoms.

Is there a known solution that obstacle?

I should point out that the exact same error is thrown with other -ff options like the build-in charmm27, amber99sb-ildn, and amber99sb force-fields.

Note: the 2BEG_model1_capped.pdb that I am using is the unmodified file from the tutorial.

Thank you.

Ivan

2

Hello Justin (@jalemkul) and András (@awacha)

I just discovered your forum conversation about the ACE residue and its port from CHARMM36 to GROMACS:

As you can see in my original post, I am facing the same problem with ACE.

From your conversation in 2021, it seems that a solution exists but it must be missing from public access.

Could you or anybody find a little time in your day to post a solution to the GROMACS community?

Thank you.

Ivan

3

GROMOS uses CA as the name for the terminal methyl carbon. Other force fields (like CHARMM) use CH3. Rename the atoms in the input coordinate file to match the expectations of the .rtp file.

4

Thanks @jalemkul for replying.

I think that I got it half-working. So far. This is what I did:

Display atom names in input coordinate file:

$ grep ' ACE ' 2BEG_model1_capped.pdb 
ATOM      1  CA  ACE A   1      18.260  29.485  42.791  1.00  0.00
ATOM      2  C   ACE A   1      19.533  30.318  42.654  1.00  0.00
ATOM      3  O   ACE A   1      20.225  30.102  41.652  1.00  0.00
ATOM    227  CA  ACE B  28      18.598  28.257  38.354  1.00  0.00
ATOM    228  C   ACE B  28      19.355  29.494  37.838  1.00  0.00
ATOM    229  O   ACE B  28      19.372  29.856  36.662  1.00  0.00
ATOM    453  CA  ACE C  55      17.981  28.049  33.725  1.00  0.00
ATOM    454  C   ACE C  55      19.122  28.901  33.156  1.00  0.00
ATOM    455  O   ACE C  55      19.864  28.403  32.305  1.00  0.00
ATOM    679  CA  ACE D  82      19.107  34.836  28.587  1.00  0.00
ATOM    680  C   ACE D  82      18.862  33.443  29.177  1.00  0.00
ATOM    681  O   ACE D  82      18.094  33.199  30.105  1.00  0.00
ATOM    905  CA  ACE E 109      16.105  31.758  24.396  1.00  0.00
ATOM    906  C   ACE E 109      17.560  32.066  24.036  1.00  0.00
ATOM    907  O   ACE E 109      18.099  31.090  23.521  1.00  0.00

Display the atom names in the .rtp file of the forcefield, for the ACE record only:

$  sed -n '/^\[ ACE \]/,/^\[ /p' ./charmm36_ljpme-jul2022.ff/aminoacids.rtp | sed '$d'
[ ACE ]
; acetylated N-terminus
  [ atoms ]
      CH3      CT3 -0.2700   1
     HH31      HA3  0.0900   1
     HH32      HA3  0.0900   1
     HH33      HA3  0.0900   1
       CA       CA  0.5100   2
        O        O -0.5100   2
  [ bonds ]
       CA   CH3
        C    +N
      CH3  HH31
      CH3  HH32
      CH3  HH33
        O     C
  [ impropers ]
       CA   CH3    +N     O
       +N    CA   +CA   +HN

As shown above, there are two atom name discrepancies

input coordinate forcefield
CA CH3
C CA
O O

So I edit my input coordinate file to change CA to CH3 and C to CA for the ACE lines only:

sed -r -e 's/CA ( * ACE )/CH3\1/' -e 's/C ( * ACE )/CA\1/' 2BEG_model1_capped.pdb > 2BEG_model1_capped_atomnamefix.pdb

which looks like this:

$ grep ' ACE ' 2BEG_model1_capped_atomnamefix.pdb
ATOM      1  CH3 ACE A   1      18.260  29.485  42.791  1.00  0.00
ATOM      2  CA  ACE A   1      19.533  30.318  42.654  1.00  0.00
ATOM      3  O   ACE A   1      20.225  30.102  41.652  1.00  0.00
ATOM    227  CH3 ACE B  28      18.598  28.257  38.354  1.00  0.00
ATOM    228  CA  ACE B  28      19.355  29.494  37.838  1.00  0.00
ATOM    229  O   ACE B  28      19.372  29.856  36.662  1.00  0.00
ATOM    453  CH3 ACE C  55      17.981  28.049  33.725  1.00  0.00
ATOM    454  CA  ACE C  55      19.122  28.901  33.156  1.00  0.00
ATOM    455  O   ACE C  55      19.864  28.403  32.305  1.00  0.00
ATOM    679  CH3 ACE D  82      19.107  34.836  28.587  1.00  0.00
ATOM    680  CA  ACE D  82      18.862  33.443  29.177  1.00  0.00
ATOM    681  O   ACE D  82      18.094  33.199  30.105  1.00  0.00
ATOM    905  CH3 ACE E 109      16.105  31.758  24.396  1.00  0.00
ATOM    906  CA  ACE E 109      17.560  32.066  24.036  1.00  0.00
ATOM    907  O   ACE E 109      18.099  31.090  23.521  1.00  0.00

But now, when I run pdb2gmx, I expose a new error:

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

Any pointer to a solution, please?

Thank you.

Ivan

5

I don’t know what’s going on with your ACE definition in the .rtp file but it is incorrect. The standard entry is:

; acetylated N-terminus
  [ atoms ]   
      CH3      CT3 -0.2700   1
     HH31      HA3  0.0900   1
     HH32      HA3  0.0900   1
     HH33      HA3  0.0900   1
        C        C  0.5100   2
        O        O -0.5100   2
  [ bonds ]  
        C   CH3
        C    +N
      CH3  HH31
      CH3  HH32
      CH3  HH33
        O     C
  [ impropers ]
        C   CH3    +N     O
       +N     C   +CA   +HN

You should not have both CA and CH3 in the same residue. The C is the carbonyl carbon and CA/CH3 is the terminal methyl. Simply substitute CA with CH3 in the coordinate file (leave the .rtp alone, and don’t rename C or O) and things should work. The error you’re getting is that the next residue in the chain wants a -C to generate a bond or CMAP or something, and if you’ve renamed it to CA, then it’s not there.