The peak of partial RDF calculated by ovito

Hello, users:
I have used the coordination analysis modifier to calculate the O-Zr of RDF in ZrO2.
Since there are two categories of Zr-O bonds in my cell and their numbers should be accurately same, so i think their corresponding peak should show the same height.
But the fact is, not.


So, why the first two red peaks have different height?

Could you please provide your input file and your coordination modifier settings? I will have a look into it.

1 Like

Of course, how nice you are.
ZrO.cif (17.6 KB)
I calculated the partial RDF by simply using the OVITO basic software, the parameters is default.
cutoff radius:4
bins:400

The value of the Radial Distribution Function (RDF) at a given r does not directly correspond to the number of neighbors at that distance. According to the RDF definition:


(CC-BY-SA Source)

You would need to reverse the normalization by dV_r \times \rho to convert from the RDF to the actual coordination number. Since dV_r does not increase linearly with r, you observe the decrease in g(r).

If you’re interested in calculating the actual number of particles in given distance, you can use the following Python modifier:

import numpy as np
from ovito.data import CutoffNeighborFinder, DataCollection
from ovito.pipeline import ModifierInterface
from traits.api import Float, Int


class NeighborCountModifier(ModifierInterface):
    num_bins = Int(400, label="Number of bins")
    r_max = Float(3.5, label="Cutoff radius")

    def modify(self, data: DataCollection, *, frame: int, **kwargs):

        bins = np.linspace(0, self.r_max, self.num_bins)

        ptypes = data.particles.particle_types

        finder = CutoffNeighborFinder(self.r_max, data)
        result = {}
        for index in range(data.particles.count):
            ptype_i = ptypes[index]
            for neigh in finder.find(index):
                if neigh.index < index:
                    continue

                ptype_j = ptypes[neigh.index]
                if ptype_i < ptype_j:
                    key = (ptype_i, ptype_j)
                else:
                    key = (ptype_j, ptype_i)

                if key not in result:
                    result[key] = np.zeros(len(bins))
                result[key][np.digitize(neigh.distance, bins)] += 1

        result_array = np.zeros((len(result), self.num_bins))
        tags = []
        for i, (k, v) in enumerate(result.items()):
            result_array[i] = v
            tags.append(
                f"{data.particles.particle_types.type_by_id(k[0]).name}-{data.particles.particle_types.type_by_id(k[1]).name}"
            )
        table = data.tables.create(
            identifier="neighbor_counts", title="Neighbor_Counts"
        )
        table.x = table.create_property("distance", data=bins)
        table.y = table.create_property("count", data=result_array.T, components=tags)

1 Like

Thank you,utt. Your reply helps me a lot giving me a deeper understanding for RDF. Wishing you a bright and shining sun today