RadialDistribution¶
- class RadialDistribution(md_trajectory, start_time=None, end_time=None, cutoff_radius=None, resolution=None, pair_selection=None, time_resolution=None, info_panel=None)¶
Class for calculating the radial distribution from an
MDTrajectory
.- Parameters:
md_trajectory (
MDTrajectory
|MoleculeConfiguration
|BulkConfiguration
|DeviceConfiguration
|SurfaceConfiguration
) – The MDTrajectory or configuration the radial distribution should be calculated for.cutoff_radius (PhysicalQuantity of type length) – Upper limit on sampled distances (must be positive). Default: Half the diagonal of the unit cell
start_time (PhysicalQuantity of type time) – The start time. Default:
0.0 * fs
end_time (PhysicalQuantity of type time) – The end time. Default: The last frame time.
resolution (PhysicalQuantity of type length) – The binning of the histogram. Default:
0.05 * Angstrom
pair_selection (sequence) – Only include contributions between this selection of atoms. Either None or a sequence containing two of the following types: Element, tag name, list of indices, or None. Default: All atoms pairs are considered.
time_resolution (PhysicalQuantity of type time) – The time interval between snapshots in the MD trajectory that are included in the analysis.
info_panel (InfoPanel (Plot2D)) – Info panel to show the calculation progress. Default: No info panel
- data()¶
Return the radial distribution function.
- distances()¶
Return the distance values associated with the radial distribution function.
Usage Examples¶
Load an MDTrajectory and calculate the RadialDistribution function between Al and O:
md_trajectory = nlread('alumina_trajectory.nc')[-1]
rdf = RadialDistribution(md_trajectory,
cutoff_radius=8.0*Angstrom,
pair_selection=[Aluminum, Oxygen])
# Get the bin_centers and the histogram of the radial distribution.
distances = rdf.distances().inUnitsOf(Angstrom)
histogram = rdf.data()
# Plot the data using pylab.
import pylab
pylab.plot(distances, histogram, label='Al-O RDF')
pylab.xlabel('r (Ang)')
pylab.ylabel('g(r)')
pylab.legend()
pylab.show()
Calculate the RadialDistribution function between the first atom in the configuration and the rest of the system:
rdf = RadialDistribution(md_trajectory,
cutoff_radius=20.0*Angstrom,
pair_selection=[[0], None])
Notes¶
The radial distribution (or pair correlation) function is calculated as
It is also possible to calculate the partial radial distribution by specifying
two selections A and B (elements or index lists) in pair_selection
.
In that case the radial distribution is calculated as
The normalization via the density ensures that for large distances the radial distribution approaches unity.
As the radial distribution assumes a homogeneous system, this analysis is not well-defined for trajectories of DeviceConfiguration or SurfaceConfiguration types. In that case only the radial distribution of the central region atoms is calculated, using infinite electrodes as boundary conditions.