Coordination number from RDF plot

GROMACS version: 2023.1
GROMACS modification: No

Hello

I am simulating ionic liquids with CO2. I want to calculate the co-ordination number of cation and anion i.e the number of CO2 surrounding cation and anion from which solubility of the gas can be evaluated. Co-ordination number (CN) can be calculated from this relation:

CN = 4πρ x ∫r2g(r)dr

NOTE: It is r-squared in the integration not r multiplied by 2. And the integration is within the limit of the first peak of the RDF curve.

From the RDF curve of anion-CO2, we used origin to evaluated integration of the first peak as seen in the attached photo. Then we multiply the results by 4πρ to get the CN.

My questions are:

  1. In the RDF plot we have g(r) against r not r-squared but in the formlua it is r-squared. How to resolve this?

  2. The g(r) has no unit, r is in nm, and ρ is in kg/m^3, then how the units can cancel out to obtain CN (unitless)?

Dear Umar,

  1. You have to multiply g(r) time r^2. You don’t need the RDF as a function of r^2.
  2. I think the density in that formula is a number density, that is, number of molecules/volume, or number of molecules/length^3. Then, r^2 x dr goes as length^3 (as a volume), and you are left with a pure number.

Also, I think if you are looking for the first neighbors you should integrate from r0=0 (so at all effects r0=0.3 is sufficient as there g(r) = 0) to the first minimum in the curve, which seems to be around 0.7, but I might be mistaken, just double check what you want to calculate and its definition.

EDIT: at r=0.3 g(r) is zero, not r.

Dear @obZehn

  1. Multiplying g(r) by r^2 and integrating direct is mathematically not right. That is why we employ a program like Origin,to do the job. But the Origin integrate the curve as it is i.e g(r) and r whereas we need integration of g(r) and r^2 as in the formula.

CN = 4πρ x ∫r^2g(r)dr

How to achieve this?

  1. Thank you for this reasonable explanation.

Yes, I want to get the first neighbors.

Regards.

What do you mean with “is mathematically not right”? What you are doing is integrating over a sphere the g(r) of your system, e.g., you are assuming the g(r) has spherical symmetry, and by looking from a distance r0=0 to r1=first minimum you are actually checking how many particles are nearby in that spherical shell, which corresponds to the first neighbors. Mathematically, that is what the 4pi is (the integration over the whole solid angle), and the r^2 is the usual radial part of the jacobian. Take a look at this and this for the spherical integration.

You can just integrate the curve with some simple code in python or similar (like this for python).

1 Like

Okay, let me re-study these references then apply the best method possible.

Thank You.
Umar