3. yaff.sampling
– Phasespace sampling¶
3.1. yaff.sampling.dof
– Abstraction layer for degrees of freedom¶
Abstraction layer for degrees of freedom
All these classes are called DOF classes, because they specify a set of degrees of freedom. These DOF classes are used for geometry/cell optimization and harmonic approximations.

class
yaff.sampling.dof.
DOF
(ff)¶ Bases:
object
Arguments:
 ff
 A force field object

ndof
¶

reset
()¶

check_delta
(x=None, eps=0.0001, zero=None)¶ Test the analytical derivatives

log
()¶

class
yaff.sampling.dof.
CartesianDOF
(ff, gpos_rms=1e05, dpos_rms=0.001, select=None)¶ Bases:
yaff.sampling.dof.DOF
Cartesian degrees of freedom
This DOF is also applicable to periodic systems. Cell parameters are not modified when this DOF is used.
Arguments:
 ff
 A force field object.
Optional arguments:
 gpos_rms, dpos_rms
Thresholds that define the convergence. If all of the actual values drop below these thresholds, the minimizer stops.
For each rms threshold, a corresponding max threshold is included automatically. The maximum of the absolute value of a component should be smaller than 3/sqrt(N) times the rms threshold, where N is the number of degrees of freedom.
 select
 A selection of atoms whose degrees of freedom are included. If not list is given, all atomic coordinates are included.
Convergence conditions:
 gpos_rms
 The rootmeansquare of the norm of the gradients of the atoms.
 dpos_rms
 The rootmeansquare of the norm of the displacements of the atoms.

fun
(x, do_gradient=False)¶ Computes the energy and optionally the gradient.
Arguments:
 x
 The degrees of freedom
Optional arguments:
 do_gradient
 When True, the gradient is also returned.

check_convergence
()¶

class
yaff.sampling.dof.
BaseCellDOF
(ff, gpos_rms=1e05, dpos_rms=0.001, grvecs_rms=1e05, drvecs_rms=0.001, do_frozen=False, freemask=None)¶ Bases:
yaff.sampling.dof.DOF
Fractional coordinates and cell parameters
Several subclasses of BaseCellDOF are implemented below. Each one considers a specific representation and subset of the cell parameters.
The following variable names are consistently used (also in subclasses):
 cellvars
 An array with all variables for the cell (specific for ja BaseCellDOF subclass).
 ncellvar
 The number of cellvars (at most 9).
 celldofs
 A selection of the elements in cellvars, based on freemask.
 ncelldof
 The number of celldofs (less than or equal to ncellvar).
 frac
 Fractional coordinates.
 x
 All degrees of freedom, i.e. celldofs and frac (in that order, frac is optional).
 The suffix 0
 Used for initial values of something.
Arguments:
 ff
 A force field object.
Optional arguments:
 gpos_rms, dpos_rms, grvecs_rms, drvecs_rms
Thresholds that define the convergence. If all of the actual values drop below these thresholds, the minimizer stops.
For each rms threshold, a corresponding max threshold is included automatically. The maximum of the absolute value of a component should be smaller than 3/sqrt(N) times the rms threshold, where N is the number of degrees of freedom.
 do_frozen
 When True, the fractional coordinates of the atoms are kept fixed.
 freemask
 When given, this must be an array of booleans indicating which cellvars are free. At least one cellvar must be free.
Convergence conditions:
 gpos_rms
 The rootmeansquare of the norm of the gradients of the atoms.
 dpos_rms
 The rootmeansquare of the norm of the displacements of the atoms.
 grvecs_rms
 The rootmeansquare of the norm of the gradients of the cell vectors.
 drvecs_rms
 The rootmeansquare of the norm of the displacements of the cell vectors.

ncellvar
¶ The number of cellvars

ncelldof
¶ The number of celldofs (free cellvars)

fun
(x, do_gradient=False)¶ Computes the energy and optionally the gradient.
Arguments:
 x
 All degrees of freedom
Optional arguments:
 do_gradient
 When True, the gradient is also returned.

check_convergence
()¶

log
()¶

class
yaff.sampling.dof.
FullCellDOF
(ff, gpos_rms=1e05, dpos_rms=0.001, grvecs_rms=1e05, drvecs_rms=0.001, do_frozen=False, freemask=None)¶ Bases:
yaff.sampling.dof.BaseCellDOF
DOF that includes all 9 components of the cell vectors
The degrees of freedom are rescaled cell vectors ordered in one row:
 3D periodic: [a_x/s, a_y/s, a_z/s, b_x/s, b_y/s, b_z/s, c_x/s, c_y/s, c_z/s] where s is the cube root of the initial cell volume such that the cell DOFs become dimensionless.
 2D periodic: [a_x/s, a_y/s, a_z/s, b_x/s, b_y/s, b_z/s] where s is the square root of the initial cell surface such that the cell DOFs become dimensionless.
 1D periodic: [a_x/s, a_y/s, a_z/s] where s is the length of the initial cell vector such that the cell DOFs become dimensionless.
Arguments:
 ff
 A force field object.
Optional arguments:
 gpos_rms, dpos_rms, grvecs_rms, drvecs_rms
Thresholds that define the convergence. If all of the actual values drop below these thresholds, the minimizer stops.
For each rms threshold, a corresponding max threshold is included automatically. The maximum of the absolute value of a component should be smaller than 3/sqrt(N) times the rms threshold, where N is the number of degrees of freedom.
 do_frozen
 When True, the fractional coordinates of the atoms are kept fixed.
 freemask
 When given, this must be an array of booleans indicating which cellvars are free. At least one cellvar must be free.
Convergence conditions:
 gpos_rms
 The rootmeansquare of the norm of the gradients of the atoms.
 dpos_rms
 The rootmeansquare of the norm of the displacements of the atoms.
 grvecs_rms
 The rootmeansquare of the norm of the gradients of the cell vectors.
 drvecs_rms
 The rootmeansquare of the norm of the displacements of the cell vectors.

class
yaff.sampling.dof.
StrainCellDOF
(ff, gpos_rms=1e05, dpos_rms=0.001, grvecs_rms=1e05, drvecs_rms=0.001, do_frozen=False, freemask=None)¶ Bases:
yaff.sampling.dof.BaseCellDOF
Eliminates rotations of the unit cell. thus six cell parameters are free.
The degrees of freedom are coefficients in symmetrix matrix transformation, A, that is applied to the initial cell vectors.
 3D periodic: [A_00, A_11, A_22, 2*A_12, 2*A_20, 2*A_01]
 2D periodic: [A_00, A_11, 2*A_01]
 1D periodic: [A_00]
Why does this work? Let R be the array with cell vectors as rows. It can always be written as a product,
R = R_0.F,where F is an arbitrary 3x3 matrix. Application of SVD to the matrix F :Yields: R = R_0.U.S.V^T = R_0.U.V^T.V.S.V^T
Then W=U.V^T is a orthonormal matrix and A=V.S.V^T is a symmetric matrix. The orthonormal matrix W is merely a rotation of the cell vectors, which can be omitted as the internal energy is invariant under such rotations. The symmetric matrix actually deforms the cell and is the part of interest.
Arguments:
 ff
 A force field object.
Optional arguments:
 gpos_rms, dpos_rms, grvecs_rms, drvecs_rms
Thresholds that define the convergence. If all of the actual values drop below these thresholds, the minimizer stops.
For each rms threshold, a corresponding max threshold is included automatically. The maximum of the absolute value of a component should be smaller than 3/sqrt(N) times the rms threshold, where N is the number of degrees of freedom.
 do_frozen
 When True, the fractional coordinates of the atoms are kept fixed.
 freemask
 When given, this must be an array of booleans indicating which cellvars are free. At least one cellvar must be free.
Convergence conditions:
 gpos_rms
 The rootmeansquare of the norm of the gradients of the atoms.
 dpos_rms
 The rootmeansquare of the norm of the displacements of the atoms.
 grvecs_rms
 The rootmeansquare of the norm of the gradients of the cell vectors.
 drvecs_rms
 The rootmeansquare of the norm of the displacements of the cell vectors.

class
yaff.sampling.dof.
AnisoCellDOF
(ff, gpos_rms=1e05, dpos_rms=0.001, grvecs_rms=1e05, drvecs_rms=0.001, do_frozen=False, freemask=None)¶ Bases:
yaff.sampling.dof.BaseCellDOF
Only the lengths of the cell vectors are free. angles are fixed.
The degrees of freedom are dimensionless scale factors for the cell lengths, using the initial cell vectors as the reference point. (This is one DOF per periodic dimension.)
Arguments:
 ff
 A force field object.
Optional arguments:
 gpos_rms, dpos_rms, grvecs_rms, drvecs_rms
Thresholds that define the convergence. If all of the actual values drop below these thresholds, the minimizer stops.
For each rms threshold, a corresponding max threshold is included automatically. The maximum of the absolute value of a component should be smaller than 3/sqrt(N) times the rms threshold, where N is the number of degrees of freedom.
 do_frozen
 When True, the fractional coordinates of the atoms are kept fixed.
 freemask
 When given, this must be an array of booleans indicating which cellvars are free. At least one cellvar must be free.
Convergence conditions:
 gpos_rms
 The rootmeansquare of the norm of the gradients of the atoms.
 dpos_rms
 The rootmeansquare of the norm of the displacements of the atoms.
 grvecs_rms
 The rootmeansquare of the norm of the gradients of the cell vectors.
 drvecs_rms
 The rootmeansquare of the norm of the displacements of the cell vectors.

class
yaff.sampling.dof.
IsoCellDOF
(ff, gpos_rms=1e05, dpos_rms=0.001, grvecs_rms=1e05, drvecs_rms=0.001, do_frozen=False, freemask=None)¶ Bases:
yaff.sampling.dof.BaseCellDOF
The cell is only allowed to undergo isotropic scaling
The only degree of freedom is an isotropic scaling factor, using the initial cell vectors as a reference.
Arguments:
 ff
 A force field object.
Optional arguments:
 gpos_rms, dpos_rms, grvecs_rms, drvecs_rms
Thresholds that define the convergence. If all of the actual values drop below these thresholds, the minimizer stops.
For each rms threshold, a corresponding max threshold is included automatically. The maximum of the absolute value of a component should be smaller than 3/sqrt(N) times the rms threshold, where N is the number of degrees of freedom.
 do_frozen
 When True, the fractional coordinates of the atoms are kept fixed.
 freemask
 When given, this must be an array of booleans indicating which cellvars are free. At least one cellvar must be free.
Convergence conditions:
 gpos_rms
 The rootmeansquare of the norm of the gradients of the atoms.
 dpos_rms
 The rootmeansquare of the norm of the displacements of the atoms.
 grvecs_rms
 The rootmeansquare of the norm of the gradients of the cell vectors.
 drvecs_rms
 The rootmeansquare of the norm of the displacements of the cell vectors.

class
yaff.sampling.dof.
FixedBCDOF
(ff, gpos_rms=1e05, dpos_rms=0.001, grvecs_rms=1e05, drvecs_rms=0.001, do_frozen=False, freemask=None)¶ Bases:
yaff.sampling.dof.BaseCellDOF
A rectangular cell that can only stretch along one axis
This cell optimization constrains the cell in the y and z direction to the original values, but allows expansion and contraction in the x direction. The system should be rotated such that the initial cell vectors look like:
a = ( ax , 0 , 0 ) b = ( 0 , by , bz ) c = ( 0 , cy , cz )
During optimization, only ax will be allowed to change.
This type of constraint can be used when looking at a structure that is periodic only in one dimension, but you have to fake a 3D structure to be able to use Ewald summation
Arguments:
 ff
 A force field object.
Optional arguments:
 gpos_rms, dpos_rms, grvecs_rms, drvecs_rms
Thresholds that define the convergence. If all of the actual values drop below these thresholds, the minimizer stops.
For each rms threshold, a corresponding max threshold is included automatically. The maximum of the absolute value of a component should be smaller than 3/sqrt(N) times the rms threshold, where N is the number of degrees of freedom.
 do_frozen
 When True, the fractional coordinates of the atoms are kept fixed.
 freemask
 When given, this must be an array of booleans indicating which cellvars are free. At least one cellvar must be free.
Convergence conditions:
 gpos_rms
 The rootmeansquare of the norm of the gradients of the atoms.
 dpos_rms
 The rootmeansquare of the norm of the displacements of the atoms.
 grvecs_rms
 The rootmeansquare of the norm of the gradients of the cell vectors.
 drvecs_rms
 The rootmeansquare of the norm of the displacements of the cell vectors.

class
yaff.sampling.dof.
FixedVolOrthoCellDOF
(ff, volume=None, gpos_rms=1e05, dpos_rms=0.001, grvecs_rms=1e05, drvecs_rms=0.001, do_frozen=False, freemask=None)¶ Bases:
yaff.sampling.dof.BaseCellDOF
Orthorombic cell optimizer with a fixed volume.
These constraints are implemented by using the following cell vectors:
a = ( s*a0*la , 0 , 0 ) b = ( 0 , s*b0*lb , 0 ) c = ( 0 , 0 , s*c0/(la*lb) )
with s = (V/V0)^(1/3)
Optional arguments (in addition to those of BaseCellDOF):
 volume
 The desired volume of the cell. (When not given, the current volume of the system is not altered.)
3.2. yaff.sampling.harmonic
– Harmonic models¶
Harmonic models

yaff.sampling.harmonic.
estimate_hessian
(dof, eps=0.0001)¶ Estimate the Hessian using the symmetric finite difference approximation.
Arguments:
 dof
 A DOF object
Optional arguments:
 eps
 The magnitude of the displacements

yaff.sampling.harmonic.
estimate_cart_hessian
(ff, eps=0.0001, select=None)¶ Estimate the Cartesian Hessian with symmetric finite differences.
Arguments:
 ff
 A force field object
Optional arguments:
 eps
 The magnitude of the Cartesian displacements
 select
 A selection of atoms for which the hessian must be computed. If not given, the entire hessian is computed.

yaff.sampling.harmonic.
estimate_elastic
(ff, eps=0.0001, do_frozen=False, ridge=0.0001)¶ Estimate the elastic constants using the symmetric finite difference approximation.
Arguments:
 ff
 A force field object
Optional arguments:
 eps
 The magnitude of the Cartesian displacements
 do_frozen
 By default this is False, which means that the changes in fractional atomic coordinates due to cell deformations are properly taken into account. When this is set to True, such displacements (other than uniform scaling) are ignored. The latter is much faster, but only correct for the simplest materials.
 ridge
 Threshold for the eigenvalues of the Cartesian Hessian. This only
matters if
do_frozen==False
.
The elastic constants are second order derivatives of the strain energy density with respect to uniform deformations. At the molecular scale, uniform deformations can be describe by a linear transformation of the atoms and the cell parameters:
x' = e . x0
When x0 corresponds to the reference point for the second order expansion, small deviations from that reference point correspond to small deviations of e from the unit matrix. The strain energy density is nothing but the energy of the system divided by its volume, minus a the energy density of the fully relaxed system. The adjective strain refers to deviations from the relaxed reference. When computing the second order derivatives, this reference point can be ignored. Assuming the reference point is the relaxed system, the strain energy density can be approximated to second order as:
u = 1/2 \sum_ijkl e_ij c_ijkl e_kl
In this equation the matrix c_ijkl contains all the elastic constants. In principle, it contains 81 values, but due to symmetry considerations, only 21 of these values are independent. [The symmetry considerations are as follows: (i) it is sufficient to consider symmetric deformations e when excluding cell rotations, (ii) the index pairs (ij) and (kl) are interchangeable.] These 21 numbers are typically represented using the (inclusive) lower diagonal of a 6x6 matrix, using the following convention to map pair indexes (ij) and (kl) onto single indexes n and m:
(ij) or (kl) n or m 11 1 22 2 22 3 23 4 31 5 12 6 (This table is also known as the compressed Voight notation.)
This routine returns the elastic constants in a symmetric 6x6 matrix, using the index conventions by Voight.
3.3. yaff.sampling.io
– Trajectory writers¶
Trajectory writers

class
yaff.sampling.io.
HDF5Writer
(f, start=0, step=1)¶ Bases:
yaff.sampling.iterative.Hook
Argument:
 f
 A h5.File object to write the trajectory to.
Optional arguments:
 start
 The first iteration at which this hook should be called.
 step
 The hook will be called every step iterations.

__call__
(iterative)¶

dump_system
(system)¶

init_trajectory
(iterative)¶

class
yaff.sampling.io.
XYZWriter
(fn_xyz, select=None, start=0, step=1)¶ Bases:
yaff.sampling.iterative.Hook
Argument:
 fn_xyz
 A filename to write the XYZ trajectory too.
Optional arguments:
 select
 A list of atom indexes that should be written to the trajectory output. If not given, all atoms are included.
 start
 The first iteration at which this hook should be called.
 step
 The hook will be called every step iterations.

__call__
(iterative)¶

class
yaff.sampling.io.
RestartWriter
(f, start=0, step=1000)¶ Bases:
yaff.sampling.iterative.Hook
Argument:
 f
 A h5.File object to write the restart information to.
Optional arguments:
 start
 The first iteration at which this hook should be called.
 step
 The hook will be called every step iterations.

init_state
(iterative)¶

__call__
(iterative)¶

dump_system
(system)¶

init_trajectory
(iterative)¶

dump_restart
(hook)¶
3.4. yaff.sampling.iterative
– Base class for iterative algorithms¶
Base class for iterative algorithms

class
yaff.sampling.iterative.
Iterative
(ff, state=None, hooks=None, counter0=0)¶ Bases:
object
Arguments:
 ff
 The ForceField instance used in the iterative algorithm
Optional arguments:
 state
 A list with state items. State items are simple objects that take or derive a property from the current state of the iterative algorithm.
 hooks
 A function (or a list of functions) that is called after every iterative.
 counter0
 The counter value associated with the initial state.

default_state
= []¶

log_name
= 'ITER'¶

initialize
()¶

call_hooks
()¶

run
(nstep=None)¶

propagate
()¶

finalize
()¶

class
yaff.sampling.iterative.
StateItem
(key)¶ Bases:
object

update
(iterative)¶

get_value
(iterative)¶

iter_attrs
(iterative)¶

copy
()¶


class
yaff.sampling.iterative.
AttributeStateItem
(key)¶ Bases:
yaff.sampling.iterative.StateItem

get_value
(iterative)¶

copy
()¶


class
yaff.sampling.iterative.
PosStateItem
¶ Bases:
yaff.sampling.iterative.StateItem

get_value
(iterative)¶


class
yaff.sampling.iterative.
DipoleStateItem
¶ Bases:
yaff.sampling.iterative.StateItem

get_value
(iterative)¶


class
yaff.sampling.iterative.
DipoleVelStateItem
¶ Bases:
yaff.sampling.iterative.StateItem

get_value
(iterative)¶


class
yaff.sampling.iterative.
VolumeStateItem
¶ Bases:
yaff.sampling.iterative.StateItem

get_value
(iterative)¶


class
yaff.sampling.iterative.
CellStateItem
¶ Bases:
yaff.sampling.iterative.StateItem

get_value
(iterative)¶


class
yaff.sampling.iterative.
EPotContribStateItem
¶ Bases:
yaff.sampling.iterative.StateItem
Keeps track of all the contributions to the potential energy.

get_value
(iterative)¶

iter_attrs
(iterative)¶


class
yaff.sampling.iterative.
EpotBondsStateItem
(do_ei=False)¶ Bases:
yaff.sampling.iterative.StateItem
Keeps track of all the Valence Bond contributions to the potential energy.
Optional arguments:
 do_ei
 If True, the electrostatic contributions of the bonded atom pair will also be tracked.

get_value
(iterative)¶

iter_atts
(iterative)¶

class
yaff.sampling.iterative.
EpotBendsStateItem
(do_ei=False)¶ Bases:
yaff.sampling.iterative.StateItem
Keeps track of all the Valence Bend contributions to the potential energy
Optional arguments:
 do_ei
 If True, the electrostatic contributions of the 13 nonbonded atom pair will also be tracked.

get_value
(iterative)¶

iter_atts
(iterative)¶

class
yaff.sampling.iterative.
EpotDihedsStateItem
(do_ei=False)¶ Bases:
yaff.sampling.iterative.StateItem
Keeps track of all the Valence Dihedral contributions to the potential energy
Optional arguments:
 do_ei
 If True, the electrostatic contributions of the 14 nonbonded atom pair will also be tracked.

get_value
(iterative)¶

iter_atts
(iterative)¶
3.5. yaff.sampling.npt
– Barostats¶
Barostats

class
yaff.sampling.npt.
TBCombination
(thermostat, barostat, start=0)¶ Bases:
yaff.sampling.verlet.VerletHook
VerletHook combining an arbitrary Thermostat and Barostat instance, which ensures these instances are called in the correct succession, and possible coupling between both is handled correctly.
Arguments:
 thermostat
 A Thermostat instance
 barostat
 A Barostat instance

name
= 'TBCombination'¶

init
(iterative)¶

pre
(iterative)¶

post
(iterative)¶

expectscall
(iterative, kind)¶

verify
()¶

class
yaff.sampling.npt.
McDonaldBarostat
(temp, press, start=0, step=1, amp=0.001)¶ Bases:
yaff.sampling.verlet.VerletHook
Warning: this code is not fully tested yet!
Arguments:
 temp
 The average temperature of the NpT ensemble
 press
 The external pressure of the NpT ensemble
Optional arguments:
 start
 The first iteration at which this hook is called
 step
 The number of iterations between two subsequent calls to this hook.
 amp
 The amplitude of the changes in the logarithm of the volume.

name
= 'McDonald'¶

kind
= 'stochastic'¶

method
= 'barostat'¶

init
(iterative)¶

pre
(iterative, chainvel0=None)¶

post
(iterative, chainvel0=None)¶

class
yaff.sampling.npt.
BerendsenBarostat
(ff, temp, press, start=0, step=1, timecon=41341.37333664682, beta=13445.401558486123, anisotropic=True, vol_constraint=False, restart=False)¶ Bases:
yaff.sampling.verlet.VerletHook
This hook implements the Berendsen barostat. The equations are derived in:
Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; Dinola, A.; Haak, J. R. J. Chem. Phys. 1984, 81, 36843690Arguments:
 ff
 A ForceField instance.
 temp
 The temperature of thermostat.
 press
 The applied pressure for the barostat.
Optional arguments:
 start
 The step at which the thermostat becomes active.
 timecon
 The time constant of the Berendsen barostat.
 beta
 The isothermal compressibility, conventionally the compressibility of liquid water
 anisotropic
 Defines whether anisotropic cell fluctuations are allowed.
 vol_constraint
 Defines whether the volume is allowed to fluctuate.

name
= 'Berendsen'¶

kind
= 'deterministic'¶

method
= 'barostat'¶

init
(iterative)¶

pre
(iterative, chainvel0=None)¶

post
(iterative, chainvel0=None)¶

class
yaff.sampling.npt.
LangevinBarostat
(ff, temp, press, start=0, step=1, timecon=41341.37333664682, anisotropic=True, vol_constraint=False)¶ Bases:
yaff.sampling.verlet.VerletHook
This hook implements the Langevin barostat. The equations are derived in:
Feller, S. E.; Zhang, Y.; Pastor, R. W.; Brooks, B. R. J. Chem. Phys. 1995, 103, 46134621Arguments:
 ff
 A ForceField instance.
 temp
 The temperature of thermostat.
 press
 The applied pressure for the barostat.
Optional arguments:
 start
 The step at which the barostat becomes active.
 timecon
 The time constant of the Langevin barostat.
 anisotropic
 Defines whether anisotropic cell fluctuations are allowed.
 vol_constraint
 Defines whether the volume is allowed to fluctuate.

name
= 'Langevin'¶

kind
= 'stochastic'¶

method
= 'barostat'¶

init
(iterative)¶

pre
(iterative, chainvel0=None)¶

post
(iterative, chainvel0=None)¶

baro
(iterative, chainvel0)¶

getR
()¶

class
yaff.sampling.npt.
MTKBarostat
(ff, temp, press, start=0, step=1, timecon=41341.37333664682, anisotropic=True, vol_constraint=False, baro_thermo=None, vel_press0=None, restart=False)¶ Bases:
yaff.sampling.verlet.VerletHook
This hook implements the MartynaTobiasKlein barostat. The equations are derived in:
Martyna, G. J.; Tobias, D. J.; Klein, M. L. J. Chem. Phys. 1994, 101, 41774189.The implementation (used here) of a symplectic integrator of this barostat is discussed in
Martyna, G. J.; Tuckerman, M. E.; Tobias, D. J.; Klein, M. L. Mol. Phys. 1996, 87, 11171157.Arguments:
 ff
 A ForceField instance.
 temp
 The temperature of thermostat.
 press
 The applied pressure for the barostat.
Optional arguments:
 start
 The step at which the barostat becomes active.
 timecon
 The time constant of the MartynaTobiasKlein barostat.
 anisotropic
 Defines whether anisotropic cell fluctuations are allowed.
 vol_constraint
 Defines whether the volume is allowed to fluctuate.
 baro_thermo
 NHCThermostat instance, coupled directly to the barostat
 vel_press0
 The initial barostat velocity tensor
 restart
 If true, the cell is not symmetrized initially

name
= 'MTTK'¶

kind
= 'deterministic'¶

method
= 'barostat'¶

init
(iterative)¶

pre
(iterative, chainvel0=None)¶

post
(iterative, chainvel0=None)¶

baro
(iterative, chainvel0)¶

add_press_cont
()¶

class
yaff.sampling.npt.
PRBarostat
(ff, temp, press, start=0, step=1, timecon=41341.37333664682, anisotropic=True, vol_constraint=False, baro_thermo=None, vel_press0=None, restart=False)¶ Bases:
yaff.sampling.verlet.VerletHook
This hook implements the ParrinelloRahman barostat for finite strains. The equations are derived in:
Parrinello, M.; Rahman, A. J. Appl. Phys. 1981, 52, 71827190.Arguments:
 ff
 A ForceField instance.
 temp
 The temperature of thermostat.
 press
 The applied pressure for the barostat.
Optional arguments:
 start
 The step at which the barostat becomes active.
 timecon
 The time constant of the TadmorMiller barostat.
 anisotropic
 Defines whether anisotropic cell fluctuations are allowed.
 vol_constraint
 Defines whether the volume is allowed to fluctuate.
 baro_thermo
 NHCThermostat instance, coupled directly to the barostat
 vel_press0
 The initial barostat velocity tensor
 restart
 If true, the cell is not symmetrized initially

name
= 'PR'¶

kind
= 'deterministic'¶

method
= 'barostat'¶

init
(iterative)¶

pre
(iterative, chainvel0=None)¶

post
(iterative, chainvel0=None)¶

baro
(iterative, chainvel0)¶

add_press_cont
()¶

class
yaff.sampling.npt.
TadmorBarostat
(ff, temp, press, start=0, step=1, timecon=41341.37333664682, anisotropic=True, vol_constraint=False, baro_thermo=None, vel_press0=None, restart=False)¶ Bases:
yaff.sampling.verlet.VerletHook
This hook implements the TadmorMiller barostat for finite strains. The equations are derived in:
Tadmor, E. B.; Miller, R. E. Modeling Materials: Continuum, Atomistic and Multiscale Techniques 2011, 520527.Arguments:
 ff
 A ForceField instance.
 temp
 The temperature of thermostat.
 press
 The applied second PiolaKirchhoff tensor for the barostat.
Optional arguments:
 start
 The step at which the barostat becomes active.
 timecon
 The time constant of the TadmorMiller barostat.
 anisotropic
 Defines whether anisotropic cell fluctuations are allowed.
 vol_constraint
 Defines whether the volume is allowed to fluctuate.
 baro_thermo
 NHCThermostat instance, coupled directly to the barostat
 vel_press0
 The initial barostat velocity tensor
 restart
 If true, the cell is not symmetrized initially

name
= 'Tadmor'¶

kind
= 'deterministic'¶

method
= 'barostat'¶

init
(iterative)¶

pre
(iterative, chainvel0=None)¶

post
(iterative, chainvel0=None)¶

baro
(iterative, chainvel0)¶

add_press_cont
()¶

class
yaff.sampling.npt.
MTKAttributeStateItem
(attr)¶ Bases:
yaff.sampling.iterative.StateItem

get_value
(iterative)¶

copy
()¶

3.6. yaff.sampling.nvt
– Thermostats¶
Thermostats

class
yaff.sampling.nvt.
AndersenThermostat
(temp, start=0, step=1, select=None, annealing=1.0)¶ Bases:
yaff.sampling.verlet.VerletHook
This is an implementation of the Andersen thermostat. The method is described in:
Andersen, H. C. J. Chem. Phys. 1980, 72, 23842393.Arguments:
 temp
 The average temperature of the NVT ensemble
Optional arguments:
 start
 The first iteration at which this hook is called
 step
 The number of iterations between two subsequent calls to this hook.
 select
 An array of atom indexes to indicate which atoms controlled by the thermostat.
 annealing
 After every call to this hook, the temperature is multiplied with this annealing factor. This effectively cools down the system.

name
= 'Andersen'¶

kind
= 'stochastic'¶

method
= 'thermostat'¶

init
(iterative)¶

pre
(iterative, G1_add=None)¶

post
(iterative, G1_add=None)¶

class
yaff.sampling.nvt.
BerendsenThermostat
(temp, start=0, timecon=4134.137333664683, restart=False)¶ Bases:
yaff.sampling.verlet.VerletHook
This is an implementation of the Berendsen thermostat. The algorithm is described in:
Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; Dinola, A.; Haak, J. R. J. Chem. Phys. 1984, 81, 36843690Arguments:
 temp
 The temperature of thermostat.
Optional arguments:
 start
 The step at which the thermostat becomes active.
 timecon
 The time constant of the Berendsen thermostat.
 restart
 Indicates whether the initalisation should be carried out.

name
= 'Berendsen'¶

kind
= 'deterministic'¶

method
= 'thermostat'¶

init
(iterative)¶

pre
(iterative, G1_add=None)¶

post
(iterative, G1_add=None)¶

class
yaff.sampling.nvt.
LangevinThermostat
(temp, start=0, timecon=4134.137333664683)¶ Bases:
yaff.sampling.verlet.VerletHook
This is an implementation of the Langevin thermostat. The algorithm is described in:
Bussi, G.; Parrinello, M. Phys. Rev. E 2007, 75, 056707Arguments:
 temp
 The temperature of thermostat.
Optional arguments:
 start
 The step at which the thermostat becomes active.
 timecon
 The time constant of the Langevin thermostat.

name
= 'Langevin'¶

kind
= 'stochastic'¶

method
= 'thermostat'¶

init
(iterative)¶

pre
(iterative, G1_add=None)¶

post
(iterative, G1_add=None)¶

thermo
(iterative)¶

class
yaff.sampling.nvt.
CSVRThermostat
(temp, start=0, timecon=4134.137333664683)¶ Bases:
yaff.sampling.verlet.VerletHook
This is an implementation of the CSVR thermostat. The equations are derived in:
Bussi, G.; Donadio, D.; Parrinello, M. J. Chem. Phys. 2007, 126, 014101The implementation (used here) is derived in
Bussi, G.; Parrinello, M. Comput. Phys. Commun. 2008, 179, 2629Arguments:
 temp
 The temperature of thermostat.
Optional arguments:
 start
 The step at which the thermostat becomes active.
 timecon
 The time constant of the CSVR thermostat.

name
= 'CSVR'¶

kind
= 'stochastic'¶

method
= 'thermostat'¶

init
(iterative)¶

pre
(iterative, G1_add=None)¶

post
(iterative, G1_add=None)¶

class
yaff.sampling.nvt.
GLEThermostat
(temp, a_p, c_p=None, start=0)¶ Bases:
yaff.sampling.verlet.VerletHook
This hook implements the coloured noise thermostat. The equations are derived in:
Ceriotti, M.; Bussi, G.; Parrinello, M J. Chem. Theory Comput. 2010, 6, 11701180.Arguments:
 temp
 The temperature of thermostat.
 a_p
 Square drift matrix, with elements fitted to the specific problem.
Optional arguments:
 c_p
 Square static covariance matrix. In equilibrium, its elements are fixed. For nonequilibrium dynamics, its elements should be fitted.
 start
 The step at which the thermostat becomes active.

name
= 'GLE'¶

kind
= 'stochastic'¶

method
= 'thermostat'¶

init
(iterative)¶

pre
(iterative, G1_add=None)¶

post
(iterative, G1_add=None)¶

thermo
(iterative)¶

class
yaff.sampling.nvt.
NHCThermostat
(temp, start=0, timecon=4134.137333664683, chainlength=3, chain_pos0=None, chain_vel0=None, restart=False)¶ Bases:
yaff.sampling.verlet.VerletHook
This hook implements the NoseHoover chain thermostat. The equations are derived in:
Martyna, G. J.; Klein, M. L.; Tuckerman, M. J. Chem. Phys. 1992, 97, 26352643.The implementation (used here) of a symplectic integrator of the NoseHoover chain thermostat is discussed in:
Martyna, G. J.; Tuckerman, M. E.; Tobias, D. J.; Klein, M. L. Mol. Phys. 1996, 87, 11171157.Arguments:
 temp
 The temperature of thermostat.
Optional arguments:
 start
 The step at which the thermostat becomes active.
 timecon
 The time constant of the NoseHoover thermostat.
 chainlength
 The number of beads in the NoseHoover chain.
 chain_pos0
 The initial thermostat chain positions
 chain_vel0
 The initial thermostat chain velocities
 restart
 Indicates whether the initalisation should be carried out

name
= 'NHC'¶

kind
= 'deterministic'¶

method
= 'thermostat'¶

init
(iterative)¶

pre
(iterative, G1_add=None)¶

post
(iterative, G1_add=None)¶

class
yaff.sampling.nvt.
NHCAttributeStateItem
(attr)¶ Bases:
yaff.sampling.iterative.StateItem

get_value
(iterative)¶

copy
()¶

3.7. yaff.sampling.opt
– Geometry/Cell optimization¶
Geometry/Cell optimization

class
yaff.sampling.opt.
OptScreenLog
(start=0, step=1)¶ Bases:
yaff.sampling.iterative.Hook

__call__
(iterative)¶


class
yaff.sampling.opt.
BaseOptimizer
(dof, state=None, hooks=None, counter0=0)¶ Bases:
yaff.sampling.iterative.Iterative
Arguments:
 dof
 A specification of the degrees of freedom. The convergence criteria are also part of this argument. This must be a DOF instance.
Optional arguments:
 state
 A list with state items. State items are simple objects that take or derive a property from the current state of the iterative algorithm.
 hooks
 A function (or a list of functions) that is called after every iterative.
 counter0
 The counter value associated with the initial state.

default_state
= [<yaff.sampling.iterative.AttributeStateItem object>, <yaff.sampling.iterative.AttributeStateItem object>, <yaff.sampling.iterative.PosStateItem object>, <yaff.sampling.iterative.DipoleStateItem object>, <yaff.sampling.iterative.VolumeStateItem object>, <yaff.sampling.iterative.CellStateItem object>, <yaff.sampling.iterative.EPotContribStateItem object>]¶

log_name
= 'XXOPT'¶

fun
(x, do_gradient=False)¶

initialize
()¶

propagate
()¶

finalize
()¶

class
yaff.sampling.opt.
CGOptimizer
(dof, state=None, hooks=None, counter0=0)¶ Bases:
yaff.sampling.opt.BaseOptimizer

log_name
= 'CGOPT'¶

initialize
()¶

propagate
()¶


class
yaff.sampling.opt.
BFGSHessianModel
(ndof, hessian0=None)¶ Bases:
yaff.sampling.opt.HessianModel
Arguments:
 ndof
 The number of degrees of freedom
Optional arguments:
 hessian0
 An initial guess for the hessian

update
(dx, dg)¶

class
yaff.sampling.opt.
SR1HessianModel
(ndof, hessian0=None)¶ Bases:
yaff.sampling.opt.HessianModel
Arguments:
 ndof
 The number of degrees of freedom
Optional arguments:
 hessian0
 An initial guess for the hessian

update
(dx, dg)¶

class
yaff.sampling.opt.
QNOptimizer
(dof, state=None, hooks=None, counter0=0, trust_radius=1.0, small_radius=1e05, too_small_radius=1e10, hessian0=None)¶ Bases:
yaff.sampling.opt.BaseOptimizer
A Quasinewton optimizer
This is just a basic implementation of the algorithm, but it has the potential to become more advanced and efficient. The following improvements will be made when time permits:
 Support for nonlinear constraints. This should be relatively easy. We need a routine that can bring the unknowns back to the constraints, and a routine to solve a constrained second order problem with linear equality/inequality constraints. These should be methods of an object that is an attribute of the dof object, which is need to give the constraint code access to the Cartesian coordinates. In the code below, some comments are added to mark where the constraint methods should be called.
 The Hessian updates and the diagonalization are currently very slow for big systems. This can be fixed with a rank1 update algorithm for the spectral decomposition.
 The optimizer would become much more efficient if redundant coordinates were used. This can be implemented efficiently by using the same machinery as the constraint code, but using the dlist and iclist concepts for the sake of efficiency.
 It is in practice not needed to keep track of the full Hessian. The LBFGS algorithm is a nice method to obtain a linear memory usage and computational cost. However, LBFGS is not compatible with the trust radius used in this class, while we want to keep the trust radius for the sake of efficiency, robustness and support for constraints. Using the rank1 updates mentioned above, it should be relatively easy to keep track of the decomposition of a subspace of the Hessian. This subspace can be defined as the basis of the last N rank1 updates. Simple assumptions about the remainder of the spectrum should be sufficient to keep the algorithm efficient.
Arguments:
 dof
 A specification of the degrees of freedom. The convergence criteria are also part of this argument. This must be a DOF instance.
Optional arguments:
 state
 A list with state items. State items are simple objects that take or derive a property from the current state of the iterative algorithm.
 hooks
 A function (or a list of functions) that is called after every iterative.
 counter0
 The counter value associated with the initial state.
 trust_radius
 The initial value for the trust radius. It is adapted by the algorithm after every step. The adapted trust radius is never allowed to increase above this initial value.
 small_radius
 If the trust radius goes below this limit, the decrease in energy is no longer essential. Instead a decrease in the norm of the gradient is used to accept/reject a step.
 too_small_radius
 If the trust radius becomes smaller than this parameter, the optimizer gives up. Insanely small trust radii are typical for potential energy surfaces that are not entirely smooth.
 hessian0
 An initial guess for the Hessian

log_name
= 'QNOPT'¶

initialize
()¶

propagate
()¶

make_step
()¶

yaff.sampling.opt.
solve_trust_radius
(grad, evals, radius, threshold=1e05)¶ Find a step in eigen space with the given radius
3.8. yaff.sampling.trajectory
– Computations on a reference trajectory¶
Computations on a reference trajectory

class
yaff.sampling.trajectory.
TrajScreenLog
(start=0, step=1)¶ Bases:
yaff.sampling.iterative.Hook

__call__
(iterative)¶


class
yaff.sampling.trajectory.
RefTrajectory
(ff, fn_traj, state=None, hooks=None, counter0=0)¶ Bases:
yaff.sampling.iterative.Iterative
Arguments:
 ff
 A ForceField instance
fn_traj
A hdf5 file name containing the trajectoryOptional arguments:
 state
 A list with state items. State items are simple objects that take or derive a property from the current state of the iterative algorithm.
 hooks
 A function (or a list of functions) that is called after every iterative.
 counter0
 The counter value associated with the initial state.

default_state
= [<yaff.sampling.iterative.AttributeStateItem object>, <yaff.sampling.iterative.AttributeStateItem object>, <yaff.sampling.iterative.PosStateItem object>, <yaff.sampling.iterative.DipoleStateItem object>, <yaff.sampling.iterative.VolumeStateItem object>, <yaff.sampling.iterative.CellStateItem object>, <yaff.sampling.iterative.EPotContribStateItem object>]¶

log_name
= 'TRAJEC'¶

initialize
()¶

propagate
()¶

finalize
()¶
3.9. yaff.sampling.utils
– Auxiliary routines for initial velocities¶
Auxiliary routines for initial velocities

yaff.sampling.utils.
get_random_vel
(temp0, scalevel0, masses, select=None)¶ Generate random atomic velocities using a MaxwellBoltzmann distribution
Arguments:
 temp0
 The temperature for the MaxwellBoltzmann distribution.
 scalevel0
 When set to True, the velocities are rescaled such that the instantaneous temperature coincides with temp0.
 masses
 An (N,) array with atomic masses.
Optional arguments:
 select
 When given, this must be an array of integer indexes. Only for these atoms (masses) initial velocities will be generated.
Returns: An (N, 3) array with random velocities. When the select option is used, the shape of the results is (M, 3), where M is the length of the select array.

yaff.sampling.utils.
remove_com_moment
(vel, masses)¶ Zero the linear centerofmass momentum.
Arguments:
 vel
 An (N, 3) array with atomic velocities. This array is modified inplace.
 masses
 An (N,) array with atomic masses
The zero linear COM momentum is achieved by subtracting translational rigid body motion from the atomic velocities.

yaff.sampling.utils.
remove_angular_moment
(pos, vel, masses)¶ Zero the global angular momentum.
Arguments:
 pos
 An (N, 3) array with atomic positions. This array is not modified.
 vel
 An (N, 3) array with atomic velocities. This array is modified inplace.
 masses
 An (N,) array with atomic masses
The zero angular momentum is achieved by subtracting angular rigid body motion from the atomic velocities. (The angular momentum is measured with respect to the center of mass to avoid that this routine reintroduces a linear COM velocity. This is also beneficial for the numerical stability.)

yaff.sampling.utils.
clean_momenta
(pos, vel, masses, cell)¶ Remove any relevant external momenta
Arguments:
 pos
 An (N, 3) array with atomic positions. This array is not modified.
 vel
 An (N, 3) array with atomic velocities. This array is modified inplace.
 masses
 An (N,) array with atomic masses
 cell
 A Cell instance describing the periodic boundary conditions.

yaff.sampling.utils.
angular_moment
(pos, vel, masses)¶ Compute the angular moment of a set of point particles
Arguments:
 pos
 An (N, 3) array with atomic positions.
 vel
 An (N, 3) array with atomic velocities.
 masses
 An (N,) array with atomic masses.
Returns: a (3,) array with the angular momentum vector.

yaff.sampling.utils.
get_ndof_internal_md
(natom, nper)¶ Return the effective number of internal degrees of freedom for MD simulations
Arguments:
 natom
 The number of atoms
 nper
 The number of periodic boundary conditions (0 for isolated systems)

yaff.sampling.utils.
cell_symmetrize
(ff, vector_list=None, tensor_list=None)¶ Symmetrizes the unit cell tensor, and updates the position vectors
Arguments:
 ff
 A ForceField instance
Optional arguments:
 vector_list
 A list of numpy vectors which should be transformed under the symmetrization. Note that the positions are already transformed automatically
 tensor_list
 A list of numpy tensors of rank 2 which should be transformed under the symmetrization.

yaff.sampling.utils.
get_random_vel_press
(mass, temp)¶ Generates symmetric tensor of barostat velocities
Arguments:
 mass
 The Barostat mass
 temp
 The temperature at which the velocities are selected

yaff.sampling.utils.
get_ndof_baro
(dim, anisotropic, vol_constraint)¶ Calculates the number of degrees of freedom associated with the cell fluctuation
Arguments:
 dim
 The dimension of the system
 anisotropic
 Boolean value determining whether anisotropic cell fluctuations are allowed
 vol_constraint
 Boolean value determining whether the cell volume can change

yaff.sampling.utils.
stabilized_cholesky_decomp
(mat)¶ Do LDL^T and transform to MM^T with negative diagonal entries of D put equal to zero Assume mat is square and symmetric (but not necessarily positive definite).
3.10. yaff.sampling.verlet
– Generic Verlet integrator¶
Generic Verlet integrator

class
yaff.sampling.verlet.
TemperatureStateItem
¶ Bases:
yaff.sampling.iterative.StateItem

get_value
(iterative)¶

iter_attrs
(iterative)¶


class
yaff.sampling.verlet.
VerletIntegrator
(ff, timestep=None, state=None, hooks=None, vel0=None, temp0=300, scalevel0=True, time0=None, ndof=None, counter0=None, restart_h5=None)¶ Bases:
yaff.sampling.iterative.Iterative
Arguments:
 ff
 A ForceField instance
Optional arguments:
 timestep
 The integration time step (in atomic units)
 state
 A list with state items. State items are simple objects that take or derive a property from the current state of the iterative algorithm.
 hooks
 A function (or a list of functions) that is called after every iterative.
 vel0
 An array with initial velocities. If not given, random velocities are sampled from the MaxwellBoltzmann distribution corresponding to the optional arguments temp0 and scalevel0
 temp0
 The (initial) temperature for the random initial velocities
 scalevel0
 If True (the default), the random velocities are rescaled such that the instantaneous temperature coincides with temp0.
 time0
 The time associated with the initial state.
 ndof
 When given, this option overrides the number of degrees of freedom determined from internal heuristics. When ndof is not given, its default value depends on the thermostat used. In most cases it is 3*natom, except for the NHC thermostat where the number if internal degrees of freedom is counted. The ndof attribute is used to derive the temperature from the kinetic energy.
 counter0
 The counter value associated with the initial state.
 restart_h5
 HDF5 object containing the restart information

default_state
= [<yaff.sampling.iterative.AttributeStateItem object>, <yaff.sampling.iterative.AttributeStateItem object>, <yaff.sampling.iterative.AttributeStateItem object>, <yaff.sampling.iterative.PosStateItem object>, <yaff.sampling.iterative.AttributeStateItem object>, <yaff.sampling.iterative.AttributeStateItem object>, <yaff.sampling.iterative.AttributeStateItem object>, <yaff.sampling.iterative.AttributeStateItem object>, <yaff.sampling.verlet.TemperatureStateItem object>, <yaff.sampling.iterative.AttributeStateItem object>, <yaff.sampling.iterative.AttributeStateItem object>, <yaff.sampling.iterative.AttributeStateItem object>, <yaff.sampling.iterative.AttributeStateItem object>, <yaff.sampling.iterative.AttributeStateItem object>, <yaff.sampling.iterative.AttributeStateItem object>, <yaff.sampling.iterative.DipoleStateItem object>, <yaff.sampling.iterative.DipoleVelStateItem object>, <yaff.sampling.iterative.VolumeStateItem object>, <yaff.sampling.iterative.CellStateItem object>, <yaff.sampling.iterative.EPotContribStateItem object>]¶

log_name
= 'VERLET'¶

initialize
()¶

propagate
()¶

compute_properties
(restart_h5=None)¶

finalize
()¶

call_verlet_hooks
(kind)¶

class
yaff.sampling.verlet.
VerletHook
(start=0, step=1)¶ Bases:
yaff.sampling.iterative.Hook
Specialized Verlet hook.
This is mainly used for the implementation of thermostats and barostats.
Optional arguments:
 start
 The first iteration at which this hook should be called.
 step
 The hook will be called every step iterations.

__call__
(iterative)¶

init
(iterative)¶

pre
(iterative)¶

post
(iterative)¶

class
yaff.sampling.verlet.
VerletScreenLog
(start=0, step=1)¶ Bases:
yaff.sampling.iterative.Hook
A screen logger for the Verlet algorithm

__call__
(iterative)¶


class
yaff.sampling.verlet.
ConsErrTracker
(restart_h5=None)¶ Bases:
object
A class that tracks the errors on the conserved quantity. Given its superior numerical accuracy, the algorithm below is used to calculate the running average. Its properties are discussed in Donald Knuth’s Art of Computer Programming, vol. 2, p. 232, 3rd edition.

update
(ekin, econs)¶

get
()¶


class
yaff.sampling.verlet.
KineticAnnealing
(annealing=0.99999, select=None, start=0, step=1)¶ Bases:
yaff.sampling.verlet.VerletHook
This annealing hook is designed to be used with a plain Verlet integrator. At every call, the velocities are rescaled with the annealing parameter.
Arguments:
 annealing
 After every call to this hook, the temperature is multiplied with this annealing factor. This effectively cools down the system.
 select
 An array mask or a list of indexes to indicate which atomic velocities should be annealed.
 start
 The first iteration at which this hook is called
 step
 The number of iterations between two subsequent calls to this hook.

init
(iterative)¶

pre
(iterative)¶

post
(iterative)¶