DiffractionPeak¶

class DiffractionPeak(q_value, cutoff_radius=None, included_atoms=None)¶

Class to calculate the intensity of the isotropic diffraction peak at the given q-value.

Parameters:

q_value (PhysicalQuantity of type inverse length) – The q-value at which the diffraction peak should be calculated.
cutoff_radius (PhysicalQuantity of type length) – The cutoff radius for the atomic distances that should be included in the scattering function.
Default: 10.0 Ang.
included_atoms (PeriodicTableElement | str | list) – If not None, this parameter specifies the atoms between which the scattering function should be calculated. Can be specified as element, list of elements, or list of indices.
Default: All atoms.

calculateCVAndDerivatives(configuration)¶

Calculate the value, derivative, and virial of the diffraction peak.

Parameters:: configuration (MoleculeConfiguration | BulkConfiguration | DeviceConfiguration | SurfaceConfiguration | DistributedConfiguration) – The configuration for which the diffraction peak should be calculated.
Returns:: A tuple containing the value of the diffraction intensity (float), the derivatives with respect to atomic positions (PhysicalQuantity of type inverse length), and the derivatives with respect to the cell components (virial) (PhysicalQuantity of type inverse volume).
Return type:: tuple

classmethod diffractionPlotData(configuration, q_values, cutoff_radius=None, included_atoms=None)¶

Helper method to create the X-ray diffraction data for a range of q-values, e.g. to create a plot.

Parameters:

configuration (MoleculeConfiguration | BulkConfiguration | DeviceConfiguration | SurfaceConfiguration) – The configuration for which the scattering plot should be calculated.
q_values (PhysicalQuantity of type inverse length) – The q-values at which the x-ray diffraction should be calculated.
cutoff_radius (PhysicalQuantity of type length) – The cutoff radius for the atomic distances that should be included in the scattering function.
Default: 10.0 Ang.
included_atoms (PeriodicTableElement | str | list) – If not None, this parameter specifies the atoms between which the scattering function should be calculated. Can be specified as element, list of elements, or list of indices.
Default: All atoms.

Returns:

The isotropic X-ray diffraction intensity at the given q-values.

Return type:

numpy.ndarray

Usage Examples¶

Use the DiffractionPeak as collective variable in a SteeredMolecularDynamics to accelerate the crystallization in an amorphous cobalt silicide configuration using the MomentTensorPotential.

# Set up the diffraction peak intensity at q=2.1/Ang as collective variable.
# The diffraction intensity is only calculated for Cobalt atoms.
diffraction_peak = DiffractionPeak(
    cutoff_radius=10.0*Ang,
    q_value=2.1/Ang,
    included_atoms=[Cobalt],
)

# Set up the Steered MD with the diffraction peak intensity as collective variable.
# The object will be used as post-step-hook in the MD simulation.
hook_function = SteeredMolecularDynamics(
    collective_variable=diffraction_peak,
    velocity=300/1000/ps,
    spring_constant=0.5*eV,
    measurement_interval=50,
)

method = NVTNoseHoover(
    initial_velocity=MaxwellBoltzmannDistribution(
        temperature=1500.0*Kelvin
    ),
    reservoir_temperature=1500.0*Kelvin,
    thermostat_timescale=500.0*fs,
)

md_trajectory = MolecularDynamics(
    configuration=last_image_2,
    constraints=[
        FixCenterOfMass()
    ],
    trajectory_filename='CoSi2_bulk_crystallization.hdf5',
    steps=1200000,
    post_step_hook=hook_function,
    measurement_hook=hook_function.measurements,
    log_interval=5000,
    method=method
)

crystallization_example.py

Calculate and plot the diffraction intensities for a SiO2 cristobalite crystal.

# Set up the q-values at which the diffraction intensity should be plotted.
n_points = 100
q_values = numpy.linspace(1.0, 6.0, n_points) / Ang

xray_intensities = DiffractionPeak.diffractionPlotData(sio2_cristobalite, q_values, cutoff_radius=12.0*Ang)
plot_model = Plot.PlotModel(Ang**-1)
line = Plot.Line(q_values, xray_intensities)
plot_model.addItem(line)
plot_model.setLimits()

plot_model.xAxis().setLabel('q')
plot_model.yAxis().setLabel('Intensity')

Plot.show(plot_model)

plot_diffraction_sio2_cristobalite.py

Notes¶

The DiffractionPeak object calculates the intensity of the isotropic X-ray diffraction (XRD) peak at a given scattering vector \(q\). The intensity of selected XRD peaks can be used as collective variable to accelerate the crystallization of bulk and interface materials in biased MD simulations [1].

The diffraction intensity is calculated as [1]:

\[I(Q) = \frac{1}{N} \sum_{i=1}^{N}\sum_{j=1}^{N} f_i(Q) f_j(Q) \frac{\sin(Q \dot R_{ij})}{Q \dot R_{ij}} W(R_{ij})\]

where \(Q\) is the magnitude of the scattering vector, \(N\) is the number of included atoms, \(W(R_{ij})\) is defined as

\[W(R) = \frac{\sin(\pi R / R_c)}{\pi R / R_c}\]

and \(f_i(Q)\) are the form factors which are calculated as described in [2] and XRayScattering.

Currently, the class can be used as collective variable in the SteeredMolecularDynamics hook object.

Before running a crystallization simulation, one should calculate the scattering plot over a range of q-vectors for the targeted crystal material (see second example above) to identify a characteristic peak for that material and use the corresponding q-value for the accelerated simulation. During the accelerated simulations the bias potential will try to change (e.g. increase) the value of the collective variable and thereby drive the configuration e.g. towards a larger value of the diffraction intensity. This normally results in a more pronounced ordering of the the atoms, which facilitates crystallization. Note that this collective variable does not introduce any bias with respect to the orientation of the crystal, as it acts only isotropically. Via the keyword included_atoms the peak intensity calculation can be restricted to a subset of atoms, e.g. only a specific element, or a group of atoms identified by their indices. This can e.g. be used if the difference in scattering intensity between disordered and crystalline structure is more pronounced for a sublattice, or if one wants to study crystallization in the presence of an interface, where the acceleration should only act on the atoms on one side of the interface.

Note, that currently this class does not support MPI domain decomposition in MolecularDynamics(), which is why the MD simulation should be run with OpenMP threading instead.