Auxiliary Modules

You can download source files to follow this tutorial from this link (26 GB).

Note

This is the same example set used in the Hopping Parameter Extraction or Mean Square Displacement section in CLI Tutorial documentation, or the Hopping Parameter Extraction section in API Tutorial documentation. If you have already completed that tutorial, you do not need to download the files again.

This tutorial provides a brief overview of VacHopPy’s auxiliary modules. For more detailed information on each module and its parameters, please consult the API Reference documentation.

First, navigate to the Example3 directory you downloaded. For this tutorial, we will only use the TRAJ_TiO2 directory and the POSCAR_TiO2 file; other files in the directory can be ignored.

import os
from vachoppy.core import Site

path_traj = 'TRAJ_TiO2'
path_structure = 'POSCAR_TiO2'
if not os.path.exists(path_traj): print(f"{path_traj} not found.")
if not os.path.exists(path_structure): print(f"{path_structure} not found.")

First, we’ll create an instance of the Site class. This object will analyze the perfect structure to define the lattice framework for our analysis.

site = Site(path_structure, 'O')

The Vibration Module

The Vibration module is dedicated to calculating the atomic vibration frequency (\(\nu^*\)) from trajectory data. Unlike the Calculator function, this module is designed to operate on a single HDF5 trajectory file.

Let’s instantiate the Vibration class, providing it with the path to a single trajectory file and our pre-defined site object.

from vachoppy.vibration import Vibration

path_traj_single = os.path.join(path_traj, 'TRAJ_2100K', 'TRAJ_O_01.h5')
vib = Vibration(path_traj_single, site)

To run the analysis, simply call the .calculate() method. This will perform the frequency calculation and print a summary of the results directly to the console.

vib.calculate(verbose=True)

====================================================
       High-Frequency Filtering Results (IQR)
====================================================
  - Cutoff Frequency              : 32.14 THz
  - Removed Outlier Frequencies   : 74 (out of 7521)
====================================================
       Vibrational Analysis Results Summary
====================================================
  - Mean Vibrational Amplitude (σ) : 0.191 Ang
  - Determined Site Radius (2 x σ) : 0.383 Ang
  - Total Vibrational Frequencies  : 7447 found
  - Mean Vibrational Frequency     : 13.324 THz
====================================================

Execution Time: 2.647 seconds
Peak RAM Usage: 0.016 GB

The summary provides key metrics, including the mean vibrational frequency and the determined site radius, which is calculated as twice the mean vibrational amplitude.

You can further visualize the results using the plot_frequencies() and plot_displacements() methods to see the distribution of frequencies and atomic displacements, respectively.

vib.plot_frequencies()

Example output:

Image

The Einstein Module

The Einstein module is used to calculate the atomic diffusivity by analyzing the mean squared displacement (MSD) of atoms and applying the Einstein relation.

While the Calculator() focuses on the discrete hopping events of vacancies, the Einstein module tracks the continuous, long-range movement of the atoms themselves. A key distinction is that the diffusivity calculated here inherently includes correlation effects. Consequently, the activation energy derived from this analysis represents the effective diffusion barrier, not the effective hopping barrier.

Similar to the Calculator(), the Einstein class can accept either a path to a single HDF5 file or a directory containing a bundle of files.

from vachoppy.einstein import Einstein

# For a bundle of trajectories
einstein_bundle = Einstein(path_traj, 'O')

# For a single trajectory file
einstein_single = Einstein(path_traj_single, 'O')

The analysis is performed by calling the .calculate() method.

When provided with a single HDF5 file, the module calculates the MSD and determines the diffusivity for that specific temperature. You can view these results using the .summary() method.

einstein_single.calculate()
einstein_single.summary()

============================================================
               MSD Analysis Summary (Single)
============================================================

-- Input Parameters --
  - Trajectory File : TRAJ_O_01.h5
  - Target Symbol   : O
  - Temperature     : 2100.0 K
  - Skip Time       : 0.00 ps
  - Segment Length  : Full

-- Fitting Range --
  - Start Time      : 1.00 ps
  - End Time        : 151.00 ps

-- Results --
  - Diffusivity (D) : 2.070e-11 m^2/s
============================================================

To visualize the underlying data, you can generate a plot of MSD vs. time using the .plot_msd() method.

einstein_single.plot_msd()

Example output:

Image

When provided with a bundle of HDF5 files from different temperatures, the Einstein module automatically performs an Arrhenius fit to calculate the pre-exponential factor and the diffusion barrier.Let’s run the calculation on the einstein_bundle object we created earlier and inspect the summary.

einstein_bundle.calculate()
einstein_bundle.summary()

Analysis complete: 100 successful, 0 failed.
Execution Time: 6.657 seconds
Peak RAM Usage: 0.211 GB
==================================================
          MSD Ensemble Analysis Summary
==================================================
  Temp (K)    Avg. Diffusivity (m^2/s)
----------  --------------------------
      1700                   8.616e-12
      1800                   1.307e-11
      1900                   1.686e-11
      2000                   2.251e-11
      2100                   2.948e-11

-- Arrhenius Fit Results --
  - Activation Energy (Ea) : 0.927 eV
  - Pre-factor (D0)        : 4.914e-09 m^2/s
  - R-squared              : 0.9970
==================================================

The summary now includes the average diffusivity calculated for each temperature, as well as the key results from the Arrhenius fit.

You can visualize the Arrhenius fit using the .plot_D() method

einstein_bundle.plot_D()

Example output:

Image

Similarly, the .plot_msd() method can be used to plot the MSD vs. time curves for all trajectories, grouped by temperature.

einstein_bundle.plot_msd()

Example output:

Image