BD_BOX manual, version 2.2, Rigid
Contents
1. ionic strength definition of boundaries values of cutoffs for intermolecular interactions specifications of the output names and types of files write frequency path to the input structure and diffusion tensor files numbers of molecules of different types in the simulated system and many others 7 RUNNING SIMULATIONS WITH BD_RIGID 20 While molecules are usually placed by the program at random positions and with random orientations inside the primary simulation cell it is also possible to define initial positions and orientations of all objects in the simulated system using a text formatted file Figure 3 The output from bd_rigid consists of a few files Figure 3 The output file out err contains the track of the BD run error messages if any and debugging information Another file in the PQR format pqr contains the definition of the molecular system under study it can be used to verify the starting configuration of molecules and also together with a trajectory file to visualize results The trajectory file either in the binary DCD or binary MOLB format dcd molb describes the time dependent behavior of the studied system MOLB trajectory contains only translations and rotations that are needed to perform transformations from the molecule fixed to the laboratory frame while DCD files contain Cartesian coordinates of atoms of all molecules in the simulated system in the laboratory frame The energy file e2. Models Within the BD_BOX rigid body framework molecules can be described either using a coarse grained representation or with fully atomistic details In both cases molecules are described as conglomerates of spherical subunits i e atoms Each kind of interaction that the molecule is involved in is assigned its own representation for example for a given molecule a set of charged spheres can be used to evaluate electrostatic interactions and a different set considering both positions and radii of spheres can be used to evaluate excluded volume interactions with other molecules in the studied system consult Figure 1 for the overall concept sphere dielectric cavity sphere excluded volume point charge Figure 1 Exemplary molecular models that are tractable within the BD_BOX rigid body framework Models consist of spheres with radii defining excluded volume interactions left Point charges can be positioned inside molecules their positions do not necessarily coincide with centers of spheres used to define excluded volume interactions left Additionally each charge can be assigned a non zero radius center Charge radii are used to define the dielectric cavity occupied by the molecule Positions and radii of spheres in the dielectric representation and in the excluded volume representation can be the same as for example in case of atomistically detailed molecular models a molecule is then represented by a set of van der Waa
3. bodies the bd_rigid module of the BD_BOX package When flexible models are employed each molecule consists of a various number of spherical subunits beads connected with deformable bonds 1 2 Bonded interactions that result in deformations of planar and dihedral angles can be also modeled Direct nonbonded interactions between molecules are evaluated using pairwise functions describing screened electrostatics in dielectric media 3 with effective charges assigned to spherical subunits 4 and Lennard Jones potential types The far field hydrodynamic effects are modeled using the config uration dependent Rotne Prager Yamakawa mobility tensor 5 6 7 or its Ewald summed form in case of periodic systems 9 Equations of motions are propagated using either the Ermak McCammon scheme or the predictor corrector IG T algorithm by Iniesta and Garcia de la Torre 11 Hydrodynamically correlated random displacements are generated either via the Cholesky factorization of the configuration dependent mobility tensor matrix 10 using the TEA HI approach proposed by Geyer and Winter 12 13 or with an approach described by Ando and Skolnick that utilizes Krylov subspaces With BD_BOX one can simulate flexible molecules in homogeneous flows 2 or external electric fields direct alternate or rotating fields With BD_BOX one can also simulate rigid bodies that are described with fully anisotropic diffusion tensors 115 Molecules treated as rigid bod
4. by the program to remove overlaps in the system in each attempt a BD step that leads to an overlap is repeated with a different random vector Equation 3 8 CONTROL FILE KEYWORDS AND COMMAND LINE PARAMETERS 30 pos_iter integer maximal number of iterations that are performed by the program to position and orient objects in the simulation cell at the begining of the BD run pos_min integer the minimal allowed surface to surface distance between objects in the simulation cell at the begining of the BD run check_overlaps string whether to check if there are overlapping objects in the simulated system whose po sitions need to be corrected yes no overlaps string defines the algorithm that is applied to look for the overlaps in the system either trees for a hierarchical algorithm or direct for a brute force algorithm overlaps_removal string defines the method that is applied to remove overlaps in the system none if the check_overlaps keyword is set to yes the program will repeat all steps that lead to overlaps with different random vectors atom overlaps will be corrected iteratively by applying appropriate forces and torques each pair of overlapping atoms is evaluated separately molecule overlaps will be corrected iteratively by applying forces and torques and evaluating all pairs of overlapping atoms within a given molecules pair simultaneously nb_list string defines a method that is used to evaluate nonbonded in
5. in case of coarse grained models 9 TOOLS 33 To run computations the user should simply execute the dmatrix command from within the directory where the input txt file is located As a result three text files should be created output txt details of computations and debugging information tensor txt 6 x 6 diffusion tensor matrix transformations txt transformations translations and rotations that should be performed to orient the molecule along main axes of the rotational diffusion tensor with the origin located at the diffusion centre This file can be fed to the orient_mstr utility NOTE when computations are performed for a coarse grained model overlaps between beads de fined in the input txt file are possible only if all beads in a model are of the same size 9 2 orient mstr Tool for orienting molecule in the molecular frame with axes coinciding with the principal axes of diffusional rotations Usage orient mstr transformations structure where transformations file resulting from the execution of the dmatrix tool structure file containing coordinates of a molecule to be transformed PDB PQR MSTR files are supported Results of the orient_mstr execution transformed coordinates of a molecule formatted according to the for mat of the input structure file are printed to the standard output 9 TOOLS 34 9 3 parsemolb Tool for converting binary molb files to text Usage parsemolb file inde
6. present in the simulated system contain elements of the diffusion tensor matrix in the plain text format Elements of D are given in dt files in the following order tt tt tt tr tr tr Df Diz Dis Dy Dia Dis tt tt tt tr tr tr Dz D32 D33 Di D3 Do tt tt tt Er tr tr Dzi Dz D33 D31 Dz D33 rt RE nt rr Gr rr Dii Di Di Ti 12 13 rt rt ne rr rr GR Dai D35 D33 DA Da D33 rt rt rt rr rr rr Dzi D53 D33 D31 D33 33 where units of Di are ps units of por are A ps and units of Dry are ps Elements of D should be given in the molecule fixed frame with origin located at the diffusion center of a molecule and axes coinciding with principal axes of rotations We provide a tool dmatrix see below that can be used to compute the diffusion tensor matrix of a molecule either based on the appropriate PDB file or a plain text file containing 7 RUNNING SIMULATIONS WITH BD_RIGID 23 its coarse grained representation 71 3 External Coordinates File Initial coordinates of objects in the studied system in the laboratory frame can be specified by the user in an external text file which should adhere to the following format O yn vy W a af a where i runs over N molecules ordering of molecules results from parameters defined in the prm file see below the vector u Uz Uy Uz describes translations that are needed to translate the diffusion center of a molecule to its initial position in the la
7. BD BOX version 2 2 user s manual bd rigid Pawet Zielinski Maciej Dtugosz Copyright C 2010 2011 2012 2013 2014 University of Warsaw Pawe Zieli ski Maciej D ugosz This program is free software you can redistribute it and or modify it under the terms of the GNU General Public License as published by the Free Software Foundation either version 3 of the License or at your option any later version This program is distributed in the hope that it will be useful but WITHOUT ANY WARRANTY without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE See the GNU General Public License for more details You should have received a copy of the GNU General Public License along with this program If not see http www gnu org licenses Contact authors Pawe Zieli ski pzielinsOcent uw edu pl Maciej D ugosz mdlugoszOcent uw edu pl Work on the BD_BOX package is supported by National Science Centre grant number N N519 646640 BD_BOX is a scalable Brownian dynamics package for UNIX LINUX platforms BD_BOX is written in C and uses modern computer architectures and technologies e MPI technology for distributed memory architectures e OpenMP technology for shared memory systems e NVIDIA CUDA framework for GPGPU e SSE vectorization for CPU Within the BD_BOX framework macromolecules can be represented either with flexible bead models the bd flex module of the BD_BOX package or as rigid
8. ETERS 26 restart string this keyword indicates that a continuation run should be performed using the string restart file default NULL string input restarts should be specified using command line parameters object string integer the name of the object molecule to be simulated and the number of objects of this kind in the system the program will look for appropriate string mstr and string dt files multiple object lines can be present in the parameter file their order determines the ordering of molecules in all output files as well as the order of objects in the external coordinate file defined with the ext coor file key benchmark string whether to benchmark evaluation of forces yes no default no 8 2 Nonbonded Interactions elec string whether to evaluate electrostatic screened Coulomb interactions yes no default yes elec_desolv string whether to evaluate desolvation cavity forces yes no default yes hdb_surf string whether to evaluate non polar hydrophobic interactions none primitive facts default none elec_cutoff float cutoff for electrostatic interactions also cavity terms default 0 hdb_cutoff float cutoff for nonpolar interactions evaluated with the FACTS model R in Equation 8 de fault 10 0 hdb_a float inner cutoff for nonpolar interactions evaluated with the primitive model a in Equation 10 default 3 1 hdb_b float outer cutoff for nonpolar interactions evaluated wi
9. Q tags is used in evaluation of electrostatic interactions Coulombic and desolvation terms equations 14 and 15 Molecular representation specified with LJ tags is used in evaluation of collisions excluded volume and hy drophobic interactions Lines that begin with the LJ tag should be formatted as follows LJ Xi Y Zi radius Ei dl d2 cO el C21 c31 Hi SASA LJ XM YM ZM radiusy EM dim d2M COM clm C2M C3M UM SASAM 7 RUNNING SIMULATIONS WITH BD RIGID 22 where x y z denote Cartesian coordinates of spheres in the molecular frame radius denotes radius used in evaluations of excluded volume and hydrophobic interactions e is the well depth of excluded volume interactions parameters d and c are for the FACTS model equations 6 7 SASA denotes the solvent accessible surface area of a given sphere in the model Equation 9 and y Equations 5 13 denotes the surface tension parameter Molecular representations given with tags Q and LJ can differ both in the number of spheres atoms and their positions We supply a tool orient mstr see below that can be used to transform rotate and translate molecular coor dinates given either as a mstr pdb or pqr files to the molecular frame that coincides with the main axes of diffusional rotations with the origin located at the diffusion center of a molecule 71 2 Tensor File Diffusion tensor files dt a separate file is needed for each type of a molecule that is
10. aluated in kcal mol e units of eps float dielectric constant of the solvent kappa float the inverse of the Debye length default 0 1 A dxfile string DX formatted file with the electrostatic potential of the molecule obtained by solving the linear Poisson Boltzmann equation input outsourcedx string DX formatted file with the input electrostic potential scaled by ext_phi_scale see below that can be used for debuging and visualization purposes output ext_phi_scale float factor for scaling multiplication values of the input electrostatic potential can be used to A convert units default value 0 5922 scaling from kT gt ggk to skip_points integer every skip_points point of the input grid in each direction x y z will be taken into account upon evaluation of electrostatic potential default 1 all points are taken into account p rfile string input PQR formatted file containing the definition of the evaluated molecule outdx string DX formatted file with the electrostatic potential generated by the derived effective charges output 9 TOOLS 37 outpgr PQR formatted output file with effective charges tol float goodness of electrostatic charges fit this keyword is used to define the stop condition for the program value of 1 0 means perfect fit sensible values of effective charges are typically obtained for tol in the range of 0 85 to 1 0 default value 0 9995 penalty float this variable
11. andom displacement vector in the BD integration scheme Equation 2 may lead to nonphysical overlaps between molecules Overlaps occur when the distance between two atoms that belong to different molecules is smaller than the sum of excluded volume radii of these atoms bd_rigid can check whether a particular simulation step leads to an overlap Next positions of molecules are corrected to remove all cases of overlaps This can be done in one of the four ways One option is that the current step of the simulation is rejected original positions of all molecules are restored and a new step is attempted with a different random vector until there are no cases of overlaps The second option is similar to the first however each of the subsequent attempts is made with a smaller time step After successful removal of all overlaps the time step is restored to its initial value As a third option we also implemented an iterative procedure that removes overlaps between atoms via translations and rotations around diffusion centers of colliding molecules A Z resulting from pseudo forces acting along the line connecting centers of overlapping atoms Meer Ager At EDM rr 20 where the scaling factor D the diffusion tensor Each pair of overlapping atoms is treated separately Translations and rotations of molecules are performed until there are no overlaps in the system The fourth option is similar in the spirit to the third one but o
12. ary conditions for computer simulations Physica A 143 441 467 1987 C W J Beenakker Ewald sum of the rotne prager tensor J Chem Phys 85 1581 1582 1986 D L Ermak and J A McCammon Brownian dynamics with hydrodynamic interactions J Chem Phys 69 1352 1360 1978 A Iniesta and J G de la Torre A second order algorithm for the simulation of the brownian dynamics of macromolecular models J Chem Phys 92 2015 2019 1990 U Winter and T Geyer Coarse grained simulations of a small peptide Effects of finite damping and hydrodynamic interactions J Chem Phys 131 104102 104107 2009 T Geyer and U Winter An o n2 approximation for hydrodynamic interactions in brownian dynamics simulations J Chem Phys 130 114905 114913 2009 REFERENCES 40 14 15 16 17 18 19 20 21 22 23 24 25 T Ando and J Skolnick Krylov subspace mthods for computing hydrodynamic interactions in Brownian dynamics simulations J Chem Phys 137 064106 2012 M X Fernandes and J G de la Torre Brownian dynamics simulations of rigid particles of arbitrary shape in external fields Biophys J 83 3039 3048 2002 P Strating Brownian dynamics simulations of hard sphere suspension Phys Rev E 59 2157 2187 1999 J Sun and H Weinstein Toward realistic modeling of dynamic processes in cell signaling Quantification of macromolecular crowding effects J Chem Phys 127 155105 155115 2007 S R McGuff
13. boratory frame and the vector O Qz Qy Qz describes rotations that are needed to orient the molecule in the laboratory frame N 4 02 05 02 The rotation matrix is defined as 23 My Mi Mi3 M Ma Ma M M3 M32 M33 8 CONTROL FILE KEYWORDS AND COMMAND LINE PARAMETERS 24 where ae 0 Q cosQ 92 11 Q2 0 0 Oz Mia oz U cosQ qna 20 Qy Mi3 02 1 cosQ y no QQ 0 Ma a 1 cosQ y sind 02 02 cosQ 0 M s Q2 QQ Qe Mo3 ee d cos q ma QrQ Qy M31 0 1 cosQ q sind 00 Qe M32 ar U cosQ qna 02 0 cosQ 0 M33 02 7 4 Simulation Control File Control file is based on a number of keywords that allow the user to specify simulations conditions Only these lines in the control file that begin with a particular keyword are recognized by the program Each line of the prm file should contain only one keyword and its value keyword value The complete list of available keywords is given in the next section 8 Control File Keywords and Command Line Parameters Below a complete list of recognized keywords and their possible values is given Values of all keywords relevant to a simulation must be specified by the user default values that are given below are these at which program variables are initialized and they are not always meaningful 8 CONTROL FILE KEYWORDS AND COMMAND LINE PARAMETERS 25 8 1 Input Outpu
14. count the fact that atoms in the molecule overlap the dielectric constant of the solvent p the dielectric constant of a molecule rij the distance between the i charge of one molecule and j atom of another molecule kK the inverse of the Debye length 5 EXTERNAL FIELDS 17 5 External Fields 5 1 Bounding Sphere A spherical surface impenetrable for molecules that encloses simulated system can be included in bd_rigid simulations The influence of such a containing sphere on the studied system is modeled by applying a central force of form gt A F r 16 Rophere je where Rsphere the radius of a spherical surface n a non negative integer A the amplitude of the force r the distance between a diffusion center of a particular molecule and the system s center of geometry The above force is applied to each of the studied molecules outside a predefined cutoff 5 2 Electric Field It is possible to simulate with bd rigid the interactions of charged molecules with external electric fields DC direct current field Ep 17 AC alternate current field E E cos wt g 18 RF rotate field E E cos wt sin wt 0 g 19 where Ep r or Ep y or Ep z Eo is the magnitude of the external field and w its frequency Except for the rotate field all other fields act along the space fixed direction 6 OVERLAPS 18 6 Overlaps The presence of the Gaussian r
15. ee and A H Elcock Atomically detailed simulations of concentrated protein solutions the effects of salt ph point mutations and protein concentration in simulations of 1000 molecule systems J Am Chem Soc 128 12098 2006 M Matsumoto and T Nishimura Mersenne twister A 623 dimensionally equidistributed uniform pseudo random number generator ACM Transactions on Modeling and Computer Simulation 8 330 1998 G E P Box and M E Muller A note on the generation of random normal deviates Annals Math Stat 29 610 611 1958 W Press B Flannery S Teukolsky and W Vetterling Numerical recipes in C the art of scientific computing Cambridge University Press 1992 B Carrasco and J G de la Torre Hydrodynamic properties of rigid particles comparison of different modeling and computational procedures Biophys J 76 3044 3057 1999 D A Beard and Schlick T Unbiased rotational moves for rigid body dynamics Biophys J 85 2973 2976 2003 U Haberthiir and A C Caflisch Facts Fast analytical continuum treatment of solvation J Comput Chem 29 701 715 2008 R R Gabdoulline and R C Wade On the contributions of diffusion and thermal activation to electron transfer between phormidium laminosum plastocyanin and cytochrome f Brownian dynamics simulations with explicit modeling of nonpolar desolvation interactions and electron transfer events J Am Chem Soc 131 9230 9238 2009 REFERENCES 41 26 B Roux and T Simonson I
16. etween atoms i and j ij is a unit vector pointing from the it atom to the neighboring j atom and Oj is a sigmoidal re 1 7 Frig lt Ra ij 8 0 if Tij gt Ri weighting function where R is a predefined radius of a sphere i e cutoff enclosing the i atom Parameters Co C1 Co C3 dy and da that are derived by fitting to exact values of SASA are defined in the PARAM22 CHARMM parameter set 27 For all atoms values of Ris are similar following the original parametrization of the FACTS model 24 4 2 2 Primitive Approach Each atom of a given molecule is assigned a constant SASA value computed in the absence of other molecules SASA The change in the SASA of a particular surface atom i of one molecule due to the proximity of surface atoms j of other molecules SASA is computed as SASA SASA X fijl rig 9 j where fis rij bara a for a lt rij oj lt b 10 Ai an b fis Tis 1 for Tij Oj lt a 11 4 INTERMOLECULAR INTERACTIONS 15 igra U for Tig 0j gt b 12 and amp is a user defined parameter and o is the van der Waals radius of the i atom Parameters a and b are also user defined see Figure 2 Energy is computed as 11000 10750 2 A 5 10500 SASA 10250 1000 10000 10250 10500 10750 11000 gagari Figure 2 Performance of the primitive approach Total values of solvent accessible surfaces calculated for 1500 differen
17. ga to sphere results in a volume of spherical layer between surfaces of spheres of radii omega_low and omega_high respectively right 10 EXAMPLES DISTRIBUTED WITH BD BOX 38 10 Examples distributed with BD BOX Below we present a list of examples that are distributed along with the BD BOX source code Examples are located in examples rigid name directories HEWL A single hen egg white lysozyme molecule All atom model SPHERES Periodic box with hard spheres at volume fraction of 30 QSPHERES Periodic box with oppositely charged spheres HEWL_pH3_100mM 420 molecules of hen egg white lysozyme at volume fraction of 10 Charges assigned at pH 3 and 100mM ionic strength 11 Final Notes We did our best to ensure that the BD_BOX code is bug free We have tested implemented features using rather simple and thus predictable models such as single spheres dumbbells chains but also more elaborate models of proteins Properties of these models derived from BD simulations such as for example translational and rotational diffusion coefficients at different conditions of temperature and viscosity chains end to end distances and their radii of gyration were validated using theoretical analytical predictions and available literature data Additionally we have also cross examined BD_BOX using various hardware platforms We invite users to send comments and questions regarding their own applications of BD BOX Reports on possible b
18. ies can be described either using a coarse grained representation or with fully atomistic details In the latter case intermolecular interactions may include electrostatic hydrophobic and Lennard Jones potentials External electric fields can also be applied to the simulated system Hydrodynamic interactions between molecules modeled within the rigid body framework are currently not supported BD_BOX simulations can be performed without or with boundaries in the latter case the containing or periodic boundary conditions can be used With BD_BOX one can effectively simulate either single molecules or multimolecular systems composed of large numbers of different species For efficient simulations of dense systems we implemented algorithms preventing the overlapping of diffusing molecules T6 17 18 This manual describes the usage of the bd rigid module of the BD BOX package The bd flex module for BD simulations of flexible molecular models is described separately CONTENTS Contents 2 Installation Rigid body Brownian dynamics 3 1 Molecular Models aaa pi a NN ne of ate a ra dO A 33 lEqliatlo s of MO tIO S oe e a a be Pate ee ee wod eS 3 4 Rotations so save kb aaa whew ea bre RR WG we e OR RER GR KOS GG 4 Intermolecular Interactions 4 1 Excluded Volume Interactions 2 2 er rer rr rr rr A 4 2 Hydrophobic InteractionSl LL LL omer rr er rr rr rr rr A 4 2 1 FACTS Approachh ect gu sg dog et e a dl a a i fa
19. inaries After downloading and unpacking the compressed BD_BOX archive gzip d bd_box ver tar gz tar xvf bd_box ver tar the user should execute the configure script from within the BD BOX directory 2 INSTALLATION 1 cd bd box ver configure configure takes instructions from Makefile in and builds Makefile and some other files that work on the user s system We provide the INSTALL file and some other files consult the README file in the bd_box ver direc tory with flags and options that can be passed to the configure script to tune the compilation parameters on a particular system Next the make command followed optionally by the make install should be executed As a result of the compilation process two binaries bd_flex and bd_rigid should be created either within the bd_box ver src directory or in the directory passed by the user to the configure script with the prefix option It is possible to prevent compilation of either one of binaries using disable flex or disable rigid flags Tools distributed with the BD_BOX package should be compiled separately Their sources along with appro priate Makefiles ale located in bd_box ver tools tool_name directories Each tool can be compiled using the make all command 2 INSTALLATION Units charge e length time ps temperature K kcal mol energy 3 RIGID BODY BROWNIAN DYNAMICS 9 3 Rigid body Brownian dynamics 3 1 Molecular
20. is used to prevent the effective charges to have unphysical extremely large magni tudes typically values below 0 01 should used default 0 001 step_attempts integer number of attempts made by the prgram to reduce the number of effective charges default 1000 initially effective charges are positioned at heavy atoms if possible this number is subsequently reduced while preserving the accuracy specified with the tol parameter seed integer seed for the random number generator 31415 omega string defines the region where the fitting is performed possible values are full sphere skin omega low float lower boundary of the fitting region relevant for all values of the omega keyword omega high float upper boundary of the fitting region relevant only for omega values skin and sphere out_filename string plain text log file default err log omega full omega sphere omega skin omega_high omega_low omega_low omega low Omega high Figure 4 Possible definitions of the region outside the molecule where its electrostatic potential is fitted upon evaluation of effective charges When the omega parameter is assigned value skin electrostatic potential is fitted in the volume bounded by surfaces obtained by inflating van der Waals surface of the molecule by omega_low and omega_high left When omega is set to full the upper bondary of the volume results from the definition of the input electrostatic grid center Setting ome
21. ls spheres atoms with central point charges right 3 RIGID BODY BROWNIAN DYNAMICS 10 3 2 Rigid body Diffusion Tensor Matrix In dilute solutions diffusional properties of molecules are determined by their size and shape as well as by the temperature and viscosity of the solvent The information required to characterize translational and rotational movements of quasi rigid Brownian particles in dilute solutions is contained in their diffusion tensors Single particle diffusion tensor D is represented by a symmetric 6 x 6 matrix that consists of four 3 x 3 blocks These blocks are related to translational tt and rotational diffusivities rr and their couplings rt tr 15 D pt D 1 Dt pr Diffusion tensors of molecules can be obtained theoretically from rigid body hydrodynamic calculations performed for detailed molecular models 22 For a spherical particle the above matrix takes a diagonal form isotropic dif fusion in which identical elements Di correspond to the translational diffusion coefficient of a sphere Similarly identical elements Di give its rotational diffusion coefficient For an arbitrarily shaped object its average trans lational and rotational diffusion coefficients are defined respectively as Dt Tr DV and D Tr D where the diffusion tensor Equation 1 is evaluated in the molecule fixed coordinate system with the origin lo cated at the particle s diffusion center and whose axes coincide wi
22. mplicit solvent models Biophys Chem 78 1 20 1999 27 A D Mackerell Jr D Bashford M Bellott R L Dunbrack Jr J D Evanseck M J Field S Fischer J Gao H Guo S Ha D Joseph McCarthy L Kuchnir K Kuczera F T K Lau C Mattos S Michnick T Ngo DT Nguyen B Prodhom W E Reiher III B Roux M Schlenkrich J C Smith R Stote J Straub M Watanabe J Wi rkiewicz Kuczera D Yin and M Karplus All atom empirical potential for molecular modeling and dynamics studies of proteins J Phys Chem B 102 3586 3616 1998 28 N A Baker D Sept S Joseph M J Holst and J A McCammon Electrostatics of nanosystems application to microtubules and the ribosome Proc Natl Acad Sci USA 98 10037 10041 2001 29 RR Gabdoulline and R C Wade Effective charges for macromolecules in solvent J Phys Chem 100 3868 3878 1996 30 X Li Y Levin and M E Fisher Cavity forces and criticality in electrolytes Europhys Lett 26 683 688 1994 31 M D ugosz and Antosiewicz J Anisotropic diffusion effects on the barnase barstar encounter kinetics J Chem Theory Comput 9 1667 1677 2013
23. nr registers different contributions to the total energy of the system during a simulation The binary restart file rst contains all the information needed to continue an interrupted BD run All files are periodically updated during a BD run with frequencies defined by the user Having prepared all the necessary files the user may run the simulation with the command bd_rigid file prm It is also possible not to use the control file at all or to override all options specified in the control file by using command line parameters that correspond to appropriate control file keywords bd rigid parameter value or bd_rigid file prm parameter value 7 RUNNING SIMULATIONS WITH BD_RIGID 21 7 1 Structure File Structure file contains representations of a molecule that are used in evaluation of intermolecular interactions In general the structure file adheres to the following format TAG list of properties Currently two tags are recognized Q and LJ Lines that begin with the tag Q define dielectric electrostatic model of a molecule Molecule is composed of dielectric spheres with central discrete charges Lines that begin with the tag Q should be formatted as follows Q x y z charge radius Q XN YN ZN Chargen radiusn where x y z denote Cartesian coordinates of a dielectric sphere in the molecular frame charge denotes its central charge and radius denotes its radius Molecular representation specified with
24. o keal 331 842 mol alpha float scaling factor for Lennard Jones interactions Equation 4 default 4 0 beta float scaling factor for desolvation cavity forces Equation 14 and 15 default 1 0 8 CONTROL FILE KEYWORDS AND COMMAND LINE PARAMETERS 28 8 3 Boundaries bc string specifies boundary conditions possible values are none an infinite system pbc periodic boundary conditions sphere an impenetrable sphere xbox float the primary simulation cell size in the x direction default 0 0 ybox float the primary simulation cell size in the y direction default 0 0 zbox float the primary simulation cell size in the z direction default 0 0 sphere_radius float the radius of an impenetrable sphere enclosing the simulated system default 0 0 sboundary string whether to use the bounding sphere potential default no another possibility yes Note that this keyword has nothing to do with the bc keyword if bc sphere is used molecules are not allow to cross a spherical surface around the studied system but no potential is used for that rather BD moves leading outside the spherical surface are simply rejected an repeated with different random vectors sboundary_A float the magnitude of the bounding sphere force Equation 16 default 0 0 sboundary_n integer the power of the radial distance dependence of the bounding sphere force Equation 16 default 0 sboundary_cutoff float the bounding sphe
25. of the particle can be computed via Cholesky decomposition of its diffusion tensor TO Brownian dynamics trajectory of an object in the laboratory frame is obtained by applying at each BD simulation step properly defined transformations translations and rotations between the molecule fixed and the laboratory coordinate frame 3 4 Rotations The following matrix is utilized to perform rotational motions 23 My Ming Mis M Ma Ma Maz M3 M32 M33 where 4 INTERMOLECULAR INTERACTIONS 12 0 Q cosQ 92 My 02 QrQ Oz Mia oz U cosQ qna 20 Qy Mi3 02 1 cosQ y no QrQ 0 Ma a 1 cosQ y sind 02 02 cosQ os M Q2 Qh Qr Mo3 d cos d q ma QrQ Qy M31 0 1 cosQ q sind 00 Qe M32 ar U cosQ qna 02 0 cosQ 0 M33 02 where the vector A Qs My 0 describes rotations around axes x y z of the molecule fixed coordinate system and 0 02 07 02 4 Intermolecular Interactions Forces and torques acting on molecules M in Equation are evaluated analytically based on derivatives of different potential functions that are described below During simulations intermolecular forces are evaluated within atom based spherical cutoffs defined by the user Minimum image convention is utilized in case of periodic systems 4 INTERMOLECULAR INTERACTIONS 13 4 1 Excluded Volume Interactions Excluded volume inte
26. ool This small utility program can be used to compute atomic SASA s using atomic coordinates and radii of atoms provided by the user in the PQR formatted input file Usage sasa tool input output N radius where input input pqr file output output file plain text N number of points used to describe the surface of an atom radius radius of the solvent probe 9 TOOLS 36 9 6 charge optimizer Utility program that can be used to derive a set of effective charges Effective charges positioned initially at heavy atoms of the molecule are derived by fitting the electrostatic potential resulting from the Debye Hiickel approxi mation to the external molecular potential obtained as a numerical solution of the Poisson Boltzmann equation see for details The total number of effective charges is optimized upon the execution of charge_optimizer Fitting of effective charges can be performed in various regions around the molecule the whole space around the molecule spherical layer centered on the molecule or the molecular skin defined by inflating van der Waals radii of molecule s atoms see Figure 4 Execution of the program is controlled via keywords given in the plain text parameter file Usage charge optimizer file prm and the following keywords can be specified in the parameter file in the keyword value manner gamma float unit conversion factor the default value of 331 342 results in electrostatic potentials ev
27. r A 4 2 2 Primitive Approach ooa 4 3 Electrostatic Interactions LL rr e 4 3 1 Screened Coulombic Interactions s smo a a a a a 4 3 2 Polarization terms 2 5 External Fields 5 1 Bounding Sphere ee b 2 Electric Field ac 4 0 2244554546544 de 462 od s r oh oh SRSA EGE a E SR dow 4 6 Overlaps imulations with bd_rigid 7 1 tructure File ss ss ss n hs L2 Tensor Eile 4 3 0 4201 5 dd rn der dh bobr w So a a i Go atak EO Bbw ok Trak one 7 3 External Coordinates File soo ses rr sr rr rr rr rr rr sr vvs vs nr 1 4 imulation Control Flel s ns CONTENTS 8 Control File Keywords and Command Line Parameters 8 1 Input Output Control 8 2 Nonbonded Interactions 8 3 Boundaries 8 4 External Electric Fields 8 5 Algorithms 8 6 Physical Conditions rr rer rr rr rr rr rr rr rr rr rr ere ser 8 7 Devices Control 9 Tools 9 2 orient mstr 9 3 parsemolb 9 4 molb2dcd 9 5 _sasa tool 9 6 charge optimizer 0 Examples distributed with BD BOX Final Notes 24 25 26 28 28 29 31 31 31 32 33 34 34 35 36 38 38 1 REQUIREMENTS 6 1 Requirements The BD_BOX package is distributed as source code A UNIX LINUX make tool is needed to build a working binary from the source code see below CPU versions of BD_BOX binaries can be run either in a serial mode or in parallel either on shared memory machines using OpenMP and MPI libraries or on architectures
28. ractions are modeled using pairwise terms of form E or fet p SE Pa 4 i rij rij where n and m non negative integers i and j indexes that run over atoms from different molecules o denotes an atomic radius elI the interaction well depth assigned to a given atom a and B user defined constants rij the distance between the two interacting atoms values of n m can be set to 12 6 12 0 8 6 8 0 4 2 Hydrophobic Interactions Hydrophobic i e nonpolar interactions can be evaluated using one of the two solvent accessible surface based approaches implemented in BD_BOX either the FACTS model or the model that we will further refer to using the therm primitive similar in the spirit to the approach described by Gabdoulline and Wade 25 4 2 1 FACTS Approach In the FACTS model the nonpolar contribution to the free energy of solvation of a molecule is modeled in relation to the amount of its surface exposed to the solvent SASA 24 N pg ronpolar DD SASA 5 i l 4 INTERMOLECULAR INTERACTIONS 14 where summation runs over all NV atoms of a molecule or molecules y denotes the empirical surface tension parameter 26 and atomic SASAs are given with 24 C1 SASA co apa A N Tun O j tj 3 1 Og yj Di V yj da a a ve sh sy N Vj N Vj Gi A F Hi jaa rj Did In the above equations V denotes the van der Waals volume of the je atom rij is the distance b
29. re force is applied outside this cutoff radius Equation 16 default 0 0 8 4 External Electric Fields E_ext string whether to switch on off the external electric field yes no default no 8 CONTROL FILE KEYWORDS AND COMMAND LINE PARAMETERS 29 E_magn float the magnitude of the external electric field default 0 0 E type string a choice between different types of the electric field possible values are AC DC or RF de fault DC E freq float the frequency where applicable of the external electric field units FH default 0 0 E_dirl string the direction of the external electric field x y or z default x E_dir2 string the second direction of the external electric field RF x y or z default y E_factor float the unit conversion factor for the external electric field to keal el 4 1 default 1 0 8 5 Algorithms bd_algorithm string this keyword defines the algorithm that is used to generate the BD trajectory of the studied system possible values are ermak_const the time step is constant during the simulation or ermak_var variable time step is used in order to prevent overlaps between simulated objects bd_steps integer the number of BD simulation steps dt float time step min_dt float minimal possible value of the time step when the variable time step algorithm is used rand_seed integer seed for the random number generator move_attempts integer the maximal number of attempts that are made
30. t Control out filename string text file containing diagnostic informations regarding the current bd_rigid run output dcd_filename string DCD formatted binary trajectory file default NULL string output save_dcd_freq integer frequency number of steps for writing to the DCD trajectory default 1 pqr_filename string PQR formatted file containing the startup configuration of the simulated system de fault NULL string output rep string this keyword is used to define which molecular representation from the mstr file tags Q and LJ should be written to the PQR formatted file possible options are j Lennard Jones spheres defined with the LJ tags or elec charges defined with the Q tags default j molb_filename string binary trajectory file mo b default NULL string output save_molb_freq integer frequency number of steps for writing to the molb trajectory default 1 enr filename string plain text energy file default NULL string output save enr_freq integer frequency number of steps for writing to the energy file default 1 rst_filename string binary restart file default NULL string output save_rst_freq integer frequency number of steps for writing to the binary restart file default 1 ext_coor_file string an optional plain text file describing an user defined configuration of the simulated sys tem default NULL string input 8 CONTROL FILE KEYWORDS AND COMMAND LINE PARAM
31. t relative configurations of two small proteins with amp 0 5 a 1 3 and b 6 4 All atom protein models were used in calculations with atomic radii taken from the PARAM22 CHARMM parameter set 27 Exact SASA values were calculated using APBS 28 N pgronpolar LL 5 SASA SASA 13 1 1 4 INTERMOLECULAR INTERACTIONS 16 4 3 Electrostatic Interactions 4 3 1 Screened Coulombic Interactions The energy of the electrostatic interaction screened by dissolved ions between two charges q and q that belong to different molecules is computed as e Brij D y EG 14 Tij where kK the inverse of the Debye length rij the distance between the two charges y the unit conversion factor an user defined parameter the dielectric constant of the solvent 4 3 2 Polarization terms Polarization desolvation terms are evaluated between the charges of one molecule and the dielectric cavity of another molecule Dielectric body of a molecule is treated as a conglomerate of spheres i e atoms 29 Energy of the interaction between the i charge qi of one molecule with the j dielectric sphere radius o of another molecule is computed as 80 Ey 1 krij e 2e ep r Evity AB ij exp erya qa 15 where y the unit conversion factor an input parameter B the user defined an input parameter constant that depends on ionic strength it also takes into ac
32. teractions brute brute force method spatial a method similar in the spirit to the cell linked lists algorithm however without explicit creation of nonbonded lists verlet Verlet lists algorithm with a twin range cutoff verlet_roff float defines the value of the outer cutoff in the twin range cutoff Verlet algorithm i e the outer cutoff is computed as a cutoff increased by verlet_roff verlet_count integer Verlet lists will be refreshed during a simulation every verlet_count steps dm_alg float defines the method used to perform thread decomposition in case of parallel runs either atomdec 9 TOOLS 31 decomposition is based on the indexes of atoms or spatialdec decomposition is performed based on Cartesian coordinates of atoms 8 6 Physical Conditions T float temperature default 298 15 kappa float the inverse of the Debye length default 0 1 eps float solvent dielectric constant default 78 54 epsin float solute dielectric constant default 1 8 7 Devices Control cuda_devices integer integer CUDA devices to be used their numbers followed by system identifiers cuda_blockinteger dimension of the thread block relevant in case of CUDA devices with the compute ca pability smaller than 2 0 9 Tools We provide a set of tools that can be used to prepare input files necessary to perform BD simulations with the BD_BOX bd rigid module These are described below 9 TOOLS 32 9 1 dmatri
33. th the primitive model a in Equation 10 default 4 35 hdb_phi float scaling factor for nopolar interactions evaluated with the primitive approach default 0 5 8 CONTROL FILE KEYWORDS AND COMMAND LINE PARAMETERS 21 lj cutoff float cutoff for Lennard Jones interactions default 10 0 setting this parameter to 1 results in cutoffs defined for each pair of atoms as 26 Ri R useful in case of the 12 6 form of the L J potential j_repul float repulsive term in the Lennard Jones potential Equation 4 r J eP possible values 8 12 default 12 pij samo lj attract float attractive term in the Lennard Jones potential Equation possible values 0 6 default 6 repul_B float the magnitude of the repulsive part in the Lennard Jones potential Equation 4 default 1 probe_radius float radius of the spherical probe used to determine surface atoms during simulations colli sions are evaluated between surface atoms default 1 4 srf_ratio float value used to determine whether an atom is buried or belong to the surface of a molecule total i e if for a given atom the ratio sur face sur face is smaller than srf ratio the atom belongs to the surface of the molecule default 0 9 points_per_sphere integer number of points used to describe the surface of an atom default 1000 gamma float unit conversion factor for electrostatic interactions Equation and 15 default conversion t
34. th the principal axes of the rotational diffusion tensor 3 3 Equations of Motions The propagation scheme for an arbitrarily shaped rigid body that is described with a fully anisotropic diffusion tensor can be written in the molecule fixed frame again by molecule fixed frame we mean that the origin of the frame is assumed to be located at the diffusion center of the molecule and its axes coincide with the principal axes of diffusional rotations as 10 15 At gt zj Az At DM R At 2 kgT where D the precomputed and constant in the simulation 6 x 6 diffusion tensor of the molecule evaluated in the molecule fixed coordinate system with the origin located at the molecule s diffusion center and axes coinciding 3 RIGID BODY BROWNIAN DYNAMICS 11 with the principal axes of particle s rotations so that the submatrix D in Equation 1 is diagonal At the time step gt Z the vector describing the position of the diffusion center 7 and orientation 4 of the molecule z 7 4 M the generalized force vector resulting for example from interactions with other molecules or external fields and having two components the total force F and the total torque T referred to the molecule s diffu sion center M F T T R At random displacement vector arising from the Brownian noise with zero mean and the variance covariance given with lt R At R At gt 2DAt 3 Random displacements
35. ugs are also welcomed References 1 K Kremer and G S Grest Dynamics of entangled linear polymer melts A molecular dynamics simulation J Chem Phys 92 5057 5087 1990 REFERENCES 39 2 3 4 5 6 7 8 9 10 11 12 13 J G de la Torre J G H Cifre A Ortega R R Schmidt M X Fernandes H E P S nchez and R Pamies Simuflex algorithms and tools for simulation of the conformation and dynamics of flexible molecules and nanoparticles in dilute solution J Chem Theory Comput 5 2606 2618 2009 G M Bell S Levine and McCartney L N Approximate methods of determining the double layer free energy of interaction between two charged colloidal spheres J of Colloid Int Sci 33 335 359 1970 M Aubouy E Trizac and L Bocquet Effective charge versus bare charge an analytical estimate for colloids in the infinite dilution limit J Phys A Math Gen 36 58355840 2003 J Rotne and S Prager Variational treatment of hydrodynamic interaction in polymers J Chem Phys 50 4831 4838 1969 H Yamakawa Transport properties of polymer chains in dilute solution hydrodynamic interaction J Chem Phys 53 436 444 1970 J G de la Torre and V A Bloomfield Hydrodynamic properties of macromolecular complexes i translation Biopolymers 16 1747 1763 1977 E R Smith K Snook and W van Megen Hydrodynamic interactions in brownian dynamics simulations i periodic bound
36. verlaps between two different molecules are evaluated simultaneously and all pairs of overlapping atoms from a given molecule pair are considered in order to determine the directions of correcting forces 7 Running Simulations with bd rigid Various input files are needed to perform BD simulations with the bd_rigid module of the BD_BOX package Figure 3 7 RUNNING SIMULATIONS WITH BD RIGID 19 INPUT molecular models control file PQR formatted structure debug info initial coordinates optional souti ou err restart optional molb dcd restart trajectory energy OUTPUT Figure 3 Schematics of the BD BOX bd rigid module Structure files mstr contain definitions of molecular models that are used to represent different types of molecules to be simulated A separate mstr file is required for each kind of a molecule that is present in the simulated system for example when monodisperse system is to be simulated a single structure file will suffice the number of molecules of each kind in the system is specified by the user in the prm file Diffusion tensor files dt contain definition of diffusion tensors of simulated molecules A separate dt file is required for each kind of a molecule that is present in the simulated system Simulation parameters are stored in the text control file prm Figure B These include for example the type of the propagation algorithm and its time step temperature
37. with distributed memory using the MPI library GPU versions of BD_BOX binaries require the CUDA Toolkit obtainable freely at http developer nvidia com Both CPU and GPU versions can be compiled with the support for either single or double precision floating point arithmetic For the generation of random Gaussian numbers BD_BOX uses functions implemented in the GNU Scientific Library GSL http www gnu org software gsl that are based on the Mersenne Twister algorithm of Mat sumoto and Nishimura 19 We also implemented a function that uses the standard system routine drand48 and the polar Box Muller transformation 20 A generic sequential algorithm implemented in BD_BOX for the Cholesky factorization follows the choldc subroutine described in 21 However it is better to compile BD_BOX with the implemented support for the LAPACK or SCALAPACK library http www netlib org lapack that offers efficient highly memory optimized routines for the Cholesky factorization LAPACK performs the vast amount of operations by exploiting BLAS Machine specific implementations of BLAS are available and such implementations can be crucial for the per formance of BD_BOX The Cholesky factorization on GPU relays on the implemented support for the MAGMA dense linear algebra library http icl cs utk edu magma which provides functionality of LAPACK 2 Installation The configure shell script written by GNU Autoconf is used to build the BD_BOX b
38. x This utility can be used to compute the rigid body diffusion tensor of a molecule In case of an atomistic model of a molecule coordinates are read from a PDB formatted file However if the user needs to compute the diffusion tensor of a coarse grained molecular model its coordinates should be specified directly in the text input file dmatrix requires an input file input txt the name of the file should not be changed that adheres to the following format input txt computations are to be performed for a molecule whose coordinates are stored in a PDB formatted file comment single line containing comments can be left empty 1 coordinates of a molecule are stored in the PDB file name of the appropriate PDB file name of the PDB file with coordinates of a molecule viscosity temperature viscosity of the solvent Poisse and temperature K solvent radius radius of solvent molecules typically set to 1 4 water input txt computations are to be performed for a coarse grained model comment single line containing comments can be left empty 0 computations are to be performed for a coarse grained model number of beads in a model number of spherical subunits atoms in a model x y z radius coordinates of a spherical subunit and its radius xy z radius xy z radius viscosity temperature viscosity of the solvent Poisse and temperature K solvent radius radius of solvent molecules typically set to 0 0A
39. x precision where file name of the molb file to be converted index index of a molecule whose coordinates are to be extracted numbering of molecules starts with 0 precision either float or double Output is printed to the standard output using the following format time Vz Vy Vz O 1 o 9 where time denotes current time vz Vy Vz are components of the translation vector from the molecule fixed to the laboratory frame and o 1 o 9 are components of the rotation matrix Current position of a given molecule in the laboratory frame can be obtained by applying rotations and then translations to the input structure of a molecule 9 4 molb2dcd Tool for converting binary molb trajectory files to binary dcd files Usage molb2dcd input type precision name where input text input file type either all or separate coordinates of molecules are either stored in a single file or separate dcd files are created for each molecule precision either double or float name input molb file 9 TOOLS 35 The text input file is used to specify molecules for which DCD trajectory will be created it contains names of PDB files and indexes of molecules of interest in the source molb file Indexes can be given as a list or as a range An exemplary text input file is given below exemplary text input file for the molb2dcd utility name0 pdb 0 namel pdb 1 10 name2 pdb 11 12 13 14 15 name3 pdb 16 9 5 sasa_t
Download Pdf Manuals
Related Search