Qassandra

General information

Qassandra (QuAntum corrections to claSSical pArameters for gas electroN DiffRAction) is an open-source software for processing the molecular dynamics (MD) or path-integral molecular dynamics (PIMD) simulations results. It was designed for the implementation of a computational method that allows to make corrections for absence of quantum effects in GED parameters obtained with classical MD and MC simulations.

The preceding program is MDVibCor created by A.Zakharov that could process the MD and PIMD trajectories.

Usage

Qassandra is a command-line application. It is executed by the following command:

<PATH TO THE EXECUTABLE FILE> [OPTIONS] [OPTIONS+INPUT FILES] [OPTIONS+VALUES] > [OUTPUT FILE] .

The side results (all interatomic distances at each frame of simulation or chosen coordinates values during the simulation) are dumped into a file that has name "[name of trajectory file]_atpairs.dat". The output for main results is the display. Therefore it should be redirected manually to the desired output file by the "> [OUTPUT FILE]" part of the command.

Options

Options that require input files

-p - MD simulation trajectory file in XYZ format (path).
-e - file containing equilibrium geometries for all desired conformers.  
-l - file containing description of internal coordinates for large amplitude motions (LAM).
-P - PIMD simulation trajectory files in XYZ format.

If the equilibrium geometry file is not presented the first configuration of MD trajectory will be considered to be equilibrium. In this case the algorithm for quantum correction for MD won't work since atomic masses were not specified. For the PIMD simulations this file should always be given.

Options that set the values

-D - a threshold value (in Å) for the distance of interatomic atomic pairs to be considered equal. Default = 0.00001 Å. 
-t - a time step (in fs) for the MD simulation. Default = 1 fs. 
-T - a MD simulation temperature (in K). Default = 300 K.
-m - a number of steps to be ignored from the start of the MD or PIMD simulation (minimal). Default = 2000.
-M - a maximal number of steps in trajectory to be used. Default = 10000.
-r - number of replicas of PIMD simulation. Default = 4.

The option "-m" deletes the part of the trajectory where the system is equilibrating. The "-t" and "-T" options are crucial for the quantum and thermal corrections calculations.

Options that trigger specific program behavior

-s - the parameters of interatomic terms for every conformer will be symmetrized. 
-S - the parameters of interatomic terms of equal conformers will be symmetrized.
-h - calculate the amplitudes as <math>l = \langle (r - r_e)^2 \rangle^{1/2} </math> (harmonic).
-H - force asymmetry parameters values (κ) to be 0 (harmonic).
-u - turns off algorithm that sorts interatomic terms by their <math>r_e</math> value in the output (unsorted).
-J - "Just do it!".

Input files

In this section the lines that should be in file start with the ">" in the line beginning. Note! The maximal length of each line that Qassandra reads is 100 symbols.

Equilibrium geometry files ("-e" option)

> <Nconf>       <- number of conformers
----------------------------------- <- this section should be given Nconf times
> <Nat>         <- number of atoms
> 17:28, 10 November 2015 (CET)~      <- comment line (it is skipped during the file processing)
*********************************** <- this section should be given Nat times
>  <Atom Name> <AtomMass> <X> <Y> <Z> 
format: [1 or 2 letters] [float number] [float number] [float number] [float number]
or(equivalently as a regexp):
'^ +[A-Z][a-z]? +[0-9]+\.[0-9]+ +-?[0-9]+\.[0-9]+ +-?[0-9]+\.[0-9]+ +-?[0-9]+\.[0-9]+'
***********************************
>                                    <- blank line
-----------------------------------

Atomic mass should be in atomic mass units (a.m.u.), and the Cartesian coordinates in Å.

MD or PIMD trajectory files ("-p" or "-P" option)

The format is basically XYZ file format. It is the default for output of CPMD or CP2K programs, for example.

----------------------------------- <- this section should be given as many times as the trajectory length is 
> <Nat>         <- number of atoms
> 17:28, 10 November 2015 (CET)~      <- comment line (it is skipped during the file processing)
*********************************** <- this section should be given Nat times
>  <Atom Name> <X> <Y> <Z> 
format: [word] [float number] [float number] [float number]
or(equivalently as a regexp):
    '^ +[a-zA-Z]+ +-?[0-9]+\.[0-9]+ +-?[0-9]+\.[0-9]+ +-?[0-9]+\.[0-9]+'
***********************************
-----------------------------------

LAM description file ("-l" option)

> <Nlam>       <- number of torsion angles (currently Qassandra supports only this type of coordinates)
----------------------------------- <- this section should be given Nlam times
>  <atom1>  <atom2>  <atom3>  <atom4>        <- four atoms making the torsion angle
    format: [integer number] [integer number] [integer number] [integer number]   
    or(equivalently as a regexp):
         '^ [0-9]+ +[0-9]+ +[0-9]+ +[0-9]+'            
*********************************** <- this section should be given Nconf times (see equilibrium geometry file)
 >  <min_value>   <max_value>    <- minimal and maximum values possible for the corresponding conformer's torsion angle 
 format: [float number] [float number]
 or(equivalently as a regexp):
    '^ +-?[0-9]+\.[0-9]+ +-?[0-9]+\.[0-9]+'
***********************************
-----------------------------------

Algorithm for quantum and thermal corrections of GED parameters obtained with MD simulations

General considerations

The MD simulations are classical, while atoms in molecules are moving according to quantum mechanics. The lighter the atom the larger is its de Broglie wavelength <math>\lambda = \frac{h}{m v}</math>, the stronger its' quantum behavior. Therefore it was noticed that for the terms containing light elements (such as the ones in the 1st and 2nd row of Periodic table) the GED parameters are badly described by the MD.

MD simulations can also be affected by the "flying ice cube effect". It is the energy drain from high energy degrees of freedom to the low energy ones, such as translations, rotations and LAM. Because of it the high energy motions, such as bond stretch, become frozen leading to even worse agreement with the quantum behavior. This effect appears as a results of some thermostats imperfections.

In order to obtain better description of the GED parameters obtained from MD simulations quantum - thermal corrections should be applied to them.

An illustration of quantum - thermal correction for pair distribution function (C-H term of ethane).