# Generation of topology for non-standard molecules

**URL:** <https://gromacs.bioexcel.eu/t/generation-of-topology-for-non-standard-molecules/2789>\
**Category:** User discussions\
**Tags:** forcefield\
**Created:** [September 9, 2021, 11:39am UTC](https://gromacs.bioexcel.eu/t/generation-of-topology-for-non-standard-molecules/2789 "2021-09-09T11:39:15Z")\
**Posts on this page:** 20\
**Page:** 1

<div class="post-metadata">

**Author:** ![Gourav](https://avatars.discourse-cdn.com/v4/letter/g/278dde/32.png) [@Gourav](https://gromacs.bioexcel.eu/u/Gourav)\
**Post date:** [September 9, 2021, 11:39am UTC](https://gromacs.bioexcel.eu/t/generation-of-topology-for-non-standard-molecules/2789/1 "2021-09-09T11:39:16Z")

</div>

GROMACS version: gmx, version 2020.1-Ubuntu-2020.1-1  
GROMACS modification: Yes/No  
Here post your question I have trying to generate my itp and top file using the pdb2gmx command. But as my structure does not have any parameters in the .rtp file of the forcefield I am getting an error.  
So I want to know how can I add a residue into the gromacs forcefield, so that I can generate the .itp and .top.  
I know the values of all the parameters required, like the force constants of each of each of the bonds, or the angles or dihedral parameters etc.  
I looked into gromacs documentation, but it is not quite clear. If someone could show an example on how to do this, then it would be very helpful.

---

<div class="post-metadata">

**Author:** ![alevilla](https://dub1.discourse-cdn.com/flex017/user_avatar/gromacs.bioexcel.eu/alevilla/32/439_2.png) [@alevilla](https://gromacs.bioexcel.eu/u/alevilla)\
**Post date:** [September 10, 2021, 12:14pm UTC](https://gromacs.bioexcel.eu/t/generation-of-topology-for-non-standard-molecules/2789/2 "2021-09-10T12:14:03Z")

</div>

Hi,  
The information may be a little scattered in GROMACS manual and example may be useful.  
I guess you have looked at [Adding a Residue to a Force Field — GROMACS 2021.3 documentation](https://manual.gromacs.org/current/how-to/topology.html#adding-a-new-residue)  
Also maybe the tool x2top or other online/offline tools may be useful to build the step (connectivity for the new residue) that you have to plug in the rtp file. You can also look at the third party category in the forum ( [Third party tools and files - GROMACS forums](https://gromacs.bioexcel.eu/c/third-party-tools-and-files/8)) if some useful tool is there.  
\Alessandra

---

<div class="post-metadata">

**Author:** ![Gourav](https://avatars.discourse-cdn.com/v4/letter/g/278dde/32.png) [@Gourav](https://gromacs.bioexcel.eu/u/Gourav)\
**Post date:** [September 10, 2021, 12:41pm UTC](https://gromacs.bioexcel.eu/t/generation-of-topology-for-non-standard-molecules/2789/3 "2021-09-10T12:41:51Z")

</div>

I tried x2top with a simple methane molecule with the values of all the parameters as given in the instructions present in the link [gmx x2top — GROMACS 5.1.1 documentation](https://manual.gromacs.org/5.1.1/onlinehelp/gmx-x2top.html)  
Here is my pdb file for methane:  
COMPND UNNAMED  
AUTHOR GENERATED BY OPEN BABEL 2.3.90  
HETATM 1 C UNL 1 -2.432 0.404 -0.000 1.00 0.00 C  
HETATM 2 H UNL 1 -1.323 0.404 0.000 1.00 0.00 H  
HETATM 3 H UNL 1 -2.802 0.634 -1.020 1.00 0.00 H  
HETATM 4 H UNL 1 -2.802 1.173 0.709 1.00 0.00 H  
HETATM 5 H UNL 1 -2.802 -0.594 0.312 1.00 0.00 H  
CONECT 1 2 3 4 5  
CONECT 2 1  
CONECT 3 1  
CONECT 4 1  
CONECT 5 1  
MASTER 0 0 0 0 0 0 0 0 5 0 5 0  
END

and this is what I got after I gave the x2top command  
gmx x2top -f methane.pdb -o out.top -r out.rtp -v -nexcl 3 -alldih no -remdih yes -name UNL -param yes -pdbq no -round yes -kb 10000 -kt 400 -kp 5

WARNING: all CONECT records are ignored  
Opening force field file /usr/share/gromacs/top/oplsaa.ff/atomname2type.n2t  
There are 23 name to type translations in file oplsaa.ff  
Generating bonds from distances…  
atom 5  
Can not find forcefield for atom C-1 with 4 bonds

* * *

Program: gmx x2top, version 2020.1-Ubuntu-2020.1-1  
Source file: src/gromacs/gmxpreprocess/x2top.cpp (line 188)

Fatal error:  
Could only find a forcefield type for 4 out of 5 atoms

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

---

<div class="post-metadata">

**Author:** ![alevilla](https://dub1.discourse-cdn.com/flex017/user_avatar/gromacs.bioexcel.eu/alevilla/32/439_2.png) [@alevilla](https://gromacs.bioexcel.eu/u/alevilla)\
**Post date:** [September 13, 2021, 8:38am UTC](https://gromacs.bioexcel.eu/t/generation-of-topology-for-non-standard-molecules/2789/4 "2021-09-13T08:38:11Z")

</div>

Hi,  
Is the molecule you want to describe is a methane or similar? Maybe it is easy to build directly an itp file and include the parameters there. Maybe someone else has already a template for methane.  
\Alessandra

---

<div class="post-metadata">

**Author:** ![Gourav](https://avatars.discourse-cdn.com/v4/letter/g/278dde/32.png) [@Gourav](https://gromacs.bioexcel.eu/u/Gourav)\
**Post date:** [September 13, 2021, 9:34am UTC](https://gromacs.bioexcel.eu/t/generation-of-topology-for-non-standard-molecules/2789/5 "2021-09-13T09:34:10Z")

</div>

No…mine is a complicated molecule. I was just trying with a simple example… Just wanted to about the procedure. One thing I want to ask is… in the rtp files for gromos, there are terms like gb\_42 or ga\_14 in the bonds and angles section. What does that mean? Here I am sharing that portion…  
[ALA]  
[atoms]  
N N -0.31000 0  
H H 0.31000 0  
CA CH1 0.00000 1  
CB CH3 0.00000 1  
C C 0.450 2  
O O -0.450 2  
[bonds]  
N H gb\_2  
N CA gb\_21  
CA CB gb\_27  
CA C gb\_27  
C O gb\_5  
C +N gb\_10  
[angles]  
; ai aj ak gromos type  
-C N H ga\_32  
-C N CA ga\_31  
H N CA ga\_18  
N CA CB ga\_13  
N CA C ga\_13  
CB CA C ga\_13  
CA C O ga\_30  
CA C +N ga\_19  
O C +N ga\_33  
[impropers]  
; ai aj ak al gromos type  
N -C CA H gi\_1  
CA N C CB gi\_2  
C CA +N O gi\_1  
[dihedrals]  
; ai aj ak al gromos type  
-CA -C N CA gd\_14  
-C N CA C gd\_39  
N CA C +N gd\_40  
as shown above…gb\_27, gd\_14… if someone could answer me even that then also I think I can try generating the topology

---

<div class="post-metadata">

**Author:** ![alevilla](https://dub1.discourse-cdn.com/flex017/user_avatar/gromacs.bioexcel.eu/alevilla/32/439_2.png) [@alevilla](https://gromacs.bioexcel.eu/u/alevilla)\
**Post date:** [September 13, 2021, 10:57am UTC](https://gromacs.bioexcel.eu/t/generation-of-topology-for-non-standard-molecules/2789/6 "2021-09-13T10:57:49Z")

</div>

Hi,  
the labels gb\_xx ga\_xx gd\_xx corresponds to the bond, angle and dihedral parameters defined in the ffbonded.itp. You find the file ffbonded.itp in the share/gromacs/top/{forcefield\_name}.ff/ffbonded.itp where {forcefield\_name} is the force field you are looking at.  
ffbonded.itp and ffnonbonded.itp are included in forcefield.itp  
I hope it helps  
Alessandra

---

<div class="post-metadata">

**Author:** ![Gourav](https://avatars.discourse-cdn.com/v4/letter/g/278dde/32.png) [@Gourav](https://gromacs.bioexcel.eu/u/Gourav)\
**Post date:** [September 14, 2021, 12:26pm UTC](https://gromacs.bioexcel.eu/t/generation-of-topology-for-non-standard-molecules/2789/7 "2021-09-14T12:26:57Z")

</div>

I tried adding the parameters in the .rtp file as well as the non bonded parameters in the corresponding ffnonbonbed.itp file but still I am facing an error. I am working with a silica slab, so my system has large number of atoms. Now, as far as I know the parameters for a unit cell should work for the entire periodic lattice, but when I am running the pdb2gmx command, only four atoms are presented in the final gro and top file, all other atoms are removed. Moreover, the bonding between these files are also getting removed. When viewed in VMD I am only getting four isolated atoms, 1 Si, 2 O atoms and 1 H atom

---

<div class="post-metadata">

**Author:** ![jalemkul](https://dub1.discourse-cdn.com/flex017/user_avatar/gromacs.bioexcel.eu/jalemkul/32/18_2.png) [@jalemkul](https://gromacs.bioexcel.eu/u/jalemkul)\
**Post date:** [September 14, 2021, 12:49pm UTC](https://gromacs.bioexcel.eu/t/generation-of-topology-for-non-standard-molecules/2789/8 "2021-09-14T12:49:22Z")

</div>

To create a topology for a surface or other crystalline material that has non-linear connectivity, you should not be using `pdb2gmx`, you should be using `x2top`.

---

<div class="post-metadata">

**Author:** ![Gourav](https://avatars.discourse-cdn.com/v4/letter/g/278dde/32.png) [@Gourav](https://gromacs.bioexcel.eu/u/Gourav)\
**Post date:** [September 14, 2021, 1:13pm UTC](https://gromacs.bioexcel.eu/t/generation-of-topology-for-non-standard-molecules/2789/9 "2021-09-14T13:13:16Z")

</div>

Ok… Is the procedure same for both? I am stating what I want to achieve here. My purpose is to generate a silica slab (which has a periodicity in it’s structure). Hence, I want to keep the terminal atoms(present at the edge of the simulation box) to have their connectivities open, as they would be generating a continuous sequence of atoms when pbc conditions are applied.  
Now, as my system has a large number of atoms I simply cannot use the external servers because they have a limitation on the number of atoms that their servers can process. Hence, I am struck with generating the topology manually.  
I have added the atoms, bonds, angles and dihedrals in the .rtp file and the LJ parameters in the nonbonded.itp file. What should I do next?

---

<div class="post-metadata">

**Author:** ![jalemkul](https://dub1.discourse-cdn.com/flex017/user_avatar/gromacs.bioexcel.eu/jalemkul/32/18_2.png) [@jalemkul](https://gromacs.bioexcel.eu/u/jalemkul)\
**Post date:** [September 14, 2021, 5:18pm UTC](https://gromacs.bioexcel.eu/t/generation-of-topology-for-non-standard-molecules/2789/10 "2021-09-14T17:18:08Z")

</div>

Explicit inclusion of the parameters in the `.rtp` file is not required. All bonded and nonbonded terms can be read from `ffbonded.itp` and `ffnonbonded.itp`.

---

<div class="post-metadata">

**Author:** ![Gourav](https://avatars.discourse-cdn.com/v4/letter/g/278dde/32.png) [@Gourav](https://gromacs.bioexcel.eu/u/Gourav)\
**Post date:** [September 14, 2021, 5:38pm UTC](https://gromacs.bioexcel.eu/t/generation-of-topology-for-non-standard-molecules/2789/11 "2021-09-14T17:38:24Z")

</div>

Ok… So I can just add the parameters in the ffbonded.itp and ffnonbonded.itp?  
I tried adding the parameters in the corresponding files and ran the pdb2gmx command but only 4 atoms were present in the final gro file and all other atoms (from the original silica slab) were removed. What was more annoying was that those 4 atoms were not connected by any bonds. All I got was 4 atoms isolated in space.  
Here I am attaching the files I added in the forcefield files:  
1)Added in the .rtp file:  
[SIL]  
[atoms]  
; name type charge chargegroup  
SI SI 1.020 0  
OW OW -0.510 1  
H H 0.255 0  
O O -0.993 1  
[bonds]  
; atom1 atom2 c0 c1  
SI OW 0.1720 2.7042e+06  
O H 0.0938 2.9551e+07  
SI O 0.1690 1.9257e+06  
[angles]  
; atom1 atom2 atom3 th0 cth  
SI O SI 166.00 64885.36  
SI OW H 119.00 2211.40  
O SI O 180.00 500.00  
OW SI O 109.50 450.00  
[dihedrals]  
; atom1 atom2 atom3 atom4 phi0 cp mult  
H OW SI O 180.00 1.00 3  
O SI O SI 0.00 1.05 3

1. Added in the ffnonbonded.itp file:  
SI 14 0.000 0.000 A 5.803853e-03 5.35358e-06  
O 8 0.000 0.000 A 2.24192e-03 1.36137e-06  
OW 8 0.000 0.000 A 2.24192e-03 1.36137e-06

2. the final gro file that I got:  
3D Atomistic  
4  
1SIL SI 1 2.169 7.014 -0.182  
1SIL OW 2 6.337 -0.752 -0.456  
1SIL H 3 0.494 5.754 -0.455  
1SIL O 4 2.952 -1.104 -0.091  
5.84310 8.11840 0.36460

I am attaching all these for your reference. The pdb file I am attaching with this message (CONVERTED TO .dat format)[topol.top](https://gromacs.bioexcel.eu/uploads/short-url/h21LZscq8rwQbKFBVtmpBaTLwQL.top) (2.4 KB)

[silica\_slab\_edited\_OW.dat](https://gromacs.bioexcel.eu/uploads/short-url/3crDNPgjDsp8MPnML5Sj5Nke9Y9.dat) (718.4 KB)

---

<div class="post-metadata">

**Author:** ![jalemkul](https://dub1.discourse-cdn.com/flex017/user_avatar/gromacs.bioexcel.eu/jalemkul/32/18_2.png) [@jalemkul](https://gromacs.bioexcel.eu/u/jalemkul)\
**Post date:** [September 15, 2021, 1:43am UTC](https://gromacs.bioexcel.eu/t/generation-of-topology-for-non-standard-molecules/2789/12 "2021-09-15T01:43:19Z")

</div>

> [@Gourav](#):
>
> Ok… So I can just add the parameters in the ffbonded.itp and ffnonbonded.itp?

Correct.

> [@Gourav](#):
>
> I tried adding the parameters in the corresponding files and ran the pdb2gmx command

Like I said before, this is not the correct approach; it will never work. Use `x2top` with a suitable `.n2t` file.

---

<div class="post-metadata">

**Author:** ![Gourav](https://avatars.discourse-cdn.com/v4/letter/g/278dde/32.png) [@Gourav](https://gromacs.bioexcel.eu/u/Gourav)\
**Post date:** [September 15, 2021, 2:55am UTC](https://gromacs.bioexcel.eu/t/generation-of-topology-for-non-standard-molecules/2789/13 "2021-09-15T02:55:11Z")

</div>

Ok. So basically the pdb2gmx command is applicable for non-periodic systems. If I were to use a simple ligand, then this would have worked, right?

And one more thing, I tried looking for the format of n2t file…This is what came up:  
H H 0.408 1.008 1 O 0.1  
O OA -0.674 15.9994 2 C 0.14 H 0.1  
C CH3 0.000 15.035 1 C 0.15  
C CH0 0.266 12.011 4 C 0.15 C 0.15 C 0.15 O 0.14

I didn’t quite understand what column 6 and onwards mean? And my research paper mentions about specific parameters for the LJ potential. How can I add these into my forcefield? I mean when I run x2top, does this take only the .n2t files or does it also takes the parameters given in the .rtp and ffbonded.itp and ffnonbonded.itp files as well?

---

<div class="post-metadata">

**Author:** ![Gourav](https://avatars.discourse-cdn.com/v4/letter/g/278dde/32.png) [@Gourav](https://gromacs.bioexcel.eu/u/Gourav)\
**Post date:** [September 15, 2021, 7:10am UTC](https://gromacs.bioexcel.eu/t/generation-of-topology-for-non-standard-molecules/2789/14 "2021-09-15T07:10:27Z")

</div>

I tried the x2top command. Here is what I wrote in the atomname2types.n2t file:  
H H 0.2550 1.008 1 OW 0.0938  
SI SI 1.0200 28.08 4 O 0.1690 OW 0.1720 O 0.1690 O 0.1690  
OW OW -0.510 15.9994 2 SI 0.1720 H 0.0938  
O O -0.993 15.9994 2 SI 0.1720 SI 0.1720

This is the error that I got… I have modified my pdb file this time. Please have a look (again I have uploaded it in slab.dat )[slab.dat](https://gromacs.bioexcel.eu/uploads/short-url/5FfPgcKwnowGbFFK00qvzi6hpQg.dat) (684.2 KB)

Program: gmx x2top, version 2021.3  
Source file: src/gromacs/gmxpreprocess/x2top.cpp (line 186)

Fatal error:  
Could only find a forcefield type for 2560 out of 4608 atoms

## For more information and tips for troubleshooting, please check the GROMACS website at [Common Errors — GROMACS webpage https://www.gromacs.org documentation](http://www.gromacs.org/Documentation/Errors)

* * *

MPI\_ABORT was invoked on rank 0 in communicator MPI\_COMM\_WORLD  
with errorcode 1.

NOTE: invoking MPI\_ABORT causes Open MPI to kill all MPI processes.  
You may or may not see output from other processes, depending on  
exactly when Open MPI kills them.

---

<div class="post-metadata">

**Author:** ![jalemkul](https://dub1.discourse-cdn.com/flex017/user_avatar/gromacs.bioexcel.eu/jalemkul/32/18_2.png) [@jalemkul](https://gromacs.bioexcel.eu/u/jalemkul)\
**Post date:** [September 15, 2021, 4:47pm UTC](https://gromacs.bioexcel.eu/t/generation-of-topology-for-non-standard-molecules/2789/15 "2021-09-15T16:47:52Z")

</div>

> [@Gourav](#):
>
> I didn’t quite understand what column 6 and onwards mean?

See [https://manual.gromacs.org/current/reference-manual/file-formats.html#n2t](https://manual.gromacs.org/current/reference-manual/file-formats.html#n2t)

> [@Gourav](#):
>
> And my research paper mentions about specific parameters for the LJ potential. How can I add these into my forcefield?

Pair-specific interactions are listed in `[nonbond_params]` in `ffnonbonded.itp`.

> [@Gourav](#):
>
> Fatal error:  
> Could only find a forcefield type for 2560 out of 4608 atoms

This means your `.n2t` file does not cover some of the bonds that need to be made.

---

<div class="post-metadata">

**Author:** ![Gourav](https://avatars.discourse-cdn.com/v4/letter/g/278dde/32.png) [@Gourav](https://gromacs.bioexcel.eu/u/Gourav)\
**Post date:** [September 15, 2021, 5:48pm UTC](https://gromacs.bioexcel.eu/t/generation-of-topology-for-non-standard-molecules/2789/16 "2021-09-15T17:48:18Z")

</div>

> [@Gourav](#):
>
> 2560 out of 4608 atom

One question I have to ask here. My system is a periodic one right. So do I have to cap the terminal Si atoms with any H or any ghost atom or sort? Because the terminal atoms still have their valencies open (as they are periodic cell and they will get connected when the pbc condition is applied).

---

<div class="post-metadata">

**Author:** ![jalemkul](https://dub1.discourse-cdn.com/flex017/user_avatar/gromacs.bioexcel.eu/jalemkul/32/18_2.png) [@jalemkul](https://gromacs.bioexcel.eu/u/jalemkul)\
**Post date:** [September 15, 2021, 6:18pm UTC](https://gromacs.bioexcel.eu/t/generation-of-topology-for-non-standard-molecules/2789/17 "2021-09-15T18:18:24Z")

</div>

`x2top` applies periodicity, so if that’s the way you want to model them in the simulation, that’s how the topology has to be created, not with any capping groups.

---

<div class="post-metadata">

**Author:** ![Gourav](https://avatars.discourse-cdn.com/v4/letter/g/278dde/32.png) [@Gourav](https://gromacs.bioexcel.eu/u/Gourav)\
**Post date:** [September 16, 2021, 6:00am UTC](https://gromacs.bioexcel.eu/t/generation-of-topology-for-non-standard-molecules/2789/18 "2021-09-16T06:00:13Z")

</div>

I tried the x2top command but I got an error. Actually I had my terminal silicon atoms uncapped, which resulted in an error, but when I capped the terminal atoms with H, the command ran successfully and I got the corresponding the .top and .rtp file. So now I am quite unsure about how the command works

---

<div class="post-metadata">

**Author:** ![jalemkul](https://dub1.discourse-cdn.com/flex017/user_avatar/gromacs.bioexcel.eu/jalemkul/32/18_2.png) [@jalemkul](https://gromacs.bioexcel.eu/u/jalemkul)\
**Post date:** [September 16, 2021, 12:32pm UTC](https://gromacs.bioexcel.eu/t/generation-of-topology-for-non-standard-molecules/2789/19 "2021-09-16T12:32:16Z")

</div>

There’s not a lot I can tell based on the information provided except that `x2top` reads through each atom in the coordinate file and attempts to assign atom type, charge, and bonded connectivity based on the information you have provided in the `.n2t` file. Distance tolerances are ±10% as is standard for all GROMACS distance searches. If it does not find a match, you get an error. Likely your uncapped structure had interatomic distances that were beyond the threshold.

---

<div class="post-metadata">

**Author:** ![Gourav](https://avatars.discourse-cdn.com/v4/letter/g/278dde/32.png) [@Gourav](https://gromacs.bioexcel.eu/u/Gourav)\
**Post date:** [September 16, 2021, 1:01pm UTC](https://gromacs.bioexcel.eu/t/generation-of-topology-for-non-standard-molecules/2789/20 "2021-09-16T13:01:27Z")

</div>

So let me get this straight, all I have to do is add my atom name and the corresponding atom types, followed by their charge, masses and connectivities (in terms of bond lengths) and I am good to go. Is that right?  
And if I have certain atoms which are not connected to other atoms the way I mentioned it, then it will throw an error. Is that right?  
Then can it be possible that the terminal atoms which have their connectivities open (as they are connected to the atoms of the next box) might be leading to the problem that I am facing?

[Next page](https://gromacs.bioexcel.eu/t/generation-of-topology-for-non-standard-molecules/2789.md?page=2)
