Output plugins are routines which can be activated in DynamO to collect measurements on the simulation while it is running. Although a lot of data can be measured from the configuration files, some properties (e.g., transport coefficients) must be measured using output plugins. There are a wide range of properties which can be measured by DynamO and these are listed below.
All data collected by active output plugins is written into a single XML file at the end of a simulation (output.xml.bz2 by default). Each output plugin has its own section in this file.
DynamO automatically loads the Misc output plugin as it samples a wide range of properties of the system which are computationally inexpensive to collect; however, all other output plugins must be enabled manually. How to enable a plugin is discussed in the following section:
Once you have chosen a plugin to load, you simply add it to the dynarun command line using the -L option. For example, to load the internal energy histogram plugin, you would run
dynarun config.xml -c 1000000 -L IntEnergyHist
And the output.xml.bz2 file generated by this dynarun command will now contain a histogram of the internal energy. Many of the output plugins can take options and can be set using a colon after the plugin name, like so:
dynarun config.xml -c 1000000 -L IntEnergyHist:BinWidth=0.1
If you need to specify multiple options, you must use a comma to delimit options. For example, the mean free time plugin has two options which can be set like so:
dynarun config.xml -c 1000000 -L MFT:BinWidth=0.5,Length=100
In the following sections, the options and output of each output plugin are listed along with details on the calculation of the properties.
Standard plugins are plugins which collect data continuously throughout a simulation. Most properties measured, such as temperature and mean free times, are collected through these types of plugins.
The Misc plugin is the default plugin and is full of a range of properties which are relatively cheap to collect. Each of these properties has a corresponding sub-tag in the output file. The output tags of the Misc plugin and how they are collected are discussed in the following subsections.
There are no options for this plugin.
This tag contains the number of particles divided by the volume of the primary image.
<Misc>
<Density val="0.5"/>
</Misc>
In non-periodic systems, this value may not have any significance as the primary image is not related to the dynamics of the system. Complex boundary effects, such as walls reducing the accessible volume of the system, are not included in this calculation.
This tag contains the excluded volume of all particles divided by the volume of the primary image. The volume of each particle is calculated from the Interaction which is specified by the particle's self-interaction.
<Misc>
<PackingFraction val="0.261799387799154"/>
</Misc>
In non-periodic systems, this value may not have any significance as the primary image is not related to the dynamics of the system. Complex boundary effects, such as walls reducing the accessible volume of the system, are not included in this calculation.
This tag contains the number of Species in the system.
<Misc>
<SpeciesCount val="1"/>
</Misc>
None.
This tag contains the number of Particles in the system.
<Misc>
<ParticleCount val="1372"/>
</Misc>
None.
This tag contains the current and average momentum of all of the particles in the system.
<Misc>
<SystemMomentum>
<Current x="-1.02140518265514e-13" y="3.42226247340705e-14" z="-9.39248678832882e-14"/>
<Average x="-1.02140518265516e-13" y="3.42226247340696e-14" z="-9.39248678832878e-14"/>
</SystemMomentum>
</Misc>
The averages in this tag are collected exactly (see the FAQ on exact averages in DynamO) and so this data is not valid when Lees-Edwards boundary conditions are applied (or currently in systems with gravity, see issue #22).
This tag contains the temperature of the particles in the system. This includes rotational degrees of freedom (if present). As with all temperature values in DynamO, the temperature reported is effectively the product, $k_B\,T$ (see the FAQ on units).
The temperature is calculated using the standard equipartition expression: \[k_B\,T=\frac{2}{N\,f}\sum_i^N E_i^{kinetic}\] where $E_i^{kinetic}$ is the kinetic energy of particle $i$ (including rotational energy) and $f$ is the degrees of freedom of the particle. $f$ is automatically set to $f=3$ in systems without rotational degrees of freedom and to $f=5$ in systems with rotational degrees of freedom.
<Misc>
<Temperature Mean="1.00000000000002" MeanSqr="0.999999999999996" Current="1.00000000000001" Min="1.00000000000001" Max="1.00000000000001"/>
</Misc>
The averages in this tag are collected exactly (see the FAQ on exact averages in DynamO) and so this data is not valid when Lees-Edwards boundary conditions are applied (or currently in systems with gravity, see issue #22).
This tag contains the interaction energy of the system. This is equal to the excess internal energy (also known as the configurational internal energy).
<Misc>
<UConfigurational Mean="-7047.30749520661" MeanSqr="49679058.1379387" Current="-7340" Min="-8232" Max="-6912"/>
</Misc>
The averages in this tag are collected exactly (see the FAQ on exact averages in DynamO) and so this data is not valid when Lees-Edwards boundary conditions are applied (or currently in systems with gravity, see issue #22).
This tag is provided for convenience. It contains the excess heat capacity, calculated from the fluctuations of the configurational internal energy using the following expression: \[\frac{C_v^{ex.}}{k_B}=\frac{\left\langle U_{conf.}^2\right\rangle - \left\langle U_{conf.}\right\rangle^2}{N^2\,k_B^2\,T^2}\] Where the $N^2$ term arises as the UConfigurational values are extensive. To calculate the full heat capacity of the system, the ideal gas contribution must be added: \[\frac{C_v}{k_B}=\frac{C_v^{ex.}}{k_B} +N\frac{f}{2}\] where $f$ is the degrees of freedom of the particle.
<Misc>
<ResidualHeatCapacity Value="80608.2639343955"/>
</Misc>
This expression and the values reported are only valid in the canonical (NVT) ensemble. The restrictions on the UConfigurational plugin also apply.
This tag contains the full pressure tensor for the system, worked out using the virial expression: \[\mathbf{P}=\mathbf{P}^{kinetic}+\mathbf{P}^{interaction}\] where the kinetic pressure is given by: \[\mathbf{P}^{kinetic}=\frac{1}{V}\sum_i^N m_i\,\mathbf{v}_i\mathbf{v}_i\] where $m_i$ is the mass of particle $i$. Please note that $\mathbf{v}_i\mathbf{v}_i$ is a dyadic product which yields a matrix result. The contribution to the pressure due to interactions is given by: \[\mathbf{P}^{interaction}=\frac{1}{V\,t_{sim}}\sum_{ij}^{event} \Delta \mathbf{p}_i\,\mathbf{r}_{ij}\] where $t_{sim}$ is the total simulation time, the summation is over each two-particle event (Interaction), $i$ and $j$ indicate the two particles involved in the event, $\Delta\mathbf{p}_i$ is the momentum impulse on particle $i$, and $\mathbf{r}_{ij}=\mathbf{r}_{i} - \mathbf{r}_{j}$ is the separation vector between the interacting particles. The hydrostatic pressure can be computed from the trace of the pressure tensor: \[p=\mathrm{tr}(\mathbf{P})/3=\left(P_{xx}+P_{yy}+P_{zz}\right)/3\]
<Misc>
<Pressure Avg="-1.92667849178502">
<Tensor>
-1.92776085293681 -0.027284740627551 -0.0330674943624325
-0.027284740627551 -1.91924243731331 -0.0444242381174982
-0.0330674943624326 -0.0444242381174982 -1.93303218510494
</Tensor>
<InteractionContribution>
-2.13900845558922 -0.0269634925340249 -0.032459636868961
-0.0269634925340249 -2.13348059155252 -0.0441115524879792
-0.032459636868961 -0.0441115524879792 -2.14406792617732
</InteractionContribution>
</Pressure>
</Misc>
The kinetic pressure averages are collected exactly (see the FAQ on exact averages in DynamO) and so only the InteractionContribution data is valid when Lees-Edwards boundary conditions are applied or in systems with gravity (see issue #22).
This tag contains information on the duration of the simulation (from the perspective of the model, not the user). See the Timing tag for information on the calculation time of the simulation.
<Misc>
<Duration Events="1059702" OneParticleEvents="0" TwoParticleEvents="1059702" VirtualEvents="1" Time="59.348997249289"/>
</Misc>
None.
This tag contains more detailed information on the types of events which have occurred during the simulation.
<Misc>
<EventCounters>
<Entry Type="Interaction" Name="Bulk" Event="BOUNCE" Count="1636222"/>
<Entry Type="System" Name="Thermostat" Event="GAUSSIAN" Count="136876"/>
<Entry Type="Global" Name="PBCSentinel" Event="VIRTUAL" Count="14"/>
...
</EventCounters>
</Misc>
None.
This tag contains information on the calculation time of the simulation. See the Duration tag for more information on the "simulated" timing.
<Misc>
<Timing Start="Mon Jan 14 20:26:37 2013 " End="Mon Jan 14 20:26:59 2013 " EventsPerSec="49262.070727909" SimTimePerSec="2.75894031681944"/>
</Misc>
None.
This tag is provided for convenience and contains the overall mean free time of a particle. This is calculated from the counts of the one-particle ($N_{1}^{events}$) and two-particle ($N_{2}^{events}$) event counts using the following formula: \[ t_{mft} = t_{sim}\,N/\left(2\,N_{2}^{events} + N_{1}^{events}\right)\]
<Misc>
<totMeanFreeTime val="0.0384196803563759"/>
</Misc>
Virtual events are not counted in this value, therfore it is not a direct measurement of the computational cost of simulating the system.
This tag contains the total number of events executed where the time was earlier than the current time (a negative time until the event occurred). This used to be a good indicator if the simulation was unstable due to accessing an invalid state; however, the introduction of the stable EDMD algorithm has made this tag largely obsolete.
<Misc>
<NegativeTimeEvents Count="0"/>
</Misc>
None.
This tag contains the best estimate DynamO can make of the memory used by DynamO at the last event executed in the simulation. Due to the way this information is obtained, it is the upper limit of the memory in use at that time, and the actual utilisation may be lower.
<Misc>
<Memusage MaxKiloBytes="6580"/>
</Misc>
None.
DynamO contains Einstein correlators for a number of transport coefficients. The definitions used here are the "mainstream" coefficients, which is the most convenient definition for molecular dynamics. \[\begin{align} {\bf X}_q &= - \frac{1}{T}\nabla T & {\bf X}_a &= -T\,\nabla \left(\frac{\mu_a}{T}\right) \end{align} \] \[\begin{align} {\bf J}_a &= \sum_b L_{ab}\,{\bf X}_b + L_{au}\,{\bf X}_u\\ {\bf J}_q &= \sum_a L_{qa}\,{\bf X}_a + L_{qq}\,{\bf X}_q\\ \end{align} \]
The Onsager reciprocity relations are $L_{ab}=L_{ba}$ and $L_{ua}=L_{au}$. As diffusion cannot result in the bulk motion of the fluid, we have $\sum_a {\bf J}_a=0$. As a result, we must have $L_{aa}=\sum_{b\neq a}L_{ab}$.
The phenomenological coefficients $L_{ab}$, $L_{aq}$, and $L_{qq}$ are distinct to the coefficients which appear in the common expressions for the fluxes, i.e., Newton's/Fourier's/Fick's law; However, they are related. As an example, the derivation of the thermal conductivity for a binary mixture is outlined below. Fourier's law is as follows: \[ q = -\lambda \nabla T \] This definition is only applicable in the absence of diffusional forces ($J_a=0$). Assuming a binary system ($a\in1,2$) and using this condition to eliminate the diffusive forces (${\bf X}_a$), we can rearrange the heat flux to read. \[\begin{align} {\bf J_q} &= \left(L_{qq} -\frac{L_{1q}^2}{L_{12}}\right) {\bf X}_q\\ &= \frac{1}{T}\left(L_{qq} -\frac{L_{1q}^2}{L_{12}}\right) \nabla T \end{align} \] Comparing this to Fourier's law allows us to define $\lambda$, but only in the absence of diffusive forces for a binary mixture.
This tag contains the Einstein correlation function which can be used to obtain estimates for the mainstream definition of the thermal conductivity of a fluid, $\mathbf{L}_{qq}$. The correlator used here is: \[\Delta t\,\mathbf{L}_{qq} (\Delta t) = \frac{1}{2\,V\,k_B\,T}\left\langle \left[\phi_{q}(t_0,\Delta t)\right]^2\right\rangle_{t_0} \] where $\Delta t$ is the correlation time, the angle brackets $\langle\rangle_{t_0}$ denote an averaging over the time origin, $t_0$, and $\phi_{q}(t_0,\Delta t)$ is the integral of the microscopic heat flux from $t_0$ to $t_0+\Delta t$. For discontinuous systems, this is most conveniently evaluated using the following expression: \[\phi_{q}(t_0, \Delta t) = \left[\int_{t_0}^{t_0+\Delta t} \sum_i^N\mathbf{v}_i\,e_i\,{\rm d}t\right] + \left[\sum_{i,j}^{events\in[t_0,t_0+\Delta t]}\mathbf{r}_{ij}\,\Delta e_i\right]\] where $e_i$ is the internal energy of particle $i$ (interaction energies should be equally split between the interacting pair). The leftmost term in square brackets contains the free streaming contribution and the sum is over all particles. The rightmost term contains the contribution due to interactions, and the sum is over all two-particle events in the time $[t_0,t_0+\Delta t]$, each involving a particle $i$ and $j$.
You should note that $\mathbf{L}_{qq}$ is a vector quantity, and so three values of the thermal conductivity are measured. These correspond to the transport in the $x$, $y$, and $z$-directions. If your system is isotropic, you may just average the three values to improve your statistics.
As with all of the Green-Kubo/Einstein relationships for the transport coefficients, the desired (hydrodynamic) values are the infinite time correlation values: \[\mathbf{L}_{qq} = \lim_{\Delta t\to\infty} \mathbf{L}_{qq}(\Delta t)\] A convenient way to extrapolate to the long time value is to plot $\Delta t\,L_{qq,x}(\Delta t)$ versus $\Delta t$. The slope of this plot is then the transport coefficient, and a straight line can be regressed to the long time section of the plot to extrapolate the value.
The notes below attempt to explain that a lot of judgement is needed
when evaluation correlation functions. There is only a "window" of
correlation times which can be used to extract the transport
coefficients.
Note 1: In periodic systems, correlation functions
should be studied only up to one sound-wave traversal time of the
periodic image, otherwise additional correlations will appear due to
the "artificial" periodicity of the system.
Note 2: At increasingly long correlation times, the
number of samples of the correlation function obtained from a single
simulation will decrease. Therefore, long time values will have poor
statistics (which should be avoided if possible).
Note 3: At short times, molecular processes such as bond
vibrations or local density fluctuations will strongly influence the
correlation function. For example, take a linear/chain polymer. At
short times, energy might rapidly diffuse up and down the polymer.
This will make it look like there's a high thermal conductivity at
short times; however, the energy will not easily transfer away from
the polymer as this requires the slower inter-polmer interactions to
occur. The energy diffusing up and down the polymer has a high
chance of returning to its starting point. This "return" of the
energy will cause the thermal conductivity to decrease sharply at
slightly longer times. You must therefore use a sufficiently long
correlation time to allow these microscopic processes to average out
so that you can capture the true long-time behaviour of the fluid.
<Misc>
<ThermalConductivity>
<Correlator>
0 0 0 0 0
1 59 3.65611579164476 3.95102401183708 3.10818958034339
2 58 6.66422300757686 10.2522327519284 5.65382677461933
...
</Correlator>
</ThermalConductivity>
</Misc>
This form is only valid in the micro-canonical ensemble (NVE). The current implementation is also only valid for systems without rotational degrees of freedom (such as the hard sphere, square-well, stepped-potential fluids), as it uses the approximation $e_i\approx m\,v^2_i/2+U_{i,config}$ (see issue #23 on the bug tracker).
This tag contains the Einstein correlation function which can be used to obtain estimates for the shear and bulk viscosity of the fluid. Please see the section on ThermalConductivity for most of the notation used here. The correlator used here is: \[\Delta t\,\mathbf{L}_{\eta\eta} (\Delta t) = \frac{1}{2\,V\,k_B\,T}\left(\left\langle \left[\phi_{\eta}(t_0,\Delta t)\right]^2\right\rangle_{t_0} - \left[\Delta t\,V\,\mathbf{P}\right]^2\right) \] where $\phi_{\eta}(t_0, \Delta t)$ is the integral of the microscopic stress flux. For discontinuous systems, this is most conveniently evaluated using the following expression: \[\phi_{\eta}(t_0, \Delta t) = \left[\int_{t_0}^{t_0+\Delta t} \sum_i^Nm_i\,\mathbf{v}_i\,\mathbf{v}_i\,{\rm d}t\right] + \left[\sum_{i,j}^{events\in[t_0,t_0+\Delta t]}\mathbf{r}_{ij}\,m_i\,\Delta \mathbf{v}_i\right]\] The leftmost term in square brackets contains the free streaming contribution and the sum is over all particles. The rightmost term contains the contribution due to interactions, and the sum is over all two-particle events in the time $[t_0,t_0+\Delta t]$, each involving a particle $i$ and $j$.
You should note that $\mathbf{L}_{\eta\eta}$ is a matrix quantity, and so there are 9 measured values. For example, the $L_{\eta\eta,xy}$ element correspons to the transport of $x$-momentum in $y$-direction and vice-versa. In isotropic systems, there are only two important phenomological coefficients related to momentum diffusion: the shear viscosity, $\eta$, and the bulk viscosity, $\kappa$. These are related to $\mathbf{L}_{\eta\eta}$ by the following expressions:
\[\eta = L_{\eta\eta,xy} = L_{\eta\eta,xz} = L_{\eta\eta,yx} = L_{\eta\eta,yz}= L_{\eta\eta,zx} = L_{\eta\eta,zy}\] \[\frac{4}{3}\eta + \kappa = L_{\eta\eta,xx} = L_{\eta\eta,yy} = L_{\eta\eta,zz} \]In isotropic systems we can then average each of the above direction to get an
See the notes in the ThermalConductivity section.
See the notes in the ThermalConductivity section.
<Misc>
<Viscosity>
<Correlator>
0 0 0 0 0 0 0 0 0 0 0
1.20049009599756 116 8.37184016842046 1.62337026137236 2.1589088807082 1.62337026137236 8.42188561600556 1.5159625309261 2.1589088807082 1.5159625309261 7.01773169575511
2.40098019199512 115 20.6408557061605 3.07915750482858 5.12305796346685 3.07915750482858 20.121862094858 3.29182185206649 5.12305796346685 3.2918218520665 16.8444542050959
...
</Correlator>
</Viscosity>
</Misc>
This correlator is valid in all molecular systems.
This tag contains the Einstein correlation function which can be used to obtain estimates for the mutual diffusion coefficients of the fluid. Please see the previous sections for most of the notation used here. The correlator used here is: \[\Delta t\,\mathbf{L}_{ab} (\Delta t) = \frac{1}{2\,V\,k_B\,T}\left\langle \phi_{a}(t_0,\Delta t)\,\phi_{b}(t_0,\Delta t)\right\rangle_{t_0} \] where $\phi_{a}(t_0, \Delta t)$ is the integral of the microscopic flux of Species $a$. For discontinuous systems, this is most conveniently evaluated using the following expression: \[\phi_{a}(t_0, \Delta t) = \int_{t_0}^{t_0+\Delta t} \left(\sum_{i\in a} m_i\,\mathbf{v}_i - c_a \sum_i^N m_i\,\mathbf{v}_i\right){\rm d}t \] The leftmost sum is only over the particles which belong to Species $a$ and the rightmost sum is over all particles. The variable $c_a$ is the mass fraction of the species given by \[c_a = \left(\sum_{i\in a} m_i\right) / \left(\sum_i^N m_i\right)\]
You should note that $\mathbf{L}_{ab}$ is a vector quantity, and so there are 3 measured values. These correspond to the transport in the $x$, $y$, and $z$-directions. If your system is isotropic, you may just average the three values to improve your statistics. For each pairing of Species in the system, there is a corresponding mutual diffusion coefficient to measure the diffusion of one Species through the other. You should note that $L_{ab}=L_{ba}$, so DynamO only collects one half of all of these pairings.
This correlator is valid in all molecular systems.
See the notes in the ThermalConductivity section.
See the notes in the ThermalConductivity section.
<Misc>
<MutualDiffusion>
<Correlator Species1="Bulk" Species2="Bulk">
0 0 0 0 0
1.20049009599756 116 4.34711745476417e-29 1.04569539768932e-28 1.01913292670329e-28
2.40098019199512 115 1.73926903951606e-28 4.20033393703228e-28 4.07395327680373e-28
...
</Correlator>
</MutualDiffusion>
</Misc>
This correlator is valid in all molecular systems.
This tag contains the Einstein correlation function which can be used to obtain estimates for the thermal diffusion coefficients of the fluid. Please see the previous sections for most of the notation used here. The correlator used here is: \[\Delta t\,\mathbf{L}_{aq} (\Delta t) = \frac{1}{2\,V\,k_B\,T}\left\langle \phi_{a}(t_0,\Delta t)\,\phi_{q}(t_0,\Delta t)\right\rangle_{t_0} \] where $\phi_{a}(t_0, \Delta t)$ is the integral of the microscopic flux of Species $a$ and $\phi_{q}(t_0, \Delta t)$ is the integral of the microscopic heat flux. These are evaluated using the expressions given in the previous sections on ThermalConductivity and MutualDiffusion.
You should note that $\mathbf{L}_{aq}$ is a vector quantity, and so there are 3 measured values. These correspond to the transport in the $x$, $y$, and $z$-directions. If your system is isotropic, you may just average the three values to improve your statistics. For each Species in the system, there is a corresponding thermal diffusion coefficient. You should note that $L_{aq}=L_{qa}$.
See the notes in the ThermalConductivity section.
See the notes in the ThermalConductivity section.
<Misc>
<ThermalDiffusion>
<Correlator Species="Bulk">
0 0 0 0 0
1.20049009599756 116 -3.8143162695225e-16 -1.34206107826795e-15 -1.47527366687156e-15
2.40098019199512 115 -1.24651298193849e-15 -5.21423877223736e-15 -5.52897720662611e-15
...
</Correlator>
</ThermalDiffusion>
</Misc>
All restrictions to the ThermalConductivity and MutualDiffusion correlators apply here.
The internal energy histogram plugin collects the exact histogram of the time the system spent in each accessible excess internal energy (see the UConfigurational plugin for more information on this value).
-L IntEnergyHist:BinWidth=0.5
There are no options for this plugin.
<EnergyHist BinWidth="1">
<HistogramWeighted TotalWeight="304.576171140717" Dimension="1" BinWidth="1" AverageVal="-7579.03906905185">
-8231 5.45152795515402e-05
-8230 4.7338829008711e-06
-8229 2.8123254892883e-06
-8228 2.01957587904452e-06
-8227 2.49918411825278e-08
...
</HistogramWeighted>
</EnergyHist>
The MSD plugin calculates the Mean Square Displacement (MSD) of different Species and Structures from their initial locations at the start of the simulation. This can be used to quickly estimate the diffusion coefficients of these objects. The MSDCorrelator plugin collects the evolution of the MSD over time if more detail is needed.
-L MSD
There are no options for this plugin.
<MSD>
<Species Name="Bulk" val="218.539491598976" diffusionCoeff="0.250540360762613"/>
</MSD>
The Trajectory plugin records the full sequence of processed events to a log file (trajectory.out). This plugin is used to debug errors in the detection, scheduling, and execution of events.
-L Trajectory
There are no options for this plugin.
Ticker-type plugins only measure properties at specified intervals of simulation time. An example of such a property is the RadialDistribution plugin.
Each ticker type plugin will take a measurement at the very start of the simulation and then at an interval $\Delta t_{ticker}$. By default, $\Delta t_{ticker}$ is set to the mean free time of the last simulation that the configuration file was run in or a value of 1 if the configuration file has not been run before.
You can set the ticker time by hand (if you are using a ticker-type plugin) using the -t option of dynarun/dynavis:
dynarun -L RadialDistribution -t 0.131 -c 1000000 config.out.xml.bz2
And the Radial Distribution will be sampled every 0.131 units of simulation time until the simulation is stopped after $10^6$ events. Please note that all ticker-type plugins are collected at the same time.
Ticker-type plugins can be used to collect statistics on a configuration by simply running a simulation of zero time. E.g.
dynarun -L RadialDistribution -c 0 config.out.xml.bz2
This will initialise a simulation, run the collection of the output plugin once and halt immediately.
The RadialDistribution plugin calculates the radial distribution function for all pairings of Species in the system. During a simulation, a histogram of the distance between two species of particles is collected. At the end of the simulation, this is converted into the radial distribution using the following expression:
\[g_{\alpha,\beta}(r)=\frac{V}{N_{ticks}\,N_\alpha\,\left(N_\beta-\delta_{\alpha,\beta}\right)\,V_{shell}(r)}\mathcal{H}_{\alpha,\beta}(r)\]where $\mathcal{H}_{\alpha,\beta}(r)$ is the histogram/count of particles of species $\beta$ at a distance of $r$ around particles of species $\alpha$, $V_{shell}(r)=\pi\left(4\,\Delta r\,r^2 + \Delta r^3 / 3\right)$ is the volume of the spherical shell of thickness $\Delta r$ centered about a distance of $r$, $\delta_{\alpha\beta}$ is the kronecker delta function, and $N_{ticks}$ is the number of samples/ticks used to collect the histogram.
-L RadialDistribution:BinWidth=0.01,Length=5
<RadialDistribution SampleCount="29">
<Species Name1="Bulk" Name2="Bulk">
0.10000000000000001 0
0.20000000000000001 0
...
</Species>
</RadialDistribution>
At discontinuities in the potential, the $g(r)$ function is also discontinuous. If you want to correctly capture these discontinuities, you must calculate the values of $g(r)$ either side of the discontinuity using event rates.
The VACF plugin calculates the Velocity AutoCorrelation Function (VACF) for all Species in the system. The velocity autocorrelation function is the following function
\[\Psi_\alpha(\Delta t)=\left\langle\mathbf{v}(t)\cdot\mathbf{v}(t+\Delta t)\right\rangle_{t,\alpha}\]As this is implemented as a ticker property, we only evaluate the VACF at discrete points in time. The implementation is as follows:
\[\begin{align} \Psi_{\alpha}(i\Delta t)&=\left\langle\mathbf{v}(t)\cdot\mathbf{v}(t+i\Delta t)\right\rangle_{t,\alpha} \\ &=N_{ticks}^{-1}\sum_{j}^{N_{ticks}}N_\alpha^{-1}\sum_k^\alpha \mathbf{v}_k(t_j)\cdot\mathbf{v}_k(t_{j+i}) \end{align}\]where $\Delta t$ is the tick interval/time, $N_{ticks}$ is the number of samples/ticks used to collect the correlator, $N_\alpha$ is the number of particles of species $\alpha$, $v_k(t_j)$ is the velocity of the $k^{th}$ particle of species $\alpha$ at the time of the $j^{th}$ tick/sample. The subscript $t,\alpha$ implies the average is over origin times $t$ and all particles of species $\alpha$. If a Topology is defined in the configuration file, the VACF of this will also be calculated and here $\alpha$ will correspond to the Topology, and $N_\alpha$ is the number of Structures in the Topology.
The VACF is usually of interest as its integral is directly related to the diffusion coefficient:
\[d\,D_\alpha = \int_0^\infty \Psi_{\alpha}(t)\,{\rm d}t\]where $d=3$ is the dimensionality of the system and $D_\alpha$ is the diffusion coefficient of species $\alpha$.
-L VACF:Length=100
<VACF>
<Particles>
<Species Name="Bulk">
0 0
0.13112226096298976 1.5495530734986163
0.26224452192597952 0.8305040521045054
...
</Species>
</Particles>
<Topology>
<Structure Name="Ring">
0 2.9826740342384427e-30
0.13112226096298976 2.9760224150436483e-30
...
</Structure>
</Topology>
</VACF>
The MSDCorrelator plugin calculates the Mean Square Displacement (MSD) as a function of time for all Species in the system. The MSD Correlator is the following function
\[\Psi(\Delta t) = \left\langle \left[{\bf r}(t+\Delta t)-{\bf r}(t)\right]^2\right\rangle_{t,\alpha}\]As this is implemented as a ticker property, we only evaluate the MSD at discrete points in time. The implementation is as follows:
\[\begin{align} \Psi_{\alpha}(i\Delta t)&=\left\langle \left[{\bf r}(t+i\Delta t)-{\bf r}(t)\right]^2\right\rangle_{t,\alpha} \\ &=N_{ticks}^{-1}\sum_{j}^{N_{ticks}}N_\alpha^{-1}\sum_k^\alpha \left[{\bf r}_k(t_j)-{\bf r}_k(t_{j+i})\right]^2 \end{align}\]where $\Delta t$ is the tick interval/time, $N_{ticks}$ is the number of samples/ticks used to collect the correlator, $N_\alpha$ is the number of particles of species $\alpha$, ${\bf r}_k(t_j)$ is the position of the $k^{th}$ particle of species $\alpha$ at the time of the $j^{th}$ tick/sample. The subscript $t,\alpha$ implies the average is over origin times $t$ and all particles of species $\alpha$. If a Topology is defined in the configuration file, the MSD of this will also be calculated and here $\alpha$ will correspond to the Topology, and $N_\alpha$ is the number of Structures in the Topology.
The MSD is usually of interest as its gradient is directly related to the diffusion coefficient:
\[D_\alpha=\lim_{\Delta t\to\infty}\frac{\left\langle \left[{\bf r}(t+\Delta t)-{\bf r}(t)\right]^2\right\rangle_{t,\alpha}}{2\,d\,t}\]where $d=3$ is the dimensionality of the system and $D_\alpha$ is the diffusion coefficient of species $\alpha$. The MSD plugin will calculate a prediction of the diffusion coefficient using the full time of the simulation.
-L MSDCorrelator:Length=100
<MSDCorrelator>
<Particles>
<Species Name="Bulk">
0 0
0.13112226096298976 1.5495530734986163
0.26224452192597952 0.8305040521045054
...
</Species>
</Particles>
<Topology>
<Structure Name="Ring">
0 2.9826740342384427e-30
0.13112226096298976 2.9760224150436483e-30
...
</Structure>
</Topology>
</MSDCorrelator>
The RadiusGyration plugin calculates and averages the radius of gyration for each Structure tag in the system. This is used as a measure of the size and shape of complex molecules such as polymers.
To calculate the radii of gyration, we attempt to account for periodic boundary conditions by walking along the chain using the previous particle's position as the origin of the coordinate system before applying boundary conditions. We first calculate the inertia tensor like so: \[ \mathbf{I} = \sum_i m_i\left[\left(\mathbf{r}_i\cdot\mathbf{r}_i\right)\mathbf{1} - \mathbf{r}_i\,\mathbf{r}_i\right] \] Then we scale out the mass and calculate the eigenvectors and values of the matrix $\mathbf{I} / \sum_i m_i$. The unit eigenvectors are the principle axes of the molecule and the eigenvalues are the radii of gyration.
This plugin also calculates the Nematic order parameter. For each molecule $i$, we assign the principle axis corresponding to the largest radius of gyration as the director $\mathbf{u}_i$. Then we construct the following tensor: \[ \mathbf{S} = \frac{1}{2\,N_{mono}}\sum_i \left(3\,\mathbf{u}_i\,\mathbf{u}_i - \mathbf{1}\right) \] We then calculate the eigenvectors and values of $\mathbf{S}$. If the system is aligned nematically, the eigenvectors indicate the alignment directions of the molecules. In particular, the largest eigenvector should correspond to the fluid director.
Assume we rotate the matrix $\mathbf{S}$ to a coordinate frame where $x$, $y$, and $z$ coincide with the eigenvectors (in increasing order). In perfectly aligned nematic, prolate fluids, we expect: \[ \tilde{\mathbf{S}} = \begin{bmatrix}-0.5 & 0 & 0\\0 & -0.5 & 0 \\ 0 & 0 & 1\end{bmatrix} \] In perfectly aligned oblate geometry we expect \[ \tilde{\mathbf{S}} = \begin{bmatrix}0.5 & 0 & 0\\0 & 0.5 & 0 \\ 0 & 0 & -1\end{bmatrix} \] If the molecules are not aligned/isotropic, we expect $\tilde{\mathbf{S}}=\mathbf{0}$.
-L RadiusGyration:BinWidthGyration=0.01,BinWidthNematic=0.001
<ChainGyration>
<Chain Name="HelixPolymer">
<GyrationRadii>
<Histogram SampleCount="1578" Dimension="1" BinWidth="0.01" AverageVal="0.42709125475285176">
0.029999999999999999 0.063371356147021538
0.070000000000000007 0.6337135614702154
...
</Histogram>
<Histogram SampleCount="1578" Dimension="1" BinWidth="0.01" AverageVal="1.6703738910012673">
0.28999999999999998 0.063371356147021538
0.32000000000000001 0.12674271229404308
...
</Histogram>
<Histogram SampleCount="1578" Dimension="1" BinWidth="0.01" AverageVal="9.7681432192648927">
2.1200000000000001 0.063371356147021538
2.2200000000000002 0.063371356147021538
...
</Histogram>
</GyrationRadii>
<NematicOrderParameter x="-0.5" y="-0.49999999999999994" z="0.99999999999999911">
<Histogram SampleCount="1578" Dimension="1" BinWidth="0.001" AverageVal="-0.5">
...
</Histogram>
<Histogram SampleCount="1578" Dimension="1" BinWidth="0.001" AverageVal="-0.5">
...
</Histogram>
<Histogram SampleCount="1578" Dimension="1" BinWidth="0.001" AverageVal="1">
...
</Histogram>
</NematicOrderParameter>
</Chain>
</ChainGyration>
Page last modified: Sunday 9th June 2024