SDA7 input file format

The input file consists of a series of grouped input parameters, delimited by the keyword GROUP

GROUP = groupname
...
END GROUP

and a list of free parameters which are not associated with a group.

All input parameters have the following form :

name-of-parameter = parameter-value

Lines begining with the character "#" or "!" are treated as comment lines and will not be parsed by the program.
Note: do not insert spaces in front of "#" if you want to comment out a line!

Example input files

The tabs below contain different examples of input files. All parameter names are linked to a short summary. A full list and description of the parameters is given below the tabs.

Select any other tab

Example input file for calculating barnase-barstar association rates.
This setup will record encounter complexes, one trajectory and the residence time for different center-to-center distances between the solutes.

GROUP = Type_Calculation
# type of simulation : sda_2proteins / sda_koff / sda_energy / sda_pmf / sdamm
type = sda_2proteins

# total number of solutes
total_solutes = 2

# total number of grids (number of groups "Solute_Grid")
total_grids = 2
END GROUP

# random number generator, 0 uses system clock
dseed = 256.0

# number of trajectories to perform
nrun = 1000

# maximum time (in ps) for each trajectory. 0 is no-limit
timemax = 0.

# size of the protein probe radius (Angstroms)
probep = 1.77

# size of the water probe radius (Angstroms)
probew = 1.40

# size of the spacing of the exclusion grid (Angstroms)
hexcl = 0.5

# threshold to determine the atoms on the surface
threshold = 0.0

# for saving the exclusion grid or the list of atoms on the surface
save_exclusion = 1
save_access = 1

# Ignore excluded volume check that prevents solute overlap. When turning this off the soft-core repulsion potential can be used instead.
ignore_exclusion = 0

# multiplicative factor for the electrostatic grid (always 0.5 for sda 7, 1.0 for previous version)
epfct = 0.5
# multiplicative factor for the electrostatic desolvation grid
edfct = 1.67
# multiplicative factor for the nonpolar grid
# hdfct = -0.013
hdfct = 0.0
# multiplicative factor for the Lennard-Jones grid
ljfct = 0.0
# multiplicative factor for the soft-core repulsion grid
lj_rep_fct = 0.

GROUP = ReactionCriteria
rxna12f = p12.rxna
computation = association
dind = 6.0
nnnons =2
nwrec = 23
END GROUP

GROUP = RateCalculation
# describe the reaction windows
win0 = 3.0
nwin = 35
dwin = 0.5
nb_contact = 4
# optional: output results of every trajectory for a post-bootstrap analysis
bootstrap = 1
# optional: compute first passage time, bootstrap is appled in this case too
fpt = 0
# optional: stop trajectories when the most strict reaction criteria is met (i.e. closest distance window and highest number of contacts.
stop_traj = 1
# optional: add an analytical correction to the Northrup, Allison and McCammon Method to account for long range interactions at the b surface.
analytical_correction = 1
END GROUP

# Specify solute 1:
GROUP = Solute_Grid
nb_solute = 1
pdb_filename = p1_noh.pdb
diffusion_trans = 0.027
diffusion_rotat = 3.92e-5
real_net_charge = 5.0
rotate = 1
surface = 0
flex = 0
epf = ep1.grd
qef = p1.echa
edf = ed1.grd
hdf = hd1.grd
#nb_lj = 3
#ljf = p1_lj.grd
#lj_repf = p1_lj.grd
END GROUP

# Specify solute 2
GROUP = Solute_Grid
nb_solute = 1
pdb_filename = p2_noh.pdb
diffusion_trans = 0.027
diffusion_rotat = 3.92e-5
real_net_charge = -6.0
rotate = 1
epf = ep2.grd
qef = p2.echa
edf = ed2.grd
hdf = hd2.grd
#nb_lj = 3
#ljf = p2_lj.grd
#lj_repf = p2_lj.grd
END GROUP

# Define the geometry of the simulation box
GROUP = Geometry
# define the geometry : sphere / box
type = sphere
# if periodic boundary condition are applied (only for box)
pbc = 0
# if the protein can escape
escape = 1
# if a surface is present, influence PDB
surface = 0
### deprecated, should use only start_position
##b 150.
# c-surface distance, only with sphere (Angstroms)
c = 300.
# b-surface distance in sphere, z distance within a box (Angstroms)
start_pos = 150.
END GROUP

GROUP = ResidenceTime
# type : distance / plan / resid_3d
type_residence = distance
# used to specify the size or the number of cells. If 2, specify both
fixed_bin_size = 0
# if fixed_bin_size = 1 or 2 : impose the size of the bin, number of cells computed from geometry
size_x = 2.0
# size_y = 2.0
# size_z = 2.0
# if fixed_bin_size = 1 or 2 : impose the size of the bin, number of cells computed from geometry
cell_x = 100
# cell_y = 100
# cell_z = 100
END GROUP

GROUP = Timestep
# the timestep is variable in sda_2proteins (and derivatives) and fixed for sdamm
variable = 1
# if variable, parameters to define the distance dependence
dt1 = 1.0
swd1 = 50
dt2 = 20.0
swd2 = 90.0
END GROUP

GROUP = Complexes
# parameters for file containing encounter complexes
# name of the output file, comment out to cancel the creation of a complexes file
fcomplexes = complexes
# 0: ascii; 1: binary
binary_complex = 0
# if 1, sum the energy terms
complex_sum_energy = 0
# maximum number of complexes
nb_complexes = 500
# maximum rmsd for replacing a stored complex by another
rmsd_min = 1.0
# if 1, the rmsd is only compared for complexes in the same trajectory.0 is recommended, compare a complex with all complexes from all trajectories
ionerun = 0

# use only with OpenMP. If commented out, the default values are applied (recommended)
# when a local complexes list is merged into the global one
merge_step = 25
#Size of the local complexes list
size_thread = 300

# parameters for the trajectory file
ftrajectories = trajectories
binary_trajectory = 0
trajectory_sum_energy = 0
# if < 0, record all runs; if 0, do not record anything; if > 0, give number of specific trajectory to record
ntraj_rec = 1
freq_print = 100
max_array_traj = 0
mem_traj_frames = 1
END GROUP

rboost = 1.0
novers = 150

Input file for docking a solute protein to a metal surface using the ProMetCS force field.

GROUP = Type_Calculation
# type of simulation : sda_2proteins / sda_koff / sda_energy /sda_pmf / sdamm
type = sda_2proteins

# total number of solutes
total_solutes = 2

# total number of grids (number of groups "Solute_Grid")
total_grids = 2
END GROUP

# random number generator, 0 uses system clock
dseed = 256.0

# number of trajectories to perform
nrun = 1000

# maximum time (in ps) for each trajectory. 0 is no-limit
timemax = 5e+04.

# size of the protein probe radius (Angstroms)
probep = 0.5

# size of the water probe radius (Angstroms)
probew = 1.40

# size of the spacing of the exclusion grid (Angstroms)
hexcl = 0.5

# threshold to determine the atoms on the surface
threshold = 0.0

# for saving the exclusion grid or the list of atoms on the surface
save_exclusion = 1
save_access = 1

# multiplicative factor for the electrostatic grid ( always 0.5 for sda 7, 1.0 for previous version)
epfct = 0.5
# multiplicative factor for the electrostatic desolvation grid
# in ProMetCS this factor is 2x that for two solute (sda_2proteins) calculations when the image charge model for electrostatics is used !
edfct = 3.34
# multiplicative factor for the nonpolar grid
hdfct = -0.0065
# multiplicative factor for the Lennard-Jones grid
ljfct = 1.0

# global variables, specific to ProMetCS
# scale_epsilon for image-charge correction
correction_image_charge = 1
# activates desolvation image charge
desolv_image_charge = 1

# define all parameters for the metal desolvation term
# these are default values and only need to be specified to override these.
# declaration of the group (even empty) is sufficient to get the default values.
GROUP = MetalDesolvation
distance_to_surface = 6.5
energy_per_water = 0.13
radius_patch = 6.0
radius_water = 1.5
END GROUP

GROUP = ReactionCriteria
rxna12f = SurfaceDock.rxna
computation = docking
nnnons =0
nwrec = 10
END GROUP

# Define the surface as solute 1
GROUP = Solute_Grid
nb_solute = 1
pdb_filename = gold3layersflex.pdb
diffusion_trans = 0.
diffusion_rotat = 0.
real_net_charge = 5.0
rotate = 0
surface = 1
flex = 0

edf = ed1.grd
hdf = hd1.grd
nb_lj = 2
# very important: only one blank space between atom types
atom_lj = Au Av
END GROUP

# Define the protein as solute 2
GROUP = Solute_Grid
nb_solute = 1
pdb_filename = p2_noh.pdb
diffusion_trans = 0.0123
diffusion_rotat = 1.36e-4
real_net_charge = -6.0
rotate = 1

# set image charge calculation for these solutes
image_charge = 1
epf = ep2.grd
qef = p2.echa
edf = ed2.grd

nb_lj = 2
ljf = p2lj1.grd p2lj2.grd
END GROUP

# Define the geometry of the simulation box
GROUP = Geometry
# define the geometry : sphere / box
type = box
# if periodic boundary conditions are applied (only for box)
pbc = 0
# if the protein can escape
escape = 1
# if a surface is present, influence PDB
surface = 1
xmin = 40.
xmax = 60.
ymin = 40.
ymax = 60.
zmin = 0.
zmax = 140.
# b-surface distance in sphere, z distance in a box
start_pos = 70.
END GROUP

GROUP = Timestep
# the timestep is variable in sda_2proteins (and derivatives) and fixed for sdamm
variable = 1
# if variable, parameters to define the distance dependence
dt1 = 0.12
swd1 = 50
dt2 = 20.0
swd2 = 90.0
END GROUP

GROUP = Complexes
# parameters for the complexes file
# name of the output file, comment out to cancel the creation of a complexes file
fcomplexes = complexes
# if 0:ascii; 1: binary
binary_complex = 0
# if 1, sum the energy terms
complex_sum_energy = 0
# maximum number of complexes
nb_complexes = 500
# maximum rmsd for a replacing a stored complex with another complex
rmsd_min = 1.0
# if 1, the rmsd is only compared with complexes in the same trajectory (0 recommended, compare with all complexes)
ionerun = 0

# use only with OpenMP. If commented out, the default values are applied (recommended)
# When a local complexes list is merged into the global one
merge_step = 25
#Size of the local complexes list
size_thread = 300

# parameters for the trajectory file
ftrajectories = trajectories
binary_trajectory = 0
trajectory_sum_energy = 0
# if < 0, record all runs; if 0, do not record anythingr;, if > 0, give number of trajectory to record
ntraj_rec = 1
freq_print = 100
max_array_traj = 0
mem_traj_frames = 1
END GROUP

rboost = 1.0
novers = 150

This sdamm example input file is similar to that in the example found in examples/lysozymes/128_flex/.
It specifies the simulation of 128 lysozyme molecules. The interaction grids for lysozyme have been computed at three different pH values (3, 6 and 9). It demonstrates that grids with different types of variation (here, protonation state) can be used for each solute.

GROUP = Type_Calculation
# type of simulation : sda_2proteins / sda_koff / sda_energy /sda_pmf / sdamm
type = sdamm

# total number of solutes
total_solutes = 128

# total number of grids (number of groups "Solute_Grid")
total_grids = 1
END GROUP

# random number generator, 0 uses system clock
dseed = 256.0

# number of trajectories to perform, always 1 with sdamm
nrun = 1

# maximum time (in ps) for the trajectory
timemax = 10000.0

# size of the protein probe radius (Angstroms)
probep = 1.77

# size of the water probe radius (Angstroms)
probew = 1.40

# threshold to determine the atoms on the surface
threshold = 0.0

# multiplicative factor for the electrostatic grid ( always 0.5 for sda 7, 1.0 for previous version)
epfct = 0.5
# multiplicative factor for the electrostatic desolvation grid computed with ionic strength set to zero
edfct = 0.36
# multiplicative factor for the nonpolar grid
hdfct = -0.013
# multiplicative factor for the Lennard-Jones grid
ljfct = 0.0
# multiplicative factor for the soft-core repulsion grid
lj_rep_fct = 0.0156

# name of the restart file, initial solute positions must be specified with sdamm
restart = restart_file

# Apply to analytical interactions
GROUP = Analytical
# Define bin size (A) of the precomputed array
h_analytic = 0.01
# Activates Debye-Huckel, for all solutes
debyeh = 1
# Defines ionic strength (in Molar)
ionic_strength = 0.100
# Activates Hydrodynamic interactions, boolean: 0/1
hydrodynamic = 1
# Radius used to define local volume for Hydrodynamic interactions (Angstroms)
lvol_rcut = 80.0
END GROUP

# 128 lysozymes, only one type of solute but it has three different structures/interactions
GROUP = Solute_Grid
nb_solute = 128
diffusion_trans = 0.01096
diffusion_rotat = 1.9e-5
real_net_charge = -6.0
rotate = 1
surface = 0
flex = 1

# stoke radius if hydrodynamics interaction and/or Debye-Huckel are activated
stoke_radius = 0.

# the following parameters apply only if flexibility is activated
# the number of conformations/structures per solute
total_conf = 3
# the method to use when choosing a new structre: metropolis, minimum, random, time
method = minimum
# the conformation/structure to choose as the initial one; -1 random ( sdamm reads from the restart file)
initial_conf = -1
# time (in ps) between trials to change conformation/structure
frequency = 500.0
std_frequency = 50.0
# filename with the list of conformations/structures (pdbs and grids)
list_conf = list_conformation.lys.rescaled.dat
END GROUP

# define the geometry of the simulation box
GROUP = Geometry
# define the geometry : sphere / box
type = box
# if periodic boundary conditions are applied (only for box)
pbc = 1
# if the protein can escape
escape = 0
# if a surface is present, influence PDB
surface = 0
# dimensions of the simulation box
xmin = 0.
xmax = 314.965
ymin = 0.
ymax = 314.965
zmin = 0.
zmax = 314.965
END GROUP

GROUP = Timestep
# the timestep is variable in sda_2proteins (and derivatives) and fixed for sdamm
variable = 0
# if variable, give parameters to define the distance dependence
dt1 = 0.4
END GROUP

GROUP = Complexes
# parameters for the trajectory file
ftrajectories = trajectories
binary_trajectory = 0
trajectory_sum_energy = 0
# if < 0, record all runs; if 0, do not record anything. If > 0, the trajectory number to record
ntraj_rec = 1
freq_print = 1000
max_array_traj = 0
mem_traj_frames = 1
END GROUP

Input file for the Potential of Mean Force (PMF) calculation of a 3 histidine chain with a metal surface using the ProMetCS force field.

GROUP = Type_Calculation
# type of simulation : sda_2proteins / sda_koff / sda_energy /sda_pmf / sdamm
type = sda_pmf

# total number of solutes
total_solutes = 2

# total number of grids (number of groups "Solute_Grid")
total_grids = 2
END GROUP

# size of the protein probe radius (Angstroms)
probep = 0.5

# size of the water probe radius (Angstroms)
probew = 1.40

# size of the spacing of the exclusion grid (Angstroms)
hexcl = 0.5

# threshold to determine the atoms on the surface
threshold = 0.0

# for saving the exclusion grid or the list of atoms on the surface
save_exclusion = 1
save_access = 1

# multiplicative factor for the electrostatic grid ( always 0.5 for sda 7, 1.0 for previous version)
epfct = 0.5
# multiplicative factor for the electrostatic desolvation grid
# in ProMetCS this factor is 2x that for two solute (sda_2proteins) calculations when the image charge model for electrostatics is used !
edfct = 3.34
# multiplicative factor for the nonpolar grid
hdfct = -0.0065
# multiplicative factor for the Lennard-Jones grid
ljfct = 1.0

# global variables, specific to ProMetCS
# scale_epsilon for image-charge correction
correction_image_charge = 1
# activates desolvation image charge
desolv_image_charge = 1

# define all parameters for the metal desolvation term
# these are default values and only need to be specified to override these.
# declaration of the group (even empty) is sufficient to get the default values.
GROUP = MetalDesolvation
distance_to_surface = 6.5
energy_per_water = 0.13
radius_patch = 6.0
radius_water = 1.5
END GROUP

# Define the surface as solute 1
GROUP = Solute_Grid
nb_solute = 1
pdb_filename = gold3layersflex.pdb
rotate = 0
surface = 1

edf = ed1.grd
hdf = hd1.grd
nb_lj = 2
# very important: only one blank space between atom types
atom_lj = Au Av
END GROUP

# Define the protein as solute 2
GROUP = Solute_Grid
nb_solute = 1
pdb_filename = p2_noh.pdb
rotate = 1

# set image charge calculation for these solutes
image_charge = 1
epf = ep2.grd
qef = p2.echa
edf = ed2.grd

nb_lj = 2
ljf = p2lj1.grd p2lj2.grd
END GROUP

# Define the geometry of the simulation box
GROUP = Geometry
# define the geometry : sphere / box
type = box
surface = 1
END GROUP

GROUP = Complexes
# parameters for the trajectories file
# name of the output file, comment out to cancel the creation of a trajectories file
ftrajectories = trajectories
END GROUP

GROUP = PMF
# Number of angles (for each of the Euler angles) for the sampling of the rotation of the rigid body structure
mentth = 12
mentfi = 24
mentom = 24
# Rectangular area of the surface to be covered
mentx = 6
menty = 6
# distance from the center of geometry of the protein to the center of the surface in z direction
pstart = 6
pfinish = 6
# incremental step along the z axis
dz = 0.1
END GROUP

In the case that a solute has more than one structure/conformation, the input pdb and grid files must be indicated in a separate file.
One file is required for each flexible SoluteGrid group.

More examples of the use of the IMPLICIT/EXPLICIT keyword are in Tutorial_SDA7/Importin_Ran/docking_flex

# Need to specify one of the 2 keywords: IMPLICIT/EXPLICIT
# IMPLICIT keyword: the filename is deduced from a template
# match the first "." beginning from the end of the string and replace X.suf by N.suf, with N incremented from first_mode to total_conf
# IMPLICIT, X.suffix will be replaced by -1.suffix, 0.suffix, 1.suffix...
IMPLICIT
first_mode = -1
# Can use all the same keywords as in a normal SoluteGrid group
GROUP = Conf1
pdb_filename = ../data_grid/mode-7-X.pdb_pqr
epf = ../data_grid/mode-7-X.grd
qef = ../data_grid/mode-7-X.echa_R
END GROUP

The atomic radii and test charges of specified atom types can be defined by the user in this file.
The file must be present in the directory where sda_flex is executed and must be named "add_atoms".
An example is provided in the src/add_atoms directory.

# List new atom types, assign radii and test charges.
# Can also override the default values.
#
# The specified types are read BEFORE the default ones (hard-coded in mod_van.f90) so they will override them.
# As soon a matchng atom type is found, the search ends.
# The order of the atom types is therefore important.
# Residue names are used for specifying charges but not for specifying radii, use * or empty string for the residue name for specifying radii.
#
### GROUP = VdW
# To add cobalt with a 1.6 Angtroms radius : (and avoid confusion with carbon!)
# CO * 1.6
# To change the default values of all carbons to 2.5 Angstroms:
# C* * 2.5
# END GROUP #
### GROUP = Test Charges
# To override all histidines, (solve problem HI2 - HIP with whatif )
# N* HIS 0.0
### N* HI* 0.5 # assigns this charge to backbone and sidechain N atoms in histidine, not recommended
# ND HI* 0.5 # overrride only ND in His
# NE HI* 0.5 # override only NE in His
# END GROUP

# The VdW radii must be indicated in the GROUP VdW
# The VdW parameters have a fixed format: 2 char, 3 spaces, 3 char., 1 space, float f6.2
GROUP = VdW
FE   HEM   3.00
MG     *   1.20
CO   FMN   4.35
END GROUP

# The values of the test charges must be indicated in the group TestCharges
# The test charges have a fixed format: 3 char, 2 spaces, 3 char., 1 space, float f6.2
GROUP = TestCharges
O1P  FMN  -1.00
O2P  FMN  -1.00
O1P  FAD  -1.00
O1A  FAD  -1.00
O1A  NAP  -1.00
O1N  NAP  -1.00
O1X  NAP  -1.00
O1A  HEM  -0.50
O2A  HEM  -0.50
O1D  HEM  -0.50
O2D  HEM  -0.50
END GROUP


Detailed description of the input parameters

GROUP = Type_Calculation:  type, total_solutes, total_grids

GROUP = Reaction Criteria:   computation, rxna12f, et_sol1, et_sol2, dind, nnnons, nwrec

Reaction criteria are used to define a "complex" between 2 solutes. Full freedom is let to the user:

The list of reaction criteria must be provided in a separate file (conventionally named with a suffix *.rxna), where the fields indicate:

For docking, there are 3 different types of reaction you can use:

For association rate calculations, all entries in the *.rxna file are considered as CNONS. But you can still record complexes and adjust the nwrec input parameter.

To get the specific format of the *.rxna file, we advise you to look at the input files in the barnase-barstar example of docking. Here are illustrated:

Or you can refer to the technical documentation, linked here in the section Reaction Criteria.
(The doxygen documentation must be generated before using this link).

Usage of dind
To avoid an artificially large number of contacts, the variable dind can be used.
On the figure below, all atoms used in the reaction criteria are represented by spheres. But only those which are at a separation greater than d_min are considered as independent (black spheres). So with this schema, there will be a maximum of 4 independent contacts.

graphics

GROUP = RateCalculation:  win0, nwin, dwin, nb_contact, bootstrap, fpt, stop_traj,analytical_correction

Calculate bimolecular association rate constants using the Northrup, Allison and McCammon method (see Northrup, Allison McCammon, (1984), J. Chem. Phys.).

GROUP = Analytical:  h_analytic, debyeh, ionic_strength, hydrodynamic

GROUP = Solute_Grid: nb_solute, pdb_filename, diffusion_trans, diffusion_rotat, real_net_charge, rotate, surface, flex, epf, qef, edf, hdf, ljf, lj_repf, nb_lj, atom_lj, list_conformation, stoke_radius, total_conf, method, initial_conf, frequency, std_frequency

GROUP = Geometry:   type, pbc, surface, escape, start_position, c_surface, min_X, max_X

Typical setup for simulations within a sphere. Trajectories are started at the b-surface (start_pos) and finish when solute 2 goes through the c-surface (c).

graphics2

GROUP = Timestep:  variable,dt1, swd1, dt2, swd2

These parameters define the variable timestep designed to accelerate simulations. As shown in the figure below, the timestep is automatically decreased to zero near the c-surface to reduce artifacts due to truncation of the trajectories. Although the profile of the time step influences the computed rates, the most important parameter is dt1. Take care not to choose a too small value for this parameter (since the computing time scales almost linearly with it). The physical parameter sqrt(dt*dm) measures an average Brownian dynamics step length. This must be less than 1 Å, which is the characteristic distance describing the roughness of any protein surface since the atomic bond lengths and van der Waals radii define this scale. Moreover, dt1 should be less than the characteristic distance over which the interaction forces vary substantially, which means that the optimal value is influenced by interaction force strength: one needs smaller time steps for stronger interacting proteins. The standard output will print additional information, like the location and the value of the maximum timestep.

graphics1

GROUP = ResidenceTime:  type_residence, fixed_bin_size, size_X, cell_X,shift_rt_X,format_rt3d,pair_solutes

GROUP = Complexes: fcomplexes,restart_complex,binary_complex, nb_complexes, rmsd_min, ionerun,complex_sum_energy merge_step, size_thread, ftrajectories,binary_trajectory, ntraj_rec, freq_print, trajectory_sum_energy,max_array_traj, mem_traj_frames,max_array_traj

GROUP = MetalDesolvation:  distance_to_surface, energy_per_water, radius_patch, radius_water

The surface desolvation energy is calculated as :

graphics9

The computed desolvated area is shown by bold lines, the solid and dashed lines corresponding to water desorption from the first and second hydration layers, respectively. The circles with crosses represent close solute atoms and that with a zero is at an intermediate distance.

graphics10


GROUP = PMF:  mentth, mentfi, mentom, mentx, menty, pstart, pfinish, dz

graphics111


Parameters not associated with a group, they apply to all solutes


[Back to Main documetation]

Imprint/Privacy