How to calculated radial distribution function using the centre of mass in a cluster of water

GROMACS version: 2021
GROMACS modification: Yes

I am using Gromacs.2021 compiled on my own pc with Intel MKL and MSVC++/Intel C++ compilers. I performed an MD simulation on a spherical cluster of ~200 H2O molecules kept at the centre of a 6 nm cubic box. The simulation ran fine.

Then I calculated the radial distribution function with gmx rdf -f md.xtc -s md.tpr -n index.ndx -ref 2 -sel 2

In the index file, 2 represents SOL (i.e. all the solvent molecules). I thought this would use the centres of masses of all the H2O molecules and then create a distribution from that. The graph I got was this: image.png - Google Drive . It does not look like what I would expect at all. There are two peaks for some reason. Is it taking the O and H atoms separately for each molecule? How can I direct it to use the centre of mass?

Any help is appreciated!

Hi,
the tool gmx rdf uses the position of atoms. To specify differently, you can use other options (e.i `-selrpos’) see more in
https://manual.gromacs.org/current/onlinehelp/gmx-rdf.html
Best regards
Alessandra

Thanks! I have managed to run H2O-H2O rdf.

However, now I want to run an rdf of all O atoms from the centre of mass of the whole unit cell. I cannot figure out how to do this with -selrpos. I have read the manual for the rdf page and selection options page.

I tried this:

gmx rdf -f md_0_1.xtc -s md_0_1.tpr -n index.ndx -ref 0 -sel 3 -selrpos whole_mol_com

But I don’t think it is giving the correct result. (In the index, group 0 is System, group 3 is all O atoms.)

Could you please help?

Hi,
The reason is that
whole_ calculates the centers for the whole residue/molecule, even if only part of it is selected, and not of the whole system (many waters molecules).

What about to calculate the rdf between one O and all the other options?

Best regards
Alessandra

Thanks!

So how do I get gmx rdf to calculate the centre of mass of the whole system. Would dyn_mol_com do the job? (I have tried looking through the gromacs manual but the documentation is minimal and difficult to understand for me)

I can’t choose one O and calculate rdf from that because I am calculating rdf over a trajectory, so I can’t choose one O because it won’t stay at the centre after a few frames.

Also do I need the -seltype option? I don’t understand what it does.

Hi,

I forgot that you have a spherical cluster of waters. My suggestion made sense for a box of water molecules.

This link may help you to understand the difference between -seltype and
-selrpos
https://manual.gromacs.org/current/onlinehelp/selections.html?highlight=selection

The option -ref combined with appropriate (without using groups in the index file) may allow you to archive what you want. For syntax see link above.

I hope it helps
Alessandra

Hi srm,

How did you calculate water-water RDF, let me know the full gmx rdf code you use for that?
I am also facing the same problem, which you faced.

Thank you in Advance!

Hi Kal,

In the end I decided not to use gmx rdf, and instead used MDAnalysis python library to read the trajectory and manually calculate the RDF (i.e. using the formula for rdf).

You can also try VMD which can do RDF plots.

Best wishes,

SRM

Hi, I ma not familiar with MDAnalysis. Can we use GROMACS data files to calculate RDF from MDAnalysis?

Cheers, Kal

Yes, MDAnalysis reads gromacs data. Check this page: https://userguide.mdanalysis.org/stable/examples/quickstart.html