# Qassandra

## General information

Qassandra (QuAntum corrections to claSSical pArameters for gas electroN DiffRAction) is an open-source software for processing results from molecular dynamics (MD) and path-integral molecular dynamics (PIMD) simulations. In contrast to the available program MDVibCor by A. Zakharov, Qassandra additionally implements a method for calculation of quantum corrections for vibrational molecular parameters obtained using classical MD and MC simulations.<bib id="VishnevskiyTCA2016" />

## 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] .

Main results are printed to stdout and should normally be redirected (> [OUTPUT FILE]) to a file. Additional information, like interatomic distances for each frame of simulation or values of selected coordinates during the simulation, are dumped into a file with name [name of trajectory file]_atpairs.dat.

## 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.
-V - MC simulation trajectory file in modified XYZ format (Viva Las Vegas!).


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 $l = \langle (r - r_e)^2 \rangle^{1/2}$ (harmonic).
-H - force asymmetry parameters values (κ) to be 0 (harmonic).
-u - turns off algorithm that sorts interatomic terms by their $r_e$ value in the output (unsorted).
-U  - print asymmetry parameters in UNEX style for "SMSModel=1" (6.0*κ).
-g  - print $r_e - r_g$ corrections instead of $r_e - r_a$.
-R  - include rotation-vibration interaction effects during MC simulations processing;
-b  - print plain bibliography for this program.
-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]+'
***********************************
-----------------------------------


### MC trajectory files ("-V" option)

The format is basically slightly modified XYZ file format.

----------------------------------- <- this section should be given as many times as the trajectory length is
> <Nat>         <- number of atoms
>  <Energy of the conformation>   [<additional weight of the configuration>]  <- second number should be given if option "-W" is used
format: [float number] [float number]
or(equivalently as a regexp):
'^ -?[0-9]+\.[0-9]+ [0-9]+\.[0-9]'
*********************************** <- 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 $\lambda = \frac{h}{m v}$, 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).