-
This section covers the parameters that are set in the towhee_input file Version 8.1.x. Each parameter is listed along with its type
(logical, character, integer, or double precision). towhee_input is the main input file for Towhee and is generally the only file
that needs to be edited on a regular basis. It has a regimented style to the input. The parameters are described here in the order
they appear in this file. Please look at one of the example files (available with the code package) for the precise file format.
- 8.1.2: Added wmultibody.F to the svn database as it was mistakenly ommitted previously. Many minor bug fixes in routines related to unneeded commas in output commands. Fixed a minor bug in the Sum et al. 2003 forcefield file that resulted in a nonsense line as an additional possible parameter for one of the bond length terms.
- 8.1.1: Fixed a bug that was resulting in compilation warnings on some platforms related to a subroutine call without enough variables in the rigid three-bond cyclic section of atd_cbmc_setup.F.
- 8.1.0: The additional thermodynamic pressure calculation associated with the pressure_thermo_freq parameter is now functional. The default setting for pressure_thermo_style is now 'range'. Fixed a problem with the dihedral parameters in the Sum et al. 2003 version of the NERD force field related to a mistake in the implementation of the Sum2003 Cosine series (see the towhee_ff manual torsion section). Rebuilt all of the NERD family of forcefields (NERD v1, NERD v2, NERD v3, and Sum et al. 2003) so they can be mixed together, although the user should be aware that the order they are listed in towhee_input will affect some of the intramolecular parameters as Towhee prioritizes the terms that appear first when searching with the molecule assembler.
- Added input parameters for the 2 box center-of-mass switch move and related parameters (pm2boxcomswitch, pm2comswboxpair, pm2comswmolpair). However, this move is not currently functional and will result in an error message if used.
Convenient links to towhee_input information:
- Nonbonded potential types are described in the classical_potential section.
- Information used to describe molecule-specific forcefield information is in the input_style section.
- Various constants declared at compile-time found in preproc.h.
# character are ignored.
Many parameters are not included in a typical towhee_input file.
Optional input file values are marked as optional parameter and if that parameter is not explicitly included in towhee_input it is instead
set specified default value. The default values are either noted as (default: x),
where x is the default value, or with (default) indicated
for parameters that have a specified list of options. Parameters that are only required based upon the values of other input parameters are either
formatted so they are indented compared to their parent parameter, or are followed by the notation only if 'another parameter' is 'some value'.
- 'Towhee' : (default) read the input parameters following the format for Towhee. This format is described in this file.
- 'LAMMPS' : read the input parameters from the lammps_input and lammps_data files. Outputs files suitable for use with Towhee.
- 'DX-1597-2-7' : (default) use the DX-1597-2-7 pseudorandom number generator described in Deng 2005 with a period of approximately 1014903.
- 'KISS99' : use a version of the Keep it Simple Stupid pseudorandom number generator originally posted to an internet newsgroup by G. Marsaglia. This has a period of approximately 109.
- 'MRG32k3a' : use the MRG32k3a pseudorandom number generator described in l'Ecuyer 1999 with a period of approximately 1057.
- 'RANLUX' : use the RANLUX psuedorandom number generator as described by
Lüscher 1994 and implemented by
James 1994 with a period of approximately 10171.
- When using this setting the following optional parameter is available
- random_luxlevel (integer) optional parameter
-
An integer value that defines the luxury level used by the pseudorandom number generator.
See James 1994 for more information about the luxury level.
- 0: Does not discard any pseudorandom numbers from the sequence.
- 1: Periodically discard 24 pseudorandom numbers from the sequence.
- 2: Periodically discard 73 pseudorandom numbers from the sequence.
- 3: (default) Periodically discard 199 pseudorandom numbers from the sequence.
- 4: Periodically discard 365 pseudorandom numbers from the sequence.
- When using this setting the following optional parameter is available
-
(default: 1302002) The 32 bit integer seed used to initialize the psuedorandom number generator
(unless performing a restart when random_allow_restart is
.TRUE.). Must be positive.
.TRUE.(default) restart the pseudorandom number generator from the sequence of integers in towhee_initial (if possible). This is done when the simulation is restarting from a towhee_initial file that was generated using Version 5.0.0 or later. If this is not possible then Towhee will initialize the pseudorandom number generator using the random_seed..FALSE.initialize the pseudorandom number generator using the random_seed.
- 'npt' : Isobaric-Isothermal Ensemble. The volume moves for each simulation box are performed as an exchange with an external pressure bath (set with pressure). The total number of molecules and the temperature are conserved. Combining this ensemble with any of the two-box molecule transfer moves turns it into the NPT-Gibbs ensemble.
- 'nvt' : Canonical Ensemble. The total volume, the total number of each type of molecule, and the temperature are conserved for entire system. Combining this ensemble with any of the two-box molecule transfer moves turns it into the NVT-Gibbs ensemble.
- 'uvt' : Grand Canonical Ensemble. The total volume and temperature are conserved. The total number of molecules in the system equilibrates with an external bath set at a specified chemical potential (set with chempot).
-
The temperature in Kelvin.
-
The external pressure in kPa.
-
The total number of molecule types in the simulation. Must be less then or equal to NTMAX.
-
The total chemical potential (includes intramolecular energy and phase space contributions) for molecules of each type (listed sequentially on a single line).
The units are in Kelvin (Free energy / Boltzmann's constant).
See the chemical potential algorithm page for more information about this complicated property.
-
The number of simulation boxes in the system. This value must be less than or equal to MAXBOX
Note that many other parameters depend upon numboxes and some Monte Carlo moves are only valid for multiple box ensembles.
-
Note that starting with version 8.0.0 this option merely changes the default values of certain towhee_input parameters, but they may still be
reassigned manually.
- 'equilibration' : designed for use when you are first setting up a simulation and still need to rapidly equilibrate all of the maximum update
values.
This option changes the default value of the following parameters based upon the value of nstep.
- printfreq = nstep/10
- blocksize = nstep/5
- moviefreq = nstep/5
- backupfreq = nstep/5
- pdb_output_freq = nstep/5
- pressure_virial_freq = 10
- trmaxdispfreq = nstep/100
- volmaxdispfreq = nstep/100
- 'production' : designed for use during production runs.
This option changes the default value of the following parameters based upon the value of nstep.
- printfreq = nstep/10
- blocksize = nstep/5
- moviefreq = nstep/20
- backupfreq = nstep/5
- pdb_output_freq = nstep/5
- pressure_virial_freq = 10
- trmaxdispfreq = nstep/10
- volmaxdispfreq = nstep/10
- 'manual' : (default) the default values are not changed.
-
(default: 0, also see controlstyle)
The step frequency for outputting information about the system to standard output.
The information is the number of Monte Carlo steps performed thus far during the run, the total energy in each box, the x-box length of each box,
the pressure of each box, and the number of molecules of each type in each box. If printfreq = 0 this information is not output.
-
(default: 0, also see controlstyle)
The size of the blocks for computing block averages.
It is best if blocksize divides cleanly into nstep.
The quantities that are averaged (in each simulation box) are the specific density, the pressure, all of the energy terms, the chemical
potential of each molecule type, number density of each molecule type, and the mole fractions. If blocksize is 0 then no block output takes place.
-
(default: 0, also see controlstyle)
The step frequency for outputting the system information to the towhee_movie file.
This file can be analyzed after the run using the analyse_movie utility routine to compute a variety
of distribution functions. This file gets quite large if you output frequently so be careful if you have a limited amount of hard disk space available.
If moviefreq is 0 this file is not written.
-
(default: 0, also see controlstyle)
The step frequency for writing a file named towhee_backup that is suitable for use as a restart file.
It overwrites the previous version of towhee_backup each time so it does not take up much space.
For more information about restart files see the Towhee restart manual.
If backupfreq is 0 this file is not written.
-
(default: 0)
The step frequency for writing a file named towhee_restart_NNN that is suitable for use as a restart file,
where NNN is the step number when this file is written. This file is identical to towhee_backup (see backupfreq),
except that it is not overwritten by subsequent restart files, and therefore allows for restarts from multiple points along a run.
For more information about restart files see the Towhee restart manual.
If restartfreq is 0 this file is not written.
- 'full' : (default) output information about the individual blocks of the block averages and information about the maximum displacement updates.
- 'blocks' : output information about the individual blocks of the block averages but not about the maximum displacement updates.
- 'updates' : does not output information about the individual blocks of the block averages but does output information about the maximum displacement updates.
- 'none' : does not output information about the individual blocks of the block averages or information about the maximum displacement updates.
-
(default: 0) The step frequency for outputting a snapshot of the simulation to a pdb file named box_xx_step_yyyyyyyyyyyyyy.pdb where xx is the box number
converted into a 2 character string and yyyyyyyyyyyyyy is the step number converted into a 14 character string.
No pdb files are output if pdb_output_freq is 0.
.FALSE.(default) do not output dft files..TRUE.output files for use with the Tramonto classical density functional theory code. This outputs dft_surfaces.dat and dft_decode.dat. See the Tramonto software link for information about these files.
.FALSE.(default) do not output LAMMPS files..TRUE.output files for use with the LAMMPS massively parallel molecular dynamics code. This outputs lammps_input and lammps_data# where the number is each of the simulation box numbers. See the LAMMPS software link for more information on how to read in these files.
.FALSE.(default) do not output DL_POLY files..TRUE.output files for use with the DL_POLY molecular dynamics code. This outputs CONFIG# and FIELD# where the number is each of the simulation box numbers. See the DL_POLY software link for more information about these files.
.FALSE.(default) do not perform Transition-Matrix Monte Carlo (TMMC)..TRUE.perform Transition-Matrix Monte Carlo (TMMC) on a single component system in the grand canonical ensemble. This requires other parameters to be specified below.- n_tmmc_min (integer) optional parameter,
only if tmmc_flag is
.TRUE.-
(default: 0) Minimum number of molecules allowed in the simulation box. Specifying a non-zero value is useful when sampling
specific molecule number ranges. Note that the corresponding n_tmmc_max is automatically set to the total number of molecules in the simulation.
- weight_freq (integer) optional parameter,
only if tmmc_flag is
.TRUE.-
(default: 0) The step frequency for updating the biasing function for Transition-Matrix Monte Carlo.
Updates are not performed if weight_freq is 0.
- c_matrix_freq (integer) optional parameter,
only if tmmc_flag is
.TRUE.-
(default: 0) The step frequency for writing to file the collection matrix for TMMC. A new output file containing the step number is created.
These files are not output if c_matrix_freq is 0.
- run_name (character) optional parameter,
only if tmmc_flag is
.TRUE.-
(default: 'tmmc_runfile')
File prefix for TMMC output.
The current estimate of the natural logarithm of the particle number probability distribution is written to a
file "run_name.tmmc_weights.dat" every time the biasing function is updated. Collection matrix info is written to
a file "run_name.c.stepnumber.dat". Accumulated semigrand potential energies are written to a file "run_name.vsg.stepnumber".
This information is useful if one is interested in the average potential energy as a function density.
- in_c_flag (logical) optional parameter,
only if tmmc_flag is
.TRUE..FALSE.: (default) start with the default collection matrix..TRUE.: read an initial collection matrix from a file and use that to initiate the TMMC simulation. This option is used to continue a previous TMMC simulation. Future versions could also use this to read in either a collection matrix or biasing function.- in_cfile (character) only if in_c_flag is
.TRUE.-
Text prepended to the collection matrix input file.
- in_cfile (character) only if in_c_flag is
- n_tmmc_min (integer) optional parameter,
only if tmmc_flag is
.FALSE.(default) histogram reweighting file are not output..TRUE.output files used for histogram reweighting. This value requires the following additional parameters- hist_label (integer) only if louthist is
.TRUE.-
An integer that is turned into a character string that creates the X portion of the named towhee_hisXY.dat file that is used to
output information for histogram reweighting.
- hist_suffix (character*1) only if louthist is
.TRUE.-
A single character that creates the Y portion of the named towhee_hisXY.dat file that is used to output information for histogram
reweighting.
- hist_nequil (integer) only if louthist is
.TRUE.-
The number of steps (cycles or moves) that are disregarded for the purposes of outputting histogram information. Set to 0 if you wish to
use all of the steps for computing the histogram, and set to a positive number in order to discard those initial steps from the histogram.
- histcalcfreq (integer) only if louthist is
.TRUE.-
The step frequency for computing the information needed for histogram reweighting.
- histdumpfreq (integer) only if louthist is
.TRUE.-
The step frequency for outputting the information needed for histogram reweighting to the various towhee_histogram files. The ratio of
histdumpfreq/histcalcfreq must be less than NDUMPHIST.
- hist_label (integer) only if louthist is
-
(default: 0, also see controlstyle)
The step frequency for computing the pressure using a deterministic algorithm that is based upon the interaction of pairs of molecules in
each simulation box. This is performed using the molecular version of the pressure virial (for continuous potentials) or the radial pressure
(for discontinuous potentials).
See the Pressure Manual for more information about the pressure calculation methods
implemented into Towhee.
Computing the pressure using these algorithms is a relatively expensive task (compared to the cost of a Monte Carlo move), especially for large systems.
This calculation is not performed if pressure_virial_freq is 0.
-
(default: 0)
The step frequency for computing the pressure by sampling the average energy change/volume change ratio to compute the "thermodynamic" pressure
each simulation box. This calculation is not yet functional (as of version 8.0.0).
This calculation is not performed if pressure_thermo_freq is 0.
-
The style of attempted volume change to use when computing the thermodynamic pressure via additional ghost volume changes.
- 'move data': attempts volume changes in the same manner as performing a single box isobaric-isothermal volume move.
- 'range': (default) attempts a ghost volume change of size Volume * Uniform[pressure_thermo_range(0), pressure_thermo_range(1)].
-
(default: -0.01d0 0.01d0)
The range for the attempted ghost volume changes used when computing the thermodynamic pressure.
-
(default: 0, also see controlstyle)
The step frequency for updating the maximum translational (atom and center-of-mass) and rotational displacements. They are periodically adjusted to
achieve the target acceptance rates (see tatraa, tatrac,
and tarot). When initially equilibrating a simulation it is a good idea to do this frequently
(every step or every 10 steps) in order to quickly converge to useful values for the maximum displacements. Once the acceptance rates are
near the target values then setting trmaxdispfreq to perform about 10 updates per simulation is probably sufficient.
Updates to these quantities are not performed if trmaxdispfreq is 0.
-
(default: 0, also see controlstyle)
The step frequency for updating the maximum volume displacements. They are periodically adjusted to achieve the target acceptance rate
(see tavol).
When initially equilibrating a simulation it is a good idea to do this frequently (every few steps) in order quickly converge to useful values
for the maximum displacements. Once the acceptance rates are near their desired values it is suggested to set volmaxdisp to perform
around 10 updates during the simulation.
Updates to the maximum volume displacements are not performed if volmaxdispfreq is set to 0.
-
(default: 0)
The number of additional trial insertions to perform in each box for at the end of every Monte Carlo step (listed sequentially for each molecule
type on a single line). This allows the measurement of chemical potential in ensembles that do not have an insertion and deletion move (such as
canonical and isobaric-isothermal).
-
.TRUE.start a new simulation and generate the positions of all of the atoms, assign initial box dimensions and maximum displacements. -
.FALSE.continue the simulation by reading in box lengths, maximum displacements, and coordinates from towhee_initial.
- 'dimensions' : the dimensions of the initial boxes are entered in order to construct the initial boxes.
- 'number density' : the total number density of molecules in the initial boxes are entered in order to compute the initial sizes of the boxes. This option generates cubic boxes.
- 'unit cell' : generates an initial structure by duplicating a unit cell. Reads information from the towhee_cell file and uses that to create an initial structure.
-
One line for each simulation box (see numboxes).
Each line contains an entry for each molecule type (see nmolty).
- 'full cbmc' : a template for this molecule type is created using configurational-bias and then replicated throughout the simulation box to generate an initial configuration. This is the easiest method for the user, but it can struggle with especially large or complicated molecules.
- 'coords' : coordinates for each atom are read from towhee_coords. This is useful if you are starting from a different file format (such as pdb), or have another code for building the initial configuration.
- 'helix cbmc' : a molecule template is generated by placing some of the backbone atoms onto a helix and then growing the rest of the atoms using CBMC.
This option requires the following parameters.
- helix_moltyp (integer)
-
An integer corresponding to the molecule type that had an initstyle of 'helix cbmc' in any simulation box.
Listed in consecutive order by molecule type.
- helix_radius (double precision)
-
The radius of the helix (units of Angstroms).
- helix_angle (double precision)
-
The pitch angle the helix makes with respect to the z-axis (units of degrees).
- helix_keytype (character*10)
- 'element' compares the helix_keyname with the character*2 parameter element that contains the 2 letter elemental code for each atom.
- 'nbname' compares the helix_keyname with the character*10 nonbond Atom Names for each atom type. This is the same parameter that is used when inputting the atom names with input_style options 'basic connectivity map' or 'advanced connectivity map'.
- 'pdbname' compares the helix_keyname with the parameter pdbname that contains the 4 character code used in the pdb format output. This is intended to work with input_style options 'polypeptide builder' or 'nucleic acid builder'.
- helix_keyname (character*10)
-
The key for finding matches of the atom with the data structures for the molecule that is being grown as a helix. You need to choose an atom name
that only appears in the backbone (e.g. 'P' for Charmm27 nucleic acids when using the element keytype, or ' CA ' for the C alpha backbone carbon
of a polypeptide when using the pdbname keytype).
- helix_conlen (double precision)
-
The distance between consecutive helix_element atoms (units of Angstroms).
- helix_phase (double precision)
-
The initial angle of the helical chain (units of degrees). Normally, this has little effect as it is just a rotation about the z-axis, but
if you are trying to set up two complementary nucleic acid chains to form a double helix then you would want their phase angles to differ
by 180 degrees.
- helix_moltyp (integer)
- 'nanotube' : The coordinates for each atom are read from towhee_nanotube. This file is generated by the Towhee code if you use the input_style for carbon nanotubes. This template is then replicated throughout the simulation box to generate an initial conformation.
- 'template' : A template for this molecule type is read from towhee_template. This template is then replicated throughout the simulation box to generate an initial configuration.
- 'partial cbmc' : The molecule is constructed in two steps. First, a partial list of atom positions, amino acid 3-letter codes, and 4-letter atom identifiers (following the pdb standard) is read from towhee_partial and matched up against the expected amino acid and atom codes from the molecule template listed in towhee_input. Then configurational-bias is used to grow all of the missing atoms to create the initial structure.
-
One line for each simulation box in the system. Each line contains a value for each molecule type.
- 'center' : put the center of mass of the molecule in the very center of the box. This is a very poor idea if you have multiples of the same molecule type in a box as they will all be right on top of each other. However, it is useful if you are trying to get multiple different molecule types in an initial configuration where one molecule is inside of another (for example: multi-walled nanotubes).
- 'none' : place the molecules exactly where their template indicates. This option makes the most sense for initstyle options like 'coords' where positions are read from a file.
- 'simple cubic' : places the first atom of the molecule onto a simple cubic lattice and the rest of the atoms in the molecule are placed relative to that atom. If you are not using an initstyle of 'coords' then this is usually the best option to generate an initial structure.
-
The initial number of each type of molecule in each box (one line per box).
-
The initial number of molecules (for initboxtype settings of 'dimensions' or 'number density') or of duplicated unit cells
(for initboxtype setting of 'unit cell') in each direction in each box.
The product of inix*iniy*iniz must be greater than or equal to the initial number of molecules in that box
(for initboxtype settings of 'dimensions' or 'number density'). While these are labeled x, y, and z
they correspond to the three coordinate vectors a, b, and c (equivalent to x, y, and z for a rectangular box).
-
The initial box dimensions (Angstroms) for the three vectors that describe the simulation box. There are nine entries (3 for each of the
3 vectors) in total for each simulation box. These are listed one vector at a time, with the three numbers which make up each vector listed on the
same line. Note that the coordinate system you choose does not have to be orthogonal, but it must follow the right hand rule. The three vectors
must also all be at least 45 degrees apart. Note that if you wish to use a rectangular box then only the diagonal elements of hmatrix will be
non-zero, and these will be equal to the box lengths in the x, y, and z dimensions.
-
The initial total number density of molecules in each simulation box. Listed as a single value per line with one line for each box as specified
by numboxes. Units are molecules per nm3.
Monte Carlo Move Section: Notes on the pm* parameters
-
The pm* parameters are used to determine which move type to perform every time we want to do a Monte Carlo move. A move is selected
by choosing a pseudorandom number between 0.0 and 1.0 and then going down the list of pm* until one is found that has a value higher than the
pseudorandom number. At least one of the parameters must be set to 1.0.
A similar procedure is performed when we want to determine which boxes or molecule types to perform the selected move upon.
These are done using the pm**pr and pm**mt arrays.
Starting with Towhee release 4.15.2, all of the leading pm* cumulative probabilities are optional, but if a cumulative probability is not listed,
then the rest of the parameters in that section are also not listed, and the cumulative probability of that move is zero.
For example, if the parameter
pmback is not included in towhee_input, then
pmbkmt must not be in the input file, and the Configurational-Bias Protein Backbone Regrowth move is never attempted.
-
(default: 0.0)
Cumulative probability of performing an isotropic volume move.
If this move is performed when ensemble is 'npt' then a single box is selected and it exchanges volume with an external
pressure bath (see pressure).
If this move is performed when ensemble is 'nvt' and numboxes > 1 then
a pair of boxes are selected and volume is exchanged between them.
- pmvlpr (double precision array)
-
Cumulative probability of performing a volume move on a particular box, or box pair. All of these parameters are listed on a single line.
If (ensemble = 'npt') then a value of pmvlpr is listed for each box.
If (ensemble = 'nvt') then a value is listed for each pair of simulation boxes where the pairs are ordered
(1, 2), (1, 3), ... (1, numboxes),
(2, 3), ... (numboxes-1, numboxes).
- rmvol (double precision)
-
The initial volume maximum displacement for all boxes or box pairs (a single initial value that is used for all boxes).
In the isobaric-isothermal ensemble this is the initial maximum volume displacement (cubic Angstroms) in each box.
In the canonical Gibbs ensemble this is the maximum displacement (logarithmic space) for each pair of boxes.
As the simulation progresses, these values are updated independently for each box, or each pair of boxes
(see volmaxdispfreq).
- tavol (double precision)
-
The target acceptance rate for the volume move. Must be a value between 0.0 and 1.0. The volume displacement array is periodically
adjusted (see volmaxdispfreq) to yield this acceptance rate.
The most popular value is 0.5, though some researchers advocate smaller values.
-
(default: 0.0)
Cumulative probability of performing an Anisotropic unit cell adjustment volume move.
If ensemble is 'npt' then a single box is selected and a single hmatrix element
is changed resulting in a volume exchange with an external pressure bath (see pressure).
If ensemble is 'nvt' and numboxes > 1: then a pair of boxes is selected and
one of those boxes is selected according to the pmcellpt
parameter and a single hmatrix element is changed in that box.
This results in a change of volume for the selected box that is countered by isotropically changing the volume in the other box.
- pmcellpr (double precision array)
-
Cumulative probability of performing a unit cell adjustment move on a particular box, or box pair. All of these parameters are listed on a single line
If (ensemble = 'npt') then a value of pmcellpr is listed for each box.
If (ensemble = 'nvt') then a value is listed for each pair of
simulation boxes where the pairs are ordered (1, 2), (1, 3), ... (1, numboxes), (2, 3),
... (numboxes-1, numboxes).
- pmcellpt (double precision array)
-
Cumulative probability of selecting the first box of the pair as the box to perform the non-isotropic volume move upon, while its partner undergoes an
isotropic volume move. This parameter is only used when (ensemble = 'nvt').
Note that you can choose to perform the non-isotropic portion of this move always on the same box and this is useful when performing a
solid-vapor equilibria calculation.
- rmcell (double precision)
-
The initial unit cell adjustment maximum displacement. In all cases, this is the maximum amount (in Angstroms) that a single element of the
hmatrix can change in a single unit cell move.
Note, the in the canonical Gibbs ensemble case it is possible for the isotropic box to undergo an hmatrix
change that is larger than this value as that box is countering the volume change caused by the non-isotropic adjustment in the
first box.
These values are independently updated for each box/pair according to volmaxdispfreq.
- tacell (double precision)
-
The target acceptance rate for the unit cell adjustment move. Must be a value between 0.0 and 1.0. The unit cell displacement (rmcell) is
periodically adjusted (see volmaxdispfreq) to yield this acceptance rate. I typically use a value of 0.5.
-
(default: 0.0)
Cumulative probability of performing a rotational-bias interbox molecule transfer move. This move takes a molecule out of one box and tries to place it
in another box. The molecule is grown using nch_nb_one attempted different orientations and positions of the center-of-mass.
- pm2rbswmt (double precision array)
-
Cumulative probability of performing a rotational-bias interbox molecule transfer move on each type of molecule in the system.
- pm2rbswpr (double precision array)
-
(default: 0.0)
Cumulative probability of performing a configurational-bias interbox molecule transfer move. This move takes a molecule out of one box and tries
to place it in another box. The molecule is grown using coupled-decoupled configurational-bias Monte Carlo.
- pm2cbswmt (double precision array)
-
Cumulative probability of performing a configurational-bias interbox molecule transfer move on each type of molecule in the system.
- pm2cbswpr (double precision array)
-
(default: 0.0)
Cumulative probability of performing a 2 box center-of-mass switch move. This move selects two boxes (a and b) and two molecule types (x and y). It then
randomly selects whether box a or b loses a molecule of type x and gains a molecule of type y, while the other box loses a molecule of type y and
gains a molecule of type x. Molecules are inserted at the same center-of-mass position where the molecule of the other type was deleted.
The coordinates relative to the center-of-mass of each inserted molecule are equal to the relative coordinates of the molecule of the same type
that was deleted.
- pm2comswboxpair (double precision array)
- pm2comswmolpair (double precision array)
-
(default: 0.0)
Cumulative probability of performing a grand-canonical configurational-bias insertion or deletion move.
- pmuvtcbmt (double precision array)
-
Cumulative probability of performing a grand-canonical configurational-bias insertion or deletion move on each type of molecule in the system.
-
Cumulative probability of performing an intrabox configurational-bias reinsertion move. This move takes a molecule out of one box and tries to place
it back into the same box. The molecule is grown using configurational-bias Monte Carlo.
- pm1cbswmt (double precision array)
-
Cumulative probability of performing an intrabox configurational-bias molecule transfer move on each type of molecule in the system.
-
(default: 0.0)
Cumulative probability of performing an intrabox exchange of the center of mass of two molecules of different types.
- pm1comswbox (double precision array)
-
Cumulative probability of performing a 1 box center-of-mass switch move in each simulation box. List one value for each simulation box.
- pm1comswpair (double precision array)
-
(default: 0.0)
Cumulative probability of performing an aggregation volume bias move of type 1, as described in
Chen and Siepmann 2000. This is useful for forming and destroying clusters in simulations with
molecules that tend to aggregate together.
- pmavb1in (double precision)
-
Probability of trying to move a molecule into an inner region for aggregation volume bias move of type 1. Must be strictly between 0.0 and 1.0.
- pmavb1mt (double precision array)
-
Cumulative probability of performing an aggregation volume bias move of type 1 where a molecule of a certain type is moved.
This is an array with one element for each molecule type.
- pmavb1ct (double precision array)
-
Cumulative probability of performing an aggregation volume bias move of type 1 where the molecule target is of a certain type. The molecule that is
moved is chosen according to pmavb1mt and then the type of molecule that is used as a reference for determining the
inner and outer regions is found using this parameter. This is a two dimensional array and uses one line of text for each type of molecule in the system.
- avb1rad (double precision)
-
The radius used to define the inner and outer volumes in the aggregation volume bias move of type 1. The distance is specified in Angstroms
and must be greater than zero, but less than or equal to rcut.
-
(default: 0.0)
Cumulative probability of performing an aggregation volume bias move of type 2, as described in
Chen and Siepmann 2001.
This is useful for forming and destroying clusters in simulations with molecules that tend to aggregate together.
- pmavb2in (double precision)
-
Probability of trying to move a molecule into an inner region for aggregation volume bias move of type 2. Must be strictly between 0.0 and 1.0.
- pmavb2mt (double precision array)
-
Cumulative probability of performing an aggregation volume bias move of type 2 where a molecule of a certain type is moved. This is an array with one
element for each molecule type.
- pmavb2ct (double precision array)
-
Cumulative probability of performing an aggregation volume bias move of type 2 where the molecule target is of a certain type. The molecule that is moved
is chosen according to pmavb2mt and then the type of molecule that is used as a reference for determining the inner and
outer regions is found using this parameter. This is a two dimensional array and uses one line of text for each type of molecule in the system.
- avb2rad (double precision)
-
The radius used to define the inner and outer volumes in the aggregation volume bias move of type 2. The distance is specified in Angstroms
and must be greater than zero, but less than or equal to rcut.
-
(default: 0.0)
Cumulative probability of performing an aggregation volume bias move of type 3, as described in
Chen and Siepmann 2001. This is useful for transferring molecules between clusters.
- pmavb3mt (double precision array)
-
Cumulative probability of performing an aggregation volume bias move of type 3 where a molecule of a certain type is moved.
This is an array with one element for each molecule type.
- pmavb3ct (double precision array)
-
Cumulative probability of performing an aggregation volume bias move of type 3 where the molecule target is of a certain type. The molecule that is moved
is chosen according to pmavb3mt and then the type of molecule that is used as a reference for determining the inner and
outer regions is found using this parameter. This is a two dimensional array and uses one line of text for each type of molecule in the system.
- avb3rad (double precision)
-
The radius used to define the inner and outer volumes in the aggregation volume bias move of type 3. The distance is specified in Angstroms
and must be greater than zero, but less than or equal to rcut.
-
(default: 0.0)
Cumulative probability of performing a molecule regrowth move on a molecule without regard to which box the molecule is currently located in. This move
chooses a molecule of the appropriate type at random, selects an atom of the molecule at random, and then regrows the molecule either entirely
(if a random number < pmall) or in all directions except for one.
The molecule is regrown using configurational-bias.
- pmcbmt (double precision array)
-
Cumulative probability of performing a molecule regrowth on each type of molecule in the system.
- pmall (double precision array)
-
Probability that a molecule regrowth move will regrow the entire molecule. This is listed for each molecule type in the simulation.
-
(default: 0.0)
Cumulative probability of performing configurational-bias fixed-endpoint regrowth of a portion of the protein backbone. This selects an atom along the
peptide backbone, chooses another backbone atom that is connected by three bonds to the first atom, and then regrows all of the atoms in between
these two atoms.
- pmbkmt (double precision array)
-
Cumulative probability of performing a backbone regrowth move on each type of molecule in the system.
-
(default: 0.0)
Cumulative probability of performing a configurational-bias regrowth move on a single side-chain of peptide. This works by compiling a list of
atoms that have a pdbname of ' CA ' and selecting one of those at the origin of a side-chain regrowth. This is functional for the normal amino
acids, plus proline and disulfide bonded cysteines.
- pmcbsidemt (double precision array)
-
Cumulative probability of performing a configurational-bias side-chain regrowth on each type of molecule in the system.
-
(default: 0.0)
Cumulative probability of performing a torsional pivot move about a random bond in the molecule. This move chooses a bond that is not entirely contained in a
single ring structure, and has at least one bond emanating from each end, and then rotates one side of the molecule about that bond.
- pmpivmt (double precision array)
-
Cumulative probability of performing a pivot move on each type of molecule in the system.
-
(default: 0.0)
Cumulative probability of performing a concerted rotation move for a sequence of 9 atoms in a molecule.
- pmcrmt (double precision array)
-
Cumulative probability of performing a concerted rotation move move on each type of molecule in the system.
-
(default: 0.0)
Cumulative probability of performing a concerted rotation move on a sequence of three peptides in a polypeptide.
This move only works for polypeptides.
- pmcrbmt (double precision array)
-
Cumulative probability of performing a protein backbone concerted rotation move on each type of molecule in the system.
-
(default: 0.0)
Cumulative probability of performing a plane shift move. This move displaces all of the molecules whose center of mass lies in a plane of width
planewidth. A new trial position for the center of the plane of atoms is generated uniformly across the available plane.
- pmplanebox (double precision array)
-
Cumulative probability of performing a plane shift in each of the simulation boxes. List one value for each simulation box.
- planewidth (double precision)
-
The width of the plane for the plane shift move. Any molecule whose center of mass is within a plane of this thickness (whose position is
chosen uniformly along one axis) will move during the plane shift move. The value of planewidth must be greater than 0.0d0 and less than the
shortest box length.
-
(default: 0.0)
Cumulative probability of performing a row shift move. This move displaces all of the molecules whose center of mass lies in a row of diameter
rowwidth. A new trial position for the center of the row of atoms is generated uniformly across the available row.
- pmrowbox (double precision array)
-
Cumulative probability of performing a row shift in each of the simulation boxes. List one value for each simulation box.
- rowwidth (double precision)
-
The width of the plan for the row shift move. Any molecule whose center of mass is within a row of this thickness (whose position is chosen
uniformly along one axis) will move during the row shift move. The value of rowwidth must be greater than 0.0d0 and less than the shortest
box length.
-
(default: 0.0)
Cumulative probability of performing a single-atom translation move on a molecule without regard to which box the molecule is currently located in. This
move chooses a molecule of the appropriate type at random, selects an atom of the molecule at random, selects a vector on a unit sphere at random,
and then attempts to displace the atom a random distance between -rmtraa and +rmtraa in that direction.
- pmtamt (double precision array)
-
Cumulativev probability of performing a single-atom translation move on each type of molecule in the system.
- rmtraa (double precision)
-
The initial Atom-translation maximum displacement (Angstroms) for all molecules types in all boxes. As the simulation progresses, these values
are updated to yield the desired acceptance rate for each molecule type in each box (see trmaxdispfreq).
- tatraa (double precision)
-
The target acceptance rate for the atom translation move. Must be a value between 0.0 and 1.0. The maximum atom translational displacement
is periodically adjusted (see trmaxdispfreq) to yield this acceptance rate.
The most common choice is a value of 0.5, though some researchers prefer smaller values.
-
(default: 0.0)
Cumulative probability of performing a composite move that consists of a random center-of-mass translation and random rotation. This move is essentially
a concatenation of the Center-of-Mass Molecule Translation Move and the Rotation about the Center-of-Mass Move.
- pmcomt (double precision array)
-
Cumulative probability of performing a composite move on each type of molecule in the system.
- rmcomtra (double precision)
-
The molecular translation displacement (angstroms) for all molecule types in all boxes. This parameter is very similar
to the rmtrac parameter for translations of molecular COM.
- rmcomrot (double precision)
-
The molecular rotation maximum displacement (radians) for all molecule types in all boxes. This parameter is very similar to the
rmrot parameter for rotations about COM.
-
(default: 0.0)
Cumulative probability of performing a center-of-mass translation move on a molecule without regard to which box the molecule is currently located in.
This move chooses a molecule of the appropriate type at random, chooses a vector on a unit sphere at random, and then attempts to displace the
entire molecule a random distance between -rmtrac and +rmtrac in that direction.
- pmtcmt (double precision array)
-
Cumulative probability of performing a center-of-mass translation move on each type of molecule in the system.
- rmtrac (double precision)
-
The initial Center-of-mass translation maximum displacement (Angstroms) for all molecule types in all boxes. As the simulation progresses,
these values are updated to yield the desired acceptance rate for each molecule type in each box
(see trmaxdispfreq).
- tatrac (double precision)
-
The target acceptance rate for the center-of-mass translation move. Must be a value between 0.0 and 1.0. The maximum center-of-mass
translational displacement is periodically adjusted (see trmaxdispfreq) to yield this acceptance rate.
A value of 0.5 is common, though some researchers prefer smaller values.
-
(default: 0.0)
Cumulative probability of performing a rotation about the center-of-mass move for a molecule without regard to the box the molecule is currently
located in. This move chooses a molecule of the appropriate type at random and then attempts to rotate the entire molecule about the x, y, and z
axes that run through the center-of-mass a random number of radians between -rmrot and
+rmrot around each of the three axes.
- pmromt (double precision array)
-
Cumulative probability of performing a rotation move on each type of molecule in the system.
- rmrot (double precision)
-
The initial molecular rotation maximum displacement (radians) for all molecule types in all boxes. As the simulation progresses, these values
are updated to yield the desired acceptance rate for each molecule type in each box (see trmaxdispfreq).
- tarot (double precision)
-
The target acceptance rate for the rotation move. Must be a value between 0.0 and 1.0. The rotation displacement is periodically
adjusted (see trmaxdispfreq) to yield this acceptance rate.
A value of 0.5 is common, though some researchers prefer smaller values.
- 'normal' : (default) A relatively small amount of output is output to towhee_cbmc.
- 'full' : A substantial amount of analysis about the predicted efficiency of configurational-bias algorithms is output to towhee_cbmc.
- 'Martin and Frischknecht 2006' : (default) Monte Carlo moves that involve the configurational-bias algorithm use the Coupled to pre-nonbond formulation of the coupled-decoupled algorithm that is described in Martin and Frischknecht 2006. This algorithm has a bond length selection, two bending angle selections, and a dihedral selection that are all decoupled from each other; while in turn being coupled to a pre-nonbond loop that is coupled to the nonbonded selection.
- 'Martin and Siepmann 1999 + Martin and Thompson 2004' : Monte Carlo moves that involve the configurational-bias algorithm use the coupled-decoupled formulation presented in the appendix of Martin and Siepmann 1999 with the addition of a decoupled bond selection first described in Martin and Thompson 2004. This algorithm has decoupled bond length selection, two decoupled angle selections, and a dihedral selection that is coupled to the nonbonded selection.
-
A convenient method for setting many parameters related to configuration-bias Monte Carlo with a single parameter.
Starting with version 8.0.0 this sets a suite of parameters while allowing individual parameters to subsequently be reset.
- 'Martin and Frischknecht' : (default) sets up the configurational-bias values according to the strategies
described in
Martin and Frischknecht 2006. This option sets the following parameters for all molecule
types in the system. For more information about these parameters see the explicit setting for this parameter.
- cbmc_fit_strategy = 'equilibrium'
- cbmc_nb_one_generation = 'uniform'
- nch_nb_one = 10
- nch_nb = 10
- nch_pre_nb =
- 1: if cbmc_formulation is 'Martin and Siepmann 1999 + Martin and Thompson 2004'
- 100: if cbmc_formulation is 'Martin and Frischknecht 2006'
- cbmc_dihedral_generation = 'ideal + autofit gaussian'
- dihedral_peak_weight_style = 'uniform'
- dihedral_sdev_multiplier = 1.0d0
- dihedral_ideal_fraction = 0.01d0
- nch_tor =
- 100: if cbmc_formulation is 'Martin and Siepmann 1999 + Martin and Thompson 2004'
- 1: if cbmc_formulation is 'Martin and Frischknecht 2006'
- nch_tor_connect =
- 100: if cbmc_formulation is 'Martin and Siepmann 1999 + Martin and Thompson 2004'
- 1: if cbmc_formulation is 'Martin and Frischknecht 2006'
- cbmc_bend_generation = 'ideal + autofit gaussian'
- bend_a_sdev_multiplier = 1.0d0
- bend_b_sdev_multiplier = 1.0d0
- bend_a_ideal_fraction = 0.01d0
- bend_b_ideal_fraction = 0.01d0
- nch_bend_a = 1
- nch_bend_b =
- 100: if cbmc_formulation is 'Martin and Siepmann 1999 + Martin and Thompson 2004'
- 1: if cbmc_formulation is 'Martin and Frischknecht 2006'
- cbmc_bond_generation = 'autofit gaussian'
- bond_sdev_multiplier = 1.0d0
- nch_vib = 1
- two_bond_fixed_endpoint_bias_style = 'analytic Boltzmann using angles'
- three_bond_fixed_endpoint_bias_style = 'analytic using max and min 2-4 distance'
- 'default autofit gaussian' : sets up all of the configurational-bias parameters to autofit gaussian settings. This option sets the following
parameters for all molecule types in the system. For more information about these parameters see the explicit setting for this parameter.
- cbmc_fit_strategy = 'equilibrium'
- cbmc_nb_one_generation = 'uniform'
- nch_nb_one = 10
- nch_nb = 10
- nch_pre_nb =
- 1: if cbmc_formulation is 'Martin and Siepmann 1999 + Martin and Thompson 2004'
- 10: if cbmc_formulation is 'Martin and Frischknecht 2006'
- cbmc_dihedral_generation = 'autofit gaussian'
- dihedral_peak_weight_style = 'uniform'
- dihedral_sdev_multiplier = 1.0d0
- nch_tor =
- 10: if cbmc_formulation is 'Martin and Siepmann 1999 + Martin and Thompson 2004'
- 1: if cbmc_formulation is 'Martin and Frischknecht 2006'
- nch_tor_connect =
- 10: if cbmc_formulation is 'Martin and Siepmann 1999 + Martin and Thompson 2004'
- 1: if cbmc_formulation is 'Martin and Frischknecht 2006'
- cbmc_bend_generation = 'autofit gaussian'
- bend_a_sdev_multiplier = 1.0d0
- bend_b_sdev_multiplier = 1.0d0
- nch_bend_a = 1
- nch_bend_b = 1
- cbmc_bond_generation = 'autofit gaussian'
- bond_sdev_multiplier = 1.0d0
- nch_vib = 1
- two_bond_fixed_endpoint_bias_style = 'autofit gaussian'
- two_bond_bias_sdev_multiplier = 1.0d0
- two_bond_bias_vibrange = 0.5d0 1.5d0
- three_bond_fixed_endpoint_bias_style = 'autofit gaussian using max and min 2-4 distance'
- three_bond_bias_sdev_multipler = 1.0d0
- 'default ideal' : sets up all of the configurational-bias parameters according to a general set of default values using the 'ideal' generation
styles. This option sets the following parameters for all molecule types in the system. For more information about these parameters see the explicit
setting for this parameter.
- cbmc_fit_strategy = 'equilibrium'
- cbmc_nb_one_generation = 'uniform'
- nch_nb_one = 10
- nch_nb = 10
- nch_pre_nb = 1
- cbmc_dihedral_generation = 'ideal'
- nch_tor = 360
- nch_tor_connect = 360
- cbmc_bend_generation = 'ideal'
- nch_bend_a = 1000
- nch_bend_b = 1000
- cbmc_bond_generation = 'r^2 with bounds'
- vibrang = 0.85 1.15
- nch_vib = 1000
- two_bond_fixed_endpoint_bias_style = 'analytic Boltzmann using angles and dihedrals'
- three_bond_fixed_endpoint_bias_style = 'analytic using max and min 2-4 distance'
- 'explicit' : no special settings for the default configurational-bias parameters.
- 'Martin and Thompson FPE 2004' : sets up all of the configurational-bias parameters to the values used in
Martin and Thompson 2004. This option sets the following parameters for all molecule types
in the system. For more information about these parameters see the explicit setting for this parameter.
- cbmc_fit_strategy = 'equilibrium'
- cbmc_nb_one_generation = 'uniform'
- nch_nb_one = 10
- nch_nb = 10
- nch_pre_nb = 1
- cbmc_dihedral_generation = 'ideal'
- nch_tor = 360
- nch_tor_connect = 360
- cbmc_bend_generation = 'ideal'
- nch_bend_a = 100
- nch_bend_b = 100
- cbmc_bond_generation = 'r^2 with bounds'
- vibrang = 0.85 1.15
- nch_vib = 1000
- two_bond_fixed_endpoint_bias_style = 'analytic Boltzmann using angles and dihedrals'
- three_bond_fixed_endpoint_bias_style = 'analytic using max and min 2-4 distance'
- 'widom ideal' : sets up all of the configurational-bias parameters according to a general set of default values using the 'ideal' generation styles
a single growth trial for all of the nonbonded selections following the (now discredited, see Martin 2013)
belief that this was the only
way to insure the chemical potential was computed correctly using the normal Widom insertion method. This option sets the following parameters for
all molecule types in the system. For more information about these parameters see the explicit setting for this parameter.
- cbmc_fit_strategy = 'equilibrium'
- cbmc_nb_one_generation = 'uniform'
- nch_nb_one = 1
- nch_nb = 1
- nch_pre_nb = 1
- cbmc_dihedral_generation = 'ideal'
- nch_tor = 360
- nch_tor_connect = 360
- cbmc_bend_generation = 'ideal'
- nch_bend_a = 1000
- nch_bend_b = 1000
- cbmc_bond_generation = 'r^2 with bounds'
- vibrang = 0.85 1.15
- nch_vib = 1000
- two_bond_fixed_endpoint_bias_style = 'none'
- three_bond_fixed_endpoint_bias_style = 'none'
- 'equilibrium' : (default) The equilibrium values are used when fitting the distributions of higher order terms for configurational-bias. For example, equilibrium bond lengths are used when fitting the bending angle distributions.
-
(default: 'uniform' for all molecule types and boxes)
One line for each simulation box and each line contains one value for each molecule type.
- 'energy bias' : the trial coordinates of the first atom inserted into this simulation box are generated by selecting a subvolume according to
the Boltzmann weighted energy of a test atom in the center of each subvolume, averaged over all possible insertion atoms in the mapped molecule.
Once a subvolume is selected then a position in that subvolume is generated uniformly. This biasing scheme was first used by
Snurr et al. 1993. Note that the energy map used in practice is a linear combination
of the Boltzmann weighted map described above (0.99 of the fraction) with a random selection of any subcube in the box (0.01 of the fraction) to
make sure that the generation probability is strictly positive in all subcubes. The insertion energies are computed using the initial positions
of all molecules in simulation box 1. If you are using this to help with insertions in a porous solid then you may want to only compute this
for the empty porous solid and use the restart of that map in future simulations.
When using this option you need the following parameters.
- mapmolty (integer)
-
The molecule type to use as the insertion probe for determining the energy biasing weights.
- lcreatemap (logical)
.TRUE.create a new towhee_map file. This file contains information about the energy profile of the porous molecule. This setting requires the following additional parameters.- cubex,cubey,cubez (integer)
-
The number of cubelets you want to have in each direction of the box. Each dimension of the box is divided by the corresponding cube*
value and for each cube an energy value is computed. This value is then used for biasing insertion/deletion moves.
- cubex,cubey,cubez (integer)
.FALSE.use a previously generated map in the towhee_map file. That file must have been generated using version 5.2.4 or later as there was a different, problematic, format prior to that version.
- mapmolty (integer)
- 'uniform' : the trial coordinates of first atom inserted into this simulation box are generated uniformly across the entire volume.
-
(default: 10, for all molecule types)
The number of trial positions that are sampled for the first atom inserted during a
configurational-bias or rotational-bias molecule exchange move
(see pm2boxrbswap, pm2boxcbswap,
and pm1boxcbswap). Listed on a single line with one value for each molecule type.
-
(default: 10, for all molecule types)
The number of trial positions that are sampled for all nonbonded energy trials except for the first atom inserted during a
configurational-bias molecule exchange move (see pm2boxcbswap
and pm1boxcbswap). This value is used for all atoms grown during a pmcb regrowth move.
Listed on a single line with one value for each molecule type.
-
(default: 1, for all molecule types)
This parameter is only meaningful when the cbmc_formulation is 'Martin and Frischknecht 2006'. The number of trials for a selection
procedure that takes place after the dihedral selection, but before the nonbond selection.
Listed on a single line with one value for each molecule type.
- 'autofit gaussian' : dihedral trials are generated according to a series of gaussian distributions that are individually fit to the Boltzmann
factor as a function of dihedral angle for each individual dihedral in the system. This fit is performed automatically at the start of each
simulation using the equilibrium bond lengths and bending angles. This bias is then removed in the acceptance rule. When using this option
you must also include the following parameter.
- dihedral_peak_weight_style see dihedral_peak_weight_style for information about this shared parameter.
- dihedral_sdev_multiplier see dihedral_sdev_multiplier for information about this shared parameter.
- 'global gaussian' : dihedral trials are generated according a series of gaussian distributions that are a function of the bondpatt
on the central two atoms. This series is hard coded into Towhee and details can be found by looking in the getcbdihed.F subroutine. This
gaussian bias is then removed in the acceptance rule. When using this option you must also include the following parameter.
- sdevtor (double precision) optional parameter,
also see cbmc_setting_style
-
(default: 20.0)
The standard deviation (units of degrees) used for each of the gaussian distributions for the dihedral angles. For best results
this number should be set to the observed distribution computed in the analyse_movie utility.
This number must be positive.
- sdevtor (double precision) optional parameter,
also see cbmc_setting_style
- 'ideal' : (default) dihedral trials are generated according to the ideal distribution. For dihedrals that means trials are generated uniformly on (-π, π).
- 'ideal + autofit gaussian' : dihedral trials are generated according to a linear combination of the ideal distribution (uniform) and a series
of gaussian distributions that are individually fit to the Boltzmann factor as a function of dihedral angle for each individual dihedral in the
system. This fit is performed automatically at the start of each simulation using the equilibrium bond lengths and bending angles. This bias
is then removed in the acceptance rule. When using this option you must also include the following parameters.
- dihedral_peak_weight_style see dihedral_peak_weight_style for information about this shared parameter.
- dihedral_sdev_multiplier see dihedral_sdev_multiplier for information about this shared parameter.
- dihedral_ideal_fraction (double precision) optional parameter,
also see cbmc_setting_style
-
(default: 0.01)
The fraction of dihedral trials that are generated using the ideal distribution. The remainder of the trials are generated using the autofit
gaussians. This number must be in the range [0.0,1.0] inclusive.
Shared parameters for cbmc_dihedral_generation.
- dihedral_peak_weight_style (character*30) optional parameter,
also see cbmc_setting_style
- 'uniform' : (default) each peak is selected with the same (uniform) probability.
- 'isolated Boltzmann' : the peaks are selected with a probability proportional to the sum of the Boltzmann weight computed during the automatic fitting process.
- dihedral_sdev_multiplier (double precision) optional parameter,
also see cbmc_setting_style
-
(default: 1.0)
The factor that is multiplied by the observed gaussian standard deviation for each peak in the automatic dihedral fit in order to create the
standard deviations that are used to generate the dihedrals during the simulation. This number must be positive and a default
value of 1.0 is suggested.
-
(default: 360, for all molecule types)
The number of trial dihedral angles sampled during configurational-bias moves. The value must be positive.
-
(default: 360, for all molecule types)
The number of trial dihedral angles that are sampled during configurational-bias moves when we have
grown the molecule such that we need to connect back up with atoms that already exist. This is needed in order to regrow cyclic molecules,
and to regrow the interiors of large molecules.
The value must be positive.
- 'autofit gaussian' : bending trials are generated according to gaussian distributions. The bending A trial distribution is created by fitting
the mean and standard deviation of a gaussian to sin(θA)e-β Ubend(θA) on
(0, π) while treating the relevant bond lengths according to the cbmc_fit_strategy.
The standard devation used to generate trials during the simulation is equal to the fitted standard deviation times the
bend_a_sdev_multiplier.
The bending B trial distributions are generated by considering the psuedo-dihedral angle (φ) between two bending angles that share the first and
second atoms. φ is sampled on (-π, π) while treating the other bond angles and bond lengths according to the
cbmc_fit_strategy and each value of φ implies a θ between the two angles and that is used to compute
the distribution of e-β Ubend(θB). This distribution is split into sections depending upon the number of
minima and then a gaussian is fit to the distribution in each subdivision.
The standard deviation used to generate trials during the simulation is equal to the fitted standard deviation times
bend_b_sdev_multiplier.
Using this option enables the following parameters.
- bend_a_sdev_multiplier see bend_a_sdev_multiplier for information about this shared parameter.
- bend_b_sdev_multiplier see bend_b_sdev_multiplier for information about this shared parameter.
- 'global gaussian' : bending trials are generated according to gaussian distributions. Bending A trials are generated according to a single
gaussian with a mean set to the equilibrium bending angle and a standard deviation set to sdevbena. Bending B trials are generated
according to one or more gaussian distributions with means set based upon the bondpatt of the central atom and a standard deviation
set to setbenb. For more details on the bending B distribution see the
getcbangle.Fsubroutine. Using this option enables the following parameters.- sdevbena (double precision) optional parameter,
also see cbmc_setting_style
-
(default: 5.0)
standard deviation used to generate all of the part A bending trials (units of degrees). Must be positive, and the sampling (and acceptance rate)
is optimized when this value closely matches the observed standard distribution of the angles as measured by the
analyse_movie utility.
- sdevbenb (double precision) optional parameter,
also see cbmc_setting_style
-
(default: 5.0)
standard deviation used to generate the part B bending trials (units of degrees). Must be positive.
- sdevbena (double precision) optional parameter,
also see cbmc_setting_style
- 'ideal' : (default) bending trials generated according to the ideal distribution. This is a normalized Sine distribution on (0, π) for bending angle type A and uniform on (-π, π) for bending angle type B.
- 'ideal + autofit gaussian' : bending trials are generated according to a linear combination of the 'ideal' and the 'autofit gaussian' distributions.
See both of those sections for details.
This option generates bending A trials as bend_a_ideal_fraction 'ideal'
+ (1.0-bend_a_ideal_fraction) 'autofit gaussian'.
This option generates bending B trials as bend_b_ideal_fraction 'ideal'
+ (1.0-bend_b_ideal_fraction) 'autofit gaussian'.
- bend_a_sdev_multiplier see bend_a_sdev_multiplier for information about this shared parameter.
- bend_b_sdev_multiplier see bend_b_sdev_multiplier for information about this shared parameter.
- bend_a_ideal_fraction (double precision) optional parameter,
also see cbmc_setting_style
-
(default: 0.01)
The fraction of bending A trials that are generated using the 'ideal' distribution. The remainder of the trials are generated
using the 'autofit gaussian' distribution. This value must be in the range [0.0,1.0] inclusive.
- bend_b_ideal_fraction (double precision) optional parameter,
also see cbmc_setting_style
-
(default: 0.01)
The fraction of bending B trials that are generated using the 'ideal' distribution. The remainder of the trials are
generated using the 'autofit gaussian' distribution. This value must be in the range [0.0,1.0] inclusive.
Shared parameters for cbmc_bend_generation.
- bend_a_sdev_multiplier (double precision) optional parameter,
also see cbmc_setting_style
-
(default: 1.0)
This value is multiplied by the observed standard deviation from performing the fit in order to create the standard deviation that is used to
generate bending A trials during the simulation. This value must be positive.
- bend_b_sdev_multiplier (double precision) optional parameter,
also see cbmc_setting_style
-
(default: 1.0)
This value is multiplied by the observed standard deviation from performing the fit in order to create the standard deviation that is used to
generate bending B trials during the simulation. This value must be positive.
-
(default: 100, for each molecule type)
The number of trial angles that are sampled during configurational-bias moves when we are selecting the
iugrow-iufrom-iuprev angle. This value must be positive. Currently suggested values are in the range of 100 to 1000 when using
cbmc_bend_generation style 'ideal'
and in the range from 1 to 10 when using cbmc_bend_generation styles 'global gaussian'
or 'autofit gaussian'.
-
(default: 100, for each molecule type)
The number of trial angles that are sampled during configurational-bias moves when we are selecting the
rotation about a cone of one of the iugrow angles relative to the others. This value must be positive. Currently suggested values are in
the range of 100 to 1000 when using cbmc_bend_generation style 'ideal' and in the range from 1 to 10 when using
cbmc_bend_generation styles 'global gaussian' or 'autofit gaussian'.
-
(default: 3.0)
The maximum allowed value for any bond length in the simulation. Any bond lengths greater than this value are given an infinite energy
regardless of the actual form of the bond length potential. This value is also used in order to normalize the bond length distribution in a
configurational-bias move and therefore will affect the measured chemical potential.
- 'autofit gaussian' : Generate trial bond lengths according to a gaussian distribution with a mean and standard deviation fitted to
r2 e-β Ubond(r). This option enables the following parameter.
- bond_sdev_multiplier (double precision) optional parameter,
also see cbmc_setting_style
-
(default: 1.0)
This value is multiplied by the observed standard deviation of the r2 e-β Ubond(r) distribution in order
to determine the standard deviation used to generate bond trials. This value must be positive.
- bond_sdev_multiplier (double precision) optional parameter,
also see cbmc_setting_style
- 'global gaussian' : Generate trial bond lengths according to a gaussian distribution with a mean set to the equilibrium bond length and a
standard deviation specified as sdevvib. When using this option you must also include the following parameter.
- sdevvib (double precision) optional parameter,
also see cbmc_setting_style
-
(default: 0.1)
The standard deviation of a gaussian distribution that is used to sample bond lengths during a
configurational-bias regrowth for a cbmc_bond_generation style of 'global gaussian'.
Units are Angstroms. For best results perform a short simulation of single-atom translation moves, analyse that data using the
analyse_movie utility, and set this value to the observed standard deviations in the bond length
distribution.
- sdevvib (double precision) optional parameter,
also see cbmc_setting_style
- 'ideal' : (default) Generate trial bond lengths according to the ideal r2 probability distribution within the range (0, max_bond_length).
- 'r^2 with bounds' :
Generate trial bond lengths according to a bounded r2 probability distribution within the ranges set by
the vibrang parameter.
This distribution is proportional to the true distribution, but has a limited sampling range while the true
distribution is of infinite extent. When using this option you also need to include the following parameter.
- vibrang (double precision, double precision) optional parameter,
also see cbmc_setting_style
-
(default: 0.85 1.15)
The range of bond lengths to sample via configurational-bias Monte Carlo. The range is expressed as a
fraction of the equilibrium bond length for the lower bound and the upper bound.
- vibrang (double precision, double precision) optional parameter,
also see cbmc_setting_style
-
(default: 1000, for each molecule type)
The number of trial bond lengths sampled during a configurational-bias move. This value must be
positive. Currently suggested values are 1000 trials for cbmc_bond_generation style 'r^2 with bounds' and a value in the range
of 1 to 10 for cbmc_bond_generation styles 'global gaussian' and 'autofit gaussian'.
- 'none' : (default) No additional biasing utilized during a configurational-bias move step that involves the growth of a new atom that is separated from an already existing atom by two bonds where the atom between those two bonds does not currently exist.
- 'analytic Boltzmann using angles' : Additional biasing is utilized during a configurational-bias move step that involves the growth of a new
atom that is separated from an already existing atom by two bonds where the atom between those two bonds does not currently exist. The biasing
has the following form where i is the atom being grown, k is an atom that already exists, and j is an atom that is bonded to i and j, but
has not yet been grown.
- pbias(rik) = Min[minbias,paBua].
- paBua
-
= pangle for rik < rij0 + rjk0
-
= pbond for rij0 + rjk0 < rik
- pangle = e[- β Uangle( θijk ) ]
-
where θijk is computed from the trial rik distance and the equilibrium bond lengths of the two missing bonds
(rij0 and rjk0).
- pbond = e[- β ( Uangle( θijk ) + Ubond( rij ) + Ubond( rjk ) ) ]
minbias is a minimum value set in the code in order to avoid division by zero. This is set in the febias.F subroutine and currently has a value of 1.0d-40.-
where θijk = &pi, rij = rij0 * rik
/ ( rij0 + rjk0 ),
and
rjk = rjk0 * rik / ( rij0 + rjk0 ),
- 'analytic Boltzmann dihedral energy sum' : Additional biasing is utilized during a configurational-bias move step that involves the growth of
a new atom that is separated from an already existing atom by two bonds where the atom between those two bonds does not currently exist.
The biasing has the following form where i is the atom being grown, h is an atom that already exists and is the atom from which i is being grown
this step, k is an atom that already exists, and j is an atom that is bonded to i and j, but has not yet been grown.
- pbias(rij) = Min[minbias,paBdes].
- paBdes
-
= pdihedral for rik < rij0 + rjk0
-
= minbias for rij0 + rjk0 < rik
- pdihedral = e[- β (Udihedral(φhijk(1)) + Udihedral(φhijk(2)) )]
-
where φhijk(1) and φhijk(2) are the two possible solutions for that dihedral given the following constraints,
- rij = rij0
- rjk = rjk0
- θhij = θhij0
- 'autofit gaussian' : Additional biasing is utilized during a configurational-bias move step that involves the growth of a new atom that is
separated from an already existing atom by two bonds, where the atom between those two bonds does not currently exist. The biasing is a
gaussian distribution that is fit to the following distribution where i is the atom being grown, k is an atom that already exists, and j is an
atom that is bonded to i and j, but has not yet been grown.
- pag = rij2 rjk2 θijk e[-β ( Ubond(rij) + Ubond(rjk) + Uangle(θijk) )]
- where this distribution is computed at the beginning of the simulation and then the mean and standard deviations are stored for use during the configurational-bias growth procedure. When using this option you must also include the following additional parameters.
- two_bond_bias_sdev_multiplier see two_bond_bias_sdev_multiplier for information about this shared parameter.
- two_bond_bias_vibrange (double precision array) optional parameter,
also see cbmc_setting_style
-
(default: 0.85 1.15)
the lower and upper bounds of sampling for the rij and rjk bond lengths that are varied at the start of
the simulation in order to determine the gaussian fit parameters. They are in units relative to the equilibrium bond lengths.
- 'self adapting gaussian using 1-3 distance' : Additional biasing is utilized during a configurational-bias move step that involves the growth
of a new atom that is separated from an already existing atom by two bonds, where the atom between those two bonds does not currently exist.
The biasing is a gaussian distribution based upon the distance between the atom being grown and the target atom that is two bonds away and
already exists. This gaussian distribution is self adapted during the course of a simulation so that it represents the observed 1-3 distance
distribution. When using this option you must also include the following additional parameters.
- two_bond_bias_sdev_multiplier see two_bond_bias_sdev_multiplier for information about this shared parameter.
- two_bond_bias_initial_value (character*50) optional parameter,
also see cbmc_setting_style
- 'autofit gaussian' : (default) the initial distribution is generated by automatically fitting a gaussian to a sampling of the rik distance. This option is recommended when starting a new simulation where there is no appropriate 'towhee_safe_initial' file available.
- 'file' : the initial distribution is read from the 'towhee_safe_initial' file in the local directory. This file is generally copied from the 'towhee_safe_final' file that is produced at the end of a simulation that employed this two_bond_fixed_endpoint_bias_style.
- two_bond_bias_compute_frequency (integer) optional parameter,
also see cbmc_setting_style
-
(default: 0)
statistics on the two bond bias distributions are taken from the simulation with a step frequency equal to this value. These statistics
are then periodically used to update the distributions with frequency
two_bond_bias_update_frequency.
A value of 0 disables the periodic computation of these distributions.
- two_bond_bias_update_frequency (integer) optional parameter,
also see cbmc_setting_style
-
(default: 0)
frequency for updating the self-adapting two bond biasing functions. The updates combine the distributions
computed with a frequency controlled by the two_bond_bias_compute_frequency parameter
with the previous version of the two bond biasing potentials to create the new two bond biasing potential.
A value of 0 disables this update.
- two_bond_bias_old_fraction (double precision) optional parameter,
also see cbmc_setting_style
-
(default: 0.5)
This factor determines the linear combination of the old distribution, and the observed distribution, to determine the new two bond bias
distribution. This number must be in the range [0.0,1.0] inclusive. Setting this value to 0.0 would completely replace the old
distribution with the new distribution, while a setting of 0.5 would combine the old and observed distributions equally in order to compute
the new distribution.
Shared parameters for two_bond_fixed_endpoint_bias_style.
- two_bond_bias_sdev_multiplier (double precision) optional parameter,
also see cbmc_setting_style
-
(default: 1.0)
this value is multiplied by the standard deviation determined from the fit to get the standard deviation used during the simulation. Must be positive.
- 'none' : No additional biasing is utilized during a configurational-bias move step that involves the growth of a new atom that is separated from an already existing atom by three bonds where the atoms between those three bonds do not currently exist.
- 'analytic using max and min 2-4 distance' : Additional biasing is utilized during a configurational-bias move step that involves the growth
of a new atom that is separated from an already existing atom by three bonds where the atoms between those three bonds do not currently exist.
The biasing has the following form where i is the atom being grown, h is an atom that already exists and is the atom from which i is being grown
this step, l is an atom that already exists, and j and k are atoms that have not yet been regrown and bridge the gap between atoms i and l.
The biasing with this option depends upon the minimum (rjlmin) and maximum (rjlmax) projected
distances between the j and l atoms given the constraints that rij and θhij are set to their equilibrium
values (rij0 and θhij0). These minimum and maximum distances are then compared with the
equilibrium distance between the j and l atoms (rjleq). This equilibrium distance is computed by setting
rjk = rjk0, rkl = rkl0,
and θjkl = θjkl0). The biasing value is computed as follows.
- paumam24d
-
= pstretch for rjleq < rjlmin
-
= 1.0 for rjlmin < rjleq < rjlmax
-
= pcompress for rjlmax < rjleq
- pstretch
-
= e[- β Uangle(θjkl)] where θjkl is computed given the constraints of
rjk0, rkl0, and rjlmin for rjlmin
< rjk0 + rkl0.
-
= e[- β ( Uangle(θjkl) + Ubond(rjk)
+ Ubond(rkl) )] where θjkl = π,
rjk = rjk0 * rjlmin
/ ( rjk0 + rkl0 ),
and
rkl = rkl0 * rjlmin
/ ( rjk0 + rkl0 )
for rjk0 + rkl0 < rjlmin.
- pcompress = e[- β Uangle(θjkl) ]
-
where θjkl is computed using the constraints rjk0, rkl0,
and rjlmax.
- paumam24d
- 'autofit gaussian using max and min 2-4 distance' : Additional biasing is utilized during a configurational-bias move step that involves
the growth of a new atom that is separated from an already existing atom by three bonds where the atoms between those three bonds do not
currently exist. The biasing has the following form where i is the atom being grown, h is an atom that already exists and is the atom from
which i is being grown this step, l is an atom that already exists, and j and k are atoms that have not yet been regrown and bridge
the gap between atoms i and l. The biasing with this option depends upon the minimum (rjlmin) and
maximum (rjlmax) projected distances between the j and l atoms given the constraints that rij
and θhij are set to their equilibrium values (rij0 and θhij0).
These extrema are then combined with the 'autofit gaussian' biasing described for the
two_bond_fixed_endpoint_bias_style option. The bias probability is set to the
integrated probability of that autofit gaussian on the limits between rjlmin and rjlmax.
When using this option you also need to include the following parameter.
- three_bond_bias_sdev_multiplier see three_bond_bias_sdev_multiplier for information about this shared parameter.
- 'self adapting gaussian using 1-4 distance' : Additional biasing is utilized during a configurational-bias move step that involves the
growth of a new atom that is separated from an already existing atom by three bonds where the atoms between those three bonds do not
currently exist. The biasing has the following form where i is the atom being grown, h is an atom that already exists and is the atom from
which i is being grown this step, l is an atom that already exists, and j and k are atoms that have not yet been regrown and bridge
the gap between atoms i and l. The biasing is a gaussian that depends only upon the ril distance. This gaussian distribution
can be set to self-adapt during the course of the simulation so that the biasing probability reflects the observed probability for the
ril distances. When using this option you also need to include the following parameters.
- three_bond_bias_sdev_multiplier see three_bond_bias_sdev_multiplier for information about this shared parameter.
- three_bond_bias_initial_value (character*50) optional parameter,
also see cbmc_setting_style
- 'autofit gaussian' : (default) the initial distribution is automatically fit to the Boltzmann weight as a function of ril distances using the equilibrium rij, rjk, and rkl bond lengths and the equilibrium θijk and θjkl angles. This option is recommended for the initial simulation when no appropriate 'towhee_safe_initial' file exists.
- 'file' : the initial distribution is read from the 'towhee_safe_initial' file in the local directory. This file is copied from the 'towhee_safe_final' file that is generated at the end of a simulation run that employed the appropriate three_bond_fixed_endpoint_bias_style.
- three_bond_bias_compute_frequency (integer) optional parameter,
also see cbmc_setting_style
-
(default: 0)
Statistics on the three bond bias distributions are taken from the simulation with a step frequency equal to this value. These
statistics are then periodically used to update the distributions, as described in the
three_bond_bias_update_frequency section.
A value of 0 disables the periodic computation of these distributions.
- three_bond_bias_update_frequency (integer) optional parameter,
also see cbmc_setting_style
-
(default: 0)
The distributions used to perform the three bond biasing are updated with this step frequency. The updates combined the distributions
computed with a frequency controlled by the three_bond_bias_compute_frequency
parameter with the previous version of the three bond biasing potentials to create the new three bond biasing potential.
A value of 0 disables this update.
- three_bond_bias_old_fraction (double precision) optional parameter,
also see cbmc_setting_style
-
(default: 0.5)
This factor determines the linear combination of the old distribution, and the observed distribution, to determine the new three bond
bias distribution. This number must be in the range [0.0,1.0] inclusive. Setting this value to 0.0 would completely replace the old
distribution with the new distribution, while a setting of 0.5 would combine the old and observed distributions equally in order to
compute the new distribution.
Shared parameters for three_bond_fixed_endpoint_bias_style
- three_bond_bias_sdev_multiplier optional parameter,
also see cbmc_setting_style
-
(default: 1.0)
This value is multiplied by the standard deviation determined from the fit in order to get the standard deviation that is used
during the simulation. Must be positive.
-
Number of force field files listed in ff_filename.
-
A list of the filenames (one per line for ffnumber lines) that contain the force field information.
On most systems list this file (complete with directory path) and then end the line.
If you are having trouble reading this on an older system it is possible that appending blank spaces to get up to 70 characters solves your problem.
-
The setting for this parameter controls the nonbonded potential type. Depending on the setting there are a number of other parameters that are
also required.
Please see the classical_potential web page for more information.
- 'none' : no electrostatic potential
- 'coulomb' : Coulomb's law for the electrostatic potential. This option requires the following parameters.
- coulombstyle (character*20) optional parameter
- 'ewald_fixed_kmax' (default) point charges with an Ewald sum that utilizes a constant number of inverse space vectors
(kmax) and a parameter electrostatic cutoff (rcelect) equal to half the current box length.
This option requires the following parameters.
- kalp (double precision) optional parameter
-
(default: 5.6d0)
Value used in the Ewald sum to compute alpha. The internal Ewald sum alpha term is equal to kalp divided by the shortest
box length.
- kmax (integer) optional parameter
-
(default: 5)
Maximum number of inverse space vectors to use in any dimension for the Ewald sum.
Maximum allowed value is contrained by MAXKMAX.
- kalp (double precision) optional parameter
- 'ewald_fixed_cutoff' point charges with an Ewald sum that utilizes a constant electrostatic cutoff (rcelect)
and adjusts the number of inverse space vectors (kmax) according to the following heuristic.
- alpha = ( 1.35 - 0.15 log[ewald_prec]) / rcelect
- kmax = ( alpha * Max[box length] / π) * (-log[ewald_prec])0.5
This option enables the following parameters.- ewald_prec (double precision)
-
Controls the precision of the Ewald summation technique. The smaller the value, the better the results (and the more expensive the
simulation). The recommended value of 1d-4 is generally adequate, while a value of 1d-5 is very good (but more expensive).
- rcelect (double precision)
-
The cutoff for electrostatic interactions computed in the "real space" portion of the Ewald sum. Decreasing this value means less
work in the "real space", but correspondingly more work in the "inverse space". Setting this equal to the general nonbonded cutoff
(see rcut in classical_potential) is recommended.
- 'minimum image' use the minimum image convention to compute the coulombic interactions between all pairs of atoms in a system. Please note that this option is implemented mainly to compute single-molecule energies for test systems and is not suggested for routine use in periodic systems. For a discussion of the perils of using cut-off methods for Coulombic interactions please see Hummer et al. 1997.
- 'ewald_fixed_kmax' (default) point charges with an Ewald sum that utilizes a constant number of inverse space vectors
(kmax) and a parameter electrostatic cutoff (rcelect) equal to half the current box length.
This option requires the following parameters.
- dielect (double precision) optional parameter
-
(default: 1.0)
The dielectric constant used when computing coulombic interactions. Generally this should be set to 1.0 as the solvated system
will act as the screening that the dielectric constant is intended to represent. If you are performing a simulation without any
solvent (for example a protein without the water) you could set this value to the dielectric constant of the missing solvent.
- lcompute_madelung (logical) optional parameter
.FALSE.(default) do not compute the Madelung constant. This is the correct option in the vast majority of simulations..TRUE.compute the Madelung constant assuming the system is a perfect crystal. This is done to check that the electrostatic potential is working properly by comparing against systems where the Madelung constant is known analytically.
- coulombstyle (character*20) optional parameter
-
(default: 0)
The number of external fields to apply in the simulation. These fields can take on a variety of forms, but are always applied relative to a
plane in one of the simulation boxes. Typical uses are for simulating the effect of a rigid surface without having to treat the surface atoms
explicitly. If this value is not specified then it is set to the default value of 0.
When there is a nonzero value of nfield, the following parameters are required for for each field.
- fieldtype (character*20)
- 'Hard Ball' : Places a hard ball of a specified diameter at a particular location in one of the boxes.
This ball interacts with each atom in the simulation and excludes all atom centers from its radius.
When using this option you must also specify the following parameters.
- hrbbox (integer)
-
The simulation box that contains this field. Allowable values range from 1 to numboxes.
- hrbpos (double precision)
-
Position of the center of the hard ball. Must be between 0.0 and the box length of each of the 3 axes. All three coordinates
are listed on a single line.
- hrbrad (double precision)
-
Radius of the hard ball. The ball does interact through the periodic boundaries.
- hrbfac (double precision)
-
Factor used to rescale the size of the hard ball. A value of 1.0 will use the hrbrad value without any adjustment.
This factor is multiplied by the hrbrad value to compute the true radius of the hard ball.
- hrb_energy_type (character*11)
- 'infinite' : hard overlap with the ball of infinite energy
- 'finite' : a specfic finite penalty for each atom that is inside of the hard ball.
Using this option requires the following parameter.
- hrb_ball_energy (double precision)
-
The energy penalty (units of K/kB) for each atom that overlaps with a hard ball.
- hrb_ball_energy (double precision)
- hrbbox (integer)
- 'Hard Wall' : Places a hard wall of a specified diameter in one of the boxes. This wall interacts with each atom in the simulation.
This option requires the following parameters.
- hrdbox (integer)
-
The simulation box that contains this field. Allowable values range from 1 to numboxes.
- hrdxyz (character*1)
- 'x' : hard wall is perpendicular to the x-axis (in the yz plane)
- 'y' : hard wall is perpendicular to the y-axis (in the xz plane)
- 'z' : hard wall is perpendicular to the z-axis (in the xy plane)
- hrdcen (double precision)
-
Position of the center of the hard wall. Must be between 0.0 and the box length of the axis that is perpendicular to the wall.
- hrdrad (double precision)
-
Radius of the hard wall. The wall interacts through the periodic boundaries.
- hrd_repulsion_style (character*11)
- 'centers' : the center of each atom in the system is the portion that interacts with the hard wall region. In other words, the hard wall interacts with the atoms as if they were point particles.
- 'hard radii' : the atoms are considered to have a radius that is determined by their nonbonded parameters. The distance between the surface of the hard wall and the closest approach of the atomic radius is used to determine the interactions between the wall and the atoms. This option is not currently supported for all choices of the classical_potential.
- hrd_energy_type (character*11)
- 'infinite' : any molecule inside of the hard wall has an infinite energy (hard overlap).
- 'finite' : any molecule inside of the hard wall has a finite energy that is specified by the hrd_wall_energy parameter.
If you use this option you must include the following parameter
- hrd_wall_energy (double precision)
-
The energy given to any atom that is inside of the hard wall (in Kelvin/kB). This option is designed to enable equilibrium of a hard
wall system as this provides an incentive for molecules to leave the hard wall area without causing a simulation ending overlap.
- hrd_wall_energy (double precision)
- hrdbox (integer)
- 'Harmonic Attractor' : Uses a harmonic potential to root certain atoms to a defined point or their initial positions.
- hafbox (integer)
-
The simulation box that contains this field. Allowable values range from 1 to numboxes.
- hafk (double precision)
-
The force constant for the harmonic potential.
- hafentries (integer)
-
The number of types or elements to which this field is applied.
- hafrefpos (character*7)
- 'Global' : use a global set of coordinates for each atom. This option requires the following parameter.
- hafglobxyz (double precision array)
-
The x,y, and z coordinates of the global position. This is listed as three numbers on the same line separated by spaces.
- hafglobxyz (double precision array)
- 'Initial' : use the initial coordinates of each atom.
- 'Global' : use a global set of coordinates for each atom. This option requires the following parameter.
- hafkey (character*7)
- 'Element' : apply this field to a specific group of atoms that based upon their element.
This option requires the following parameters.
- hafmolec (integer)
-
The molecule type for this field.
- hafelement (character*2)
-
The element type for this field.
- hafmolec (integer)
- 'FFtype' : apply this field to a specific group of atoms that are all the same molecule and nonbond type.
This option requires the following parameters.
- hafmolec (integer)
-
The molecule type for this field.
- hafname (character*10)
-
The nonbond type for this field.
- hafmolec (integer)
- 'Element' : apply this field to a specific group of atoms that based upon their element.
This option requires the following parameters.
- hafbox (integer)
- 'Hooper Umbrella' : Place a Hooper Umbrella field (see Hooper et al. 2000) in a
simulation box. This is a 4th power energy function based on displacement along a single axis.
- v(d) = umba * [ (d - umbcenter) / umbcenter ]4
With this option you must also specify the following parameters.- umbbox (integer)
-
The simulation box that contains this field. Allowable values range from 1 to numboxes.
- umbxyz (character*1)
- 'x' : Field is perpendicular to the x-axis (in the yz plane)
- 'y' : Field is perpendicular to the y-axis (in the xz plane)
- 'z' : Field is perpendicular to the z-axis (in the xy plane)
- umbcenter (double precision)
-
The zero energy point of the field, listed as a distance in Angstroms along the axis specified in umbxyz.
- umba (double precision)
-
The energy constant in units of K/kB.
- 'LJ 9-3 Wall' : Places a 9-3 Lennard-Jones wall in one of the simulation boxes. The wall potential has the following functional form.
- v(d) = [ 2/3 π εwf σwf3 ρwall ] * [ 2/15 (σwf/d)9 - (σwf/d)3 ]
With this option you must also specify the following parameters.- ljfbox (integer)
-
The simulation box that contains this field. Allowable values range from 1 to numboxes.
- ljfxyz (character*1)
- 'x' : Lennard-Jones wall is perpendicular to the x-axis (in the yz plane)
- 'y' : Lennard-Jones wall is perpendicular to the y-axis (in the xz plane)
- 'z' : Lennard-Jones wall is perpendicular to the z-axis (in the xy plane)
- ljfcen (double precision)
-
Position of the center of the Lennard-Jones wall. Must be between 0.0 and the box length of the axis that is perpendicular to the wall.
- ljfdir (integer)
- -1: Atoms only interact with the "left" face of this wall. This does not extend through the periodic boundary.
- 1: Atoms only interact with the "right" face of this wall. This does not extend through the periodic boundary.
- ljfcut (double precision)
-
The distance beyond which the wall-atom interactions are not computed and assumed to be zero.
- ljfshift (logical)
.TRUE.: shift the lj wall potential to be zero at the cutoff..FALSE.: do not shift the potential.
- ljfrho (double precision)
-
The number density of atoms in the integrated wall potential (units of atoms per cubic Angstrom).
- ljfntypes (integer)
-
The number of atom types in the system that interact with the wall. Any atom type not listed here will have zero interaction with
the wall. For each value of ljfntypes you must list the following parameters.
- ljfname (character*6)
-
The nonbonded name of the atom. This is compared with the nonbonded atom names that are either listed in the appropriate
input_style option for each molecule, or those directly listed in the appropriate towhee_ff_* files.
- ljfsig (double precision)
-
σ parameter for the interaction between this atom and the wall atoms. Units are Angstroms.
- ljfeps (double precision)
-
ε parameter for the interaction between this atom and the wall atoms. Units are K/kB.
- ljfname (character*6)
- 'Steele Wall' : Places a 10-4 Lennard-Jones wall in one of the simulation boxes. The wall potential has the following functional form.
- v(z) = εw [ 2/5 (σsf/z)10 - (σsf/z)4 - σsf4 / [ 3 Δ ( z + 0.61 Δ )3 ] ]
where- εw = 2 π εsf ρs σsf2 Δ
This potential is attributed to Steele 1973, but I found that reference a bit confusing so I implemented the equations as presented in Lastoskie et al. 1993 and the parameter names here follow their notation.- With this option you must also specify the following parameters.
- steele box (integer)
-
The simulation box that contains this field. Allowable values range from 1 to numboxes.
- steele xyz (character*1)
- 'x' : wall is perpendicular to the x-axis (in the yz plane)
- 'y' : wall is perpendicular to the y-axis (in the xz plane)
- 'z' : wall is perpendicular to the z-axis (in the xy plane)
- steele surface (double precision)
-
Position of the surface of the wall. Must be between 0.0 and the box length of the axis that is perpendicular to the wall.
- steele dir (integer)
- -1: Atoms only interact with the "left/bottom" face of this wall. This does not extend through the periodic boundary.
- 1: Atoms only interact with the "right/top" face of this wall. This does not extend through the periodic boundary.
- steele cutoff (double precision)
-
The distance beyond which the wall-atom interactions are not computed and assumed to be zero.
- steele shift (logical)
.TRUE.: shift the wall potential to be zero at the cutoff..FALSE.: do not shift the potential.
- steele delta (double precision)
-
Δ parameter for the spacing between the layers in the solid. Units are in Angstroms.
- steele rho_s (double precision)
-
ρs parameter for the density of the atom in the solid. Units are in atoms per cubic Angstrom.
- steele ntype (integer)
-
The number of atom types in the system that interact with the wall. Any atom type not listed here will have zero interaction with the
wall. For each type you must list the following parameters.
- steele name (character*6)
-
The nonbonded name of the atom. This is compared with the nonbonded atom names that are either listed in the appropriate
input_style option for each molecule, or those directly listed in the appropriate towhee_ff_* files.
- sigma_sf (double precision)
-
σ parameter for the interaction between this atom and the wall atoms. Units are Angstroms.
- epsilon_sf (double precision)
-
ε parameter for the interaction between this atom and the wall atoms. Units are K/kB.
- steele name (character*6)
- 'Hard Ball' : Places a hard ball of a specified diameter at a particular location in one of the boxes.
This ball interacts with each atom in the simulation and excludes all atom centers from its radius.
When using this option you must also specify the following parameters.
- 'none' : (default) no additional solvation potential.
- 'internal' : use one of the solvation potentials internal to Towhee. This option requires the following parameter.
- solvation_type (character*20)
- 'EEF1' : solvation using the Charmm19-EEF1 potential.
- solvation_type (character*20)
- 'external' : use a solvation potential via an external program that has been compiled as a library and linked into Towhee.
- solvation_type (character*20)
- 'tramonto' : solvation using the classical density functional theory code Tramonto to compute a solvation free energy.
- solvation_type (character*20)
- 'advanced connectivity map' : advanced connectivity map is a powerful, general way to set up molecules
- 'basic connectivity map' : basic connectivity map
- 'explicit' : explicit declaration of all terms.
- 'nanotube builder' : nanotube builder designed for carbon nanotubes
- 'nucleic acid builder' : nucleic acid builder designed for DNA
- 'polypeptide builder' : polypeptide builder designed for proteins