analyze_md - Trebonius91/utils4VASP GitHub Wiki

This program evaluates MD trajectories of surface slabs and bulk systems

The analysis is always done on an XDATCAR file.

To start the program, simply type the name and a list of possible commands appears. Since this program is quite big in scope, the commands are grouped into main commands activating a certain preprocessing or analysis option and additional commands that are used to fine-tune the chosen evaluations.

The main commands are the following:

  • -part_number=[number] If the trajectory shall be evaluated as a whole or if all analysis tasks shall be done on N separate parts on a trajectory, which can be useful to study the time-dependent behavior of some properties, to see if the system has evolved in a certain way or switched its phase during the simulation time. This keyword is optional, if its not given, the trajectory will be evaluated as one part. If all is analyzed in one part, the results are written to the folder analysis, else, they are written in folders analysis_part1 to analysis_partN. The handling of all following keywords is independent of the choice of parts.
  • -rdf Activates the calculation of radial distribution function (RDF) calculation. RDFs are calculated for each element combination in the system and stored into a folder RDFs within the analysis folder. For example, for a system containing Ga and Cu atoms, the files RDF_Ga-Ga.dat, RDF_Cu-Cu.dat, RDF_Ga-Cu.dat and RDF_Cu-Ga.dat are written. The last two files contain the same RDF profile. Further, integrated RDF profiles are written to file, that contain the average number of elements in a sphere around any atom of another element up to different sphere radii. Of these, also four files are written int the RDFs folder, for example RDF_int_Ga-Cu.dat. Here, files with different orderings of elements are usually not identical, because it makes a difference, how many Cu atoms are around a Ga atom and how many Ga atoms are around a Cu atoms, if the total numbers of them are different.
  • -dens_elems The one-dimensional densities of all elements in the systems are plotted either along a coordinate axis, for example to investigate the surface or interface structures of slabs or membranes, or around the center of mass of one or several elements. The details of this analysis need to be defined with the -dens_mode keyword.
  • -diffusion The three-dimensional self-diffusion coefficients of all elements present in the slab are calculated from the trajectory. The calculation is done via the mean-square displacement (MSD) calculation and averaged over all atoms of the respective element. For surface slabs, also 2D diffusion coefficients can be calculated, by setting the keyword -diff_2d (see below).
  • -dt=[value] The MD time step in fs between two frames of the trajectory. This can be identical to the time step of the MD simulation or it is a multiple of it, if only every N'th frame was written during the dynamics (ML_OUTBLOCK command in VASP).
  • -track_atoms=[list of numbers] The time-dependent x-,y- and z-coordinates of the chosen atoms during the dynamics (indices in the trajectory file, beginning from 1) are written in plot files in the analysis folder(s), one file for each atom. The files are called track_atom[number].dat. Example: -track_atoms=2,28,107 would produce three files, track_atom2.dat, track_atom28.dat and track_atom107.dat. Each line contains the step number (or time, if the -dt keyword is given, see above) and the x-, y- and z-coordinates of the respective atom.
  • dens_cube 3D spatially resolved atom densities are calculated for all elements in the system and stored as Gaussian Cube files in the analysis folder(s), for example into density_Ga.cube. The contained densities can be visualized, e.g., with VMD. Open the file with vmd, then go to Representations -> Create Rep, then in the new representation to Draw style -> Drawing Method -> Isosurface. Change Draw to Solid Surface and modify Isovalue with the mouse bar until the shape looks good.
  • -dft_element=[element] If templates for DFT single point calculations (for example CLS or Bader chares) with the atoms of the chosen elements at several different z-values (near the surface or rather into the bulk etc) shall be made, give the symbol of the element here, for example -dft_element=Pt. After evaluation of the trajectory, a folder dft_input is generated, with subfolders named round1 to roundN, where N is the total number of dft_rounds, see below. In each of these subfolders, a POSCAR file named POSCAR_slice1 to POSCAR_sliceM (where M is the number of dft_slices, see below) is written which can be used as POSCAR file for a following DFT single point calculation.
  • -tension The surface tension of a surface slab is calculated, averaged over the whole trajectory. The calculated value is printed on the shell after the calculation. For this, the system is assumed to have vacuum introduced along the z-axis! Else, useless values will result.
  • -stability Determine if any chemical bond breaks during the MD trajectory and if yes, which bond is it at what time. In the first frame, chemical bonds are detected based on the van der Waals radii of the involved elements. Then, in the remaining dynamics it is checked if any of the detected bonds got elongated by more than a factor of 1.8 which equals decomposition of the bond. If this happens, the indices and elements of the connected atoms and the number of the MD frame are written out. This can be helpful for automatically determining the stability of a system or to check what the most common decomposition/reaction pathway for a system could be.

Besides these main commands, a couple of additional commands can be used to fine-tune the evaluation:

  • -frame_first=[number] Determines the first frame of the XDATCAR file to be considered for the analysis. If for example -frame_first=100 is chosen, the first 99 frames are skipped. Default value is 1.
  • -frame_last=[number] Determines the last frame of the XDATCAR file to be considered for the analysis. Default value is the total number of frames in XDATCAR.
  • -write_traj If in addition to the analysis, all trajectory frames shall be written to the file trajectory.xyz, which might be helpful for visualization together with analysis.
  • -rdf_bins=[number] Number of bins along the atom-atom distance, for which RDF values are calculated. Default is 201.
  • -rdf_cutoff=[value] The maximum distance (in Angstroms) for which the RDF is calculated. It always starts from 0 Angstroms, up to the given value. Default is 8 Angstroms.
  • -diff_collect If -diffusion is activated, it will be calculated by treating all atoms of each element as a single particle, such that, for example, only the dynamics of the center of mass of a particle consisting of Ga atoms, and not its intrinsic diffusion dynamics, is considered.
  • -diff_2d If diffusions of atoms on a surface or near a surface of a surface slab shall be calculated, where no relevant movements along the z-axis are possible, set this option. Then, the diffusion coefficient formula is adapted to represent a two-dimensional system explicitly.
  • -dens_bins=[number] How many data points shall be written if 1D element densities are called with -dens_elems. The density profiles are written to file dens_elems.dat. This parameter varies the number of bins along the z-axis, for which the density is calculated. Default is 501.
  • -dens_mode=[string] Chooses the axis or sphere along which the element densities shall be evaluated with -dens_elems. One can choose a coordinate axis with x, y or z, or the sphere around the center of mass of one or several elements, with sphere:[list of elements]. In the latter case, the element densities are evaluated similar to RDF profiles in the spheres around a given atom, normalized to the sphere areas. In this case, however, no atom is the central point but the center of mass of one or several of the elements in the system. This can be used, e.g., to describe the averaged or time-dependent (with -part_number) structure of a metal nanoparticle being in contact with a gaseous or liquid surround. For example: -dens_mode=elements:Pt,Cu calculates the densities of all elements around the center of mass of all Pt and Cu atoms. If a coordinate axis is chosen, the concentrations of all elements beneath the surface of a slab are given, if a slab exists with a vacuum region.
  • -dens_depth=[value] If -dens_elems is activated and -dens_mode is set to x, y or z, element concentrations beneath the surface of a slab (if existent) are calculated. This option determines how deep beneath the slab surface the integration area shall range. Default value is 5 Angstroms.
  • -dft_slices=[number] For the DFT element single point calculations, the number of slices along the z-axis, from which a structure with at least one atom of the respective element shall be picked. In practice, the first frame fulfilling the condition of an atom of the chosen element in the frame is picked. If no frame with an atom in the respective area is found in the whole trajectory, no folder is generated for this frame. The default number is 100.
  • -dft_rounds=[number] As explained above, the first frame with an element in the respective slice is chosen as example screenshot for the respective slice. In order to get a better averaging over the whole simulation time, the trajectory (or rather each part of it, after division according to -part_number) can be subdivided into a number of sub-parts, for each part the same analysis is repeated (first frame with an atom in the slice etc.). The default value is 1.