Re: [gmx-users] Finding centroid for a bunch of residues
Hi there, This is the question out of gromacs..but I need it urgently.. and I hope this is the only place where I can get such expert to solve my query... while trying to restrict my MDRUN for a particular site of the protein molecule I want to visualize the site and find out the centroid for the particular site.. Is there any visualization tool that can do the same .. I mean Is there any molecular visualization tool that can help in finding out the ...centroid between a group of resuidues ? Pymol has a script that computes the center of mass for a selection, and optionally, moves this center of mass to the origin: http://pymolwiki.org/index.php/COM Haven't used it myself, so I can't vouch for its effectiveness. Steve Kirk -- Dr. Steven R. Kirk [EMAIL PROTECTED], [EMAIL PROTECTED] Dept. of Technology, Mathematics Computer Science (P)+46 520 223215 University West (F)+46 520 223299 Trollhattan 461 86 SWEDEN http://beacon.webhop.org ___ gmx-users mailing listgmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/search before posting! Please don't post (un)subscribe requests to the list. Use the www interface or send it to [EMAIL PROTECTED] Can't post? Read http://www.gromacs.org/mailing_lists/users.php
[gmx-users] Checkpointing GROMACS jobs
Hello, I have been using GROMACS for some very long (in wall clock terms) simulations, and am curious as to how other users on this list solve the problem of checkpointing long MD runs. It's a problem because of the tendency of computational nodes in large HPC facilities (the more processors, the more prevalent the problem, it seems) to keel over near the end of a very time consuming run. Intermittent disk and scheduler faults can also trigger such conditions. Checkpointing at the operating system level is very system-specific, and occasionally compilers can produce executable 'dump' files that continue from where your program left off, but I'm thinking that someone must have automated this process directly using conventionally-compiled GROMACS executables. Of course, it is possible to do an exact continuation from a crashed run using .edr and trajectory (.trr) files by generating a new .tpr from the last trajectory frame that had both position and velocity data. This seems to be, by necessity, an entirely interactive process (unless someone out there has a cool auto-restart script ..). I am thinking more in terms of 'proactive' checkpointing for long jobs, by the following process: A script parses the desired .mdp file describing the user's MD run of T timesteps, then asks the user how many sections (N) to split the run into. The script will then auto-generate a shell script containing all the necessary GROMACS commands to: * Generate a new .mdp file almost identical to the original, but with the number of timesteps set to T/N. * Run N successive mdrun commands, where the output .trr and .edr files from each short run using the modified .mdp file are used, to generate an 'exact restart' .tpr file for the next 'mdrun' command, with the appropriate continuation flag set. * Log (to a file) how many of the N partial runs have been completed, in such a way that if the shell script containing the commands is restarted, it will jump to the correct point in the sequence, restarting from the most recently completed partial run. Has anyone else already solved this problem, or have a method implementing some of the desirable properties above that I can then extend to do exactly the things described above? -- Dr. Steven R. Kirk [EMAIL PROTECTED], [EMAIL PROTECTED] Dept. of Technology, Mathematics Computer Science (P)+46 520 223215 University West (F)+46 520 223299 P.O. Box 957 Trollhattan 461 29 SWEDEN http://beacon.webhop.org ___ gmx-users mailing listgmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/search before posting! Please don't post (un)subscribe requests to the list. Use the www interface or send it to [EMAIL PROTECTED] Can't post? Read http://www.gromacs.org/mailing_lists/users.php
Re: [gmx-users] Re: Targeted MD
Mark Abraham [EMAIL PROTECTED] wrote Subject: Re: [gmx-users] Re: Targeted MD To: Discussion list for GROMACS users gmx-users@gromacs.org Message-ID: [EMAIL PROTECTED] Content-Type: text/plain; charset=ISO-8859-1; format=flowed wei-xin xu wrote: Some hints on practices that generally *not a good idea* to use: * Do not use separate thermostats for different components of your system. Some molecular dynamics thermostats only work well in the thermodynamic limit. If you use one thermostat for, say, a small molecule, another for protein, and another for water, you are likely introducing errors and artifacts that are hard to predict. In particular, do not couple ions in aqueous solvent differently from that * solvent. Sorry that I do not actually understand here. The link I copied above shows that better not to couple ions in aqueous solvent differently from that solvent. Maybe not separately but differently (mean different temperature)? differently is intended to mean in a separate group. I'll reword my wiki sentence. The original poster showed an .mdp file where tc-grps = Protein ; SOL CL UNK (or something like that). Now actually, the semicolon starts a comment, so he's only thermostatting the protein. That's a bad idea because it will lead to net heat flow from the protein to the rest of the system. Even if there were no semicolon, there are probably a few thousand solvent molecules and a handful of chloride ions. Temperature is defined from the average kinetic energy, and the average kinetic energy of a handful of ions in thermal contact with many other atoms will have large fluctuations, and this will lead to the thermostat doing lots of corrections, for lots of heat flow in and out of the system. So treating solvent+ions+other_small_stuff as one group for T-coupling purposes is a good idea, and the standard group Non-Protein serves well here. So a usual tc-grps line has Protein Non-Protein for a protein simulation. Mark A supplementary question. The tc-grps line can take predefined standard group names such as 'System', 'Protein' and 'Non-Protein'. Does the 'Protein' group need to actually BE a protein, or are 'Protein' and 'Non-Protein' really synonyms for 'PresumablyBigMoleculeOfInterestProteinOrNot' and 'EverythingElse' ? Thanks! Steve Kirk -- Dr. Steven R. Kirk [EMAIL PROTECTED], [EMAIL PROTECTED] Dept. of Technology, Mathematics Computer Science (P)+46 520 223215 University West (F)+46 520 223299 P.O. Box 957 Trollhattan 461 29 SWEDEN http://beacon.webhop.org ___ gmx-users mailing listgmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/search before posting! Please don't post (un)subscribe requests to the list. Use the www interface or send it to [EMAIL PROTECTED] Can't post? Read http://www.gromacs.org/mailing_lists/users.php
[gmx-users] Re: Diagnosing an explosion
Apologies for replying to myself - maybe I can sharpen up some of my original questions below: Reply below I am having a problem with my BD simulation either crashing with a Range error or locking up the mdrun process completely. I am running Gromacs 3.3.1 (must use this version ATM because HPC provider has not upgraded yet) on a quad core Intel machine, w/ Fedora Core 7. I have 1000 'particles' in a cubic box (configuration here: http://datavet.hv.se/personal/SK/publicfiles/colloid.pdb where the particles are placed randomly with no two particle centres closer than 4 nm. The particles are defined in the topology here: ; Complete topology file ; ; particles.itp [ defaults ] 1 1 no 1.0 1.0 [ atomtypes ] AR 28699.81134 0.00 A 1.0 1.0 Looks like you are using dimensionless eps and sigma while using otherwise MD units (e.g. T). Check chapter 2 of the manual. You want to multiply your T by boltz (0.008314 or something like that) probably. Aha! Thank you David. My goal was to use MD units in everything, since my tabulated data file was extracted from another GROMACS simulation and is hence already in MD units (kJ/mol vs. nm), and I would prefer to specify temperatures in K. My reason for specifying 1.0 for the last two parameters was that I thought this would allow me to use, unscaled, the tabulated potential (which is already in MD units). If I want to do this, what alternative values should I use for the last two parameters on the [ atomtypes ] line above (given that I am really only using the g(r) and g''(r) columns to hold my potential data: the h(r) and h''(r), which I don't care about, are all zeroes) ? If I change the combination type to 2 or 3 in the [ defaults ] line, will this remove the unit scaling problem ? I have checked Chapter 2 in the manual and I am still unsure how to handle this. I took your advice (which matches with manual Section 2.4 on *reduced* units) and scaled my temperature to 2.494353 = 300K * 0.008314. This did not give a crash, but increasing the timestep and decreasing bd_fric made the problem return. Scaling the temperature value in this way did not fix the 'equilibration around 2x requested temperature' problem mentioned below - g_energy still shows the temperature settling at a numerical value around 4.988. Can anyone suggest a reason for this strange thermostat behaviour ? Many thanks, Steve Kirk -- Dr. Steven R. Kirk [EMAIL PROTECTED], [EMAIL PROTECTED] Dept. of Technology, Mathematics Computer Science (P)+46 520 223215 University West (F)+46 520 223299 P.O. Box 957 Trollhattan 461 29 SWEDEN http://taconet.webhop.org ___ gmx-users mailing listgmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/search before posting! Please don't post (un)subscribe requests to the list. Use the www interface or send it to [EMAIL PROTECTED] Can't post? Read http://www.gromacs.org/mailing_lists/users.php
[gmx-users] Re: diagnosing an explosion
Reply below I am having a problem with my BD simulation either crashing with a Range error or locking up the mdrun process completely. I am running Gromacs 3.3.1 (must use this version ATM because HPC provider has not upgraded yet) on a quad core Intel machine, w/ Fedora Core 7. I have 1000 'particles' in a cubic box (configuration here: http://datavet.hv.se/personal/SK/publicfiles/colloid.pdb where the particles are placed randomly with no two particle centres closer than 4 nm. The particles are defined in the topology here: ; Complete topology file ; ; particles.itp [ defaults ] 1 1 no 1.0 1.0 [ atomtypes ] AR 28699.81134 0.00 A 1.0 1.0 Looks like you are using dimensionless eps and sigma while using otherwise MD units (e.g. T). Check chapter 2 of the manual. You want to multiply your T by boltz (0.008314 or something like that) probably. Aha! Thank you David. My goal was to use MD units in everything, since my tabulated data file was extracted from another GROMACS simulation and is hence already in MD units (kJ/mol vs. nm), and I would prefer to specify temperatures in K. My reason for specifying 1.0 for the last two parameters was that I thought this would allow me to use, unscaled, the tabulated potential (which is already in MD units). If I want to do this, what alternative values should I use for the last two parameters on the [ atomtypes ] line above (given that I am really only using the g(r) and g''(r) columns to hold my potential data: the h(r) and h''(r), which I don't care about, are all zeroes) ? I have checked Chapter 2 in the manual and I am still unsure how to handle this. I took your advice (which matches with manual Section 2.4 on *reduced* units) and scaled my temperature to 2.494353 = 300K * 0.008314. This did not give a crash, but increasing the timestep and decreasing bd_fric made the problem return. Scaling the temperature value in this way did not fix the 'equilibration around 2x requested temperature' problem mentioned below - g_energy still shows the temperature settling at a numerical value around 4.988. Many thanks for any further suggestions! [ bondtypes ] ; molecule.itp [ nonbond_params ] AR AR 1 1.0 1.0 [ moleculetype ] AR 0 [ atoms ] 1 AR 1 AR AR 1 0 28699.81134 ; colloid.top ; #include particles.itp ; #include molecule.itp [ system ] Aggregation simulation [ molecules ] AR 1000 The particles interact via the tabulated non-bonded potential held here: http://datavet.hv.se/personal/SK/publicfiles/table.xvg and the mdrun is controlled by the mdp file settings: ; Input file ; title = Aggregation cpp = /lib/cpp constraints = none integrator = bd dt = 0.002; ps ! nsteps = 100 ; nstcomm = 1 nstxout = 10 nstvout = 10 nstfout = 0 nstlog = 1 nstenergy = 1 nstlist = 5 ns_type = grid coulombtype = User vdwtype = User rlist = 5.036 rcoulomb= 5.036 rvdw= 5.036 table-extension = 0.3 comm-mode = None bd_fric = 5000 ; Berendsen temperature coupling is on in two groups Tcoupl = berendsen tc-grps = System tau_t = 0.1 ref_t = 300 ; Energy monitoring energygrps = AR ; Isotropic pressure coupling is now off Pcoupl = no Pcoupltype = isotropic tau_p = 0.5 compressibility = 4.5e-5 ref_p = 1.0 ; Generate velocites is on at 300 K. gen_vel = yes gen_temp= 300.0 gen_seed= -1 Two strange things happen in this simulation - between the first and second step, the temperature approximately doubles to 600K, and remains approximately at this value for almost as long as the simulation runs. Why does the thermostat stabilize at almost exactly 2x the requested temperature? The second strange event is the explosion (I have looked at the output from g_energy up to this point, and the total and potential energies drop nicely up to the point where the simulation crashes, either with a range error or a frozen process. From md.log: Step Time Lambda 59840 119.680010.0 Energies (kJ/mol) LJ (SR) Coulomb (SR) PotentialKinetic En. Total Energy -2.18045e+040.0e+00 -2.18045e+047.47528e+03 -1.43292e+04 Temperature Pressure (bar) 5.99376e+021.67395e-01 Step Time Lambda 59841 119.682010.0 Energies (kJ/mol) LJ (SR) Coulomb (SR) PotentialKinetic En. Total Energy -2.85992e+290.0e+00 -2.85992e+299.18128e+25 -2.85900e+29 Temperature Pressure (bar) 7.36165e+24
[gmx-users] Diagnosing an explosion
Hello, I am having a problem with my BD simulation either crashing with a Range error or locking up the mdrun process completely. I am running Gromacs 3.3.1 (must use this version ATM because HPC provider has not upgraded yet) on a quad core Intel machine, w/ Fedora Core 7. I have 1000 'particles' in a cubic box (configuration here: http://datavet.hv.se/personal/SK/publicfiles/colloid.pdb where the particles are placed randomly with no two particle centres closer than 4 nm. The particles are defined in the topology here: ; Complete topology file ; ; particles.itp [ defaults ] 1 1 no 1.0 1.0 [ atomtypes ] AR 28699.81134 0.00 A 1.0 1.0 [ bondtypes ] ; molecule.itp [ nonbond_params ] AR AR 1 1.0 1.0 [ moleculetype ] AR 0 [ atoms ] 1 AR 1 AR AR 1 0 28699.81134 ; colloid.top ; #include particles.itp ; #include molecule.itp [ system ] Aggregation simulation [ molecules ] AR 1000 The particles interact via the tabulated non-bonded potential held here: http://datavet.hv.se/personal/SK/publicfiles/table.xvg and the mdrun is controlled by the mdp file settings: ; Input file ; title = Aggregation cpp = /lib/cpp constraints = none integrator = bd dt = 0.002; ps ! nsteps = 100 ; nstcomm = 1 nstxout = 10 nstvout = 10 nstfout = 0 nstlog = 1 nstenergy = 1 nstlist = 5 ns_type = grid coulombtype = User vdwtype = User rlist = 5.036 rcoulomb= 5.036 rvdw= 5.036 table-extension = 0.3 comm-mode = None bd_fric = 5000 ; Berendsen temperature coupling is on in two groups Tcoupl = berendsen tc-grps = System tau_t = 0.1 ref_t = 300 ; Energy monitoring energygrps = AR ; Isotropic pressure coupling is now off Pcoupl = no Pcoupltype = isotropic tau_p = 0.5 compressibility = 4.5e-5 ref_p = 1.0 ; Generate velocites is on at 300 K. gen_vel = yes gen_temp= 300.0 gen_seed= -1 Two strange things happen in this simulation - between the first and second step, the temperature approximately doubles to 600K, and remains approximately at this value for almost as long as the simulation runs. Why does the thermostat stabilize at almost exactly 2x the requested temperature? The second strange event is the explosion (I have looked at the output from g_energy up to this point, and the total and potential energies drop nicely up to the point where the simulation crashes, either with a range error or a frozen process. From md.log: Step Time Lambda 59840 119.680010.0 Energies (kJ/mol) LJ (SR) Coulomb (SR) PotentialKinetic En. Total Energy -2.18045e+040.0e+00 -2.18045e+047.47528e+03 -1.43292e+04 Temperature Pressure (bar) 5.99376e+021.67395e-01 Step Time Lambda 59841 119.682010.0 Energies (kJ/mol) LJ (SR) Coulomb (SR) PotentialKinetic En. Total Energy -2.85992e+290.0e+00 -2.85992e+299.18128e+25 -2.85900e+29 Temperature Pressure (bar) 7.36165e+241.50577e+21 Step Time Lambda 59842 119.684010.0 Energies (kJ/mol) LJ (SR) Coulomb (SR) PotentialKinetic En. Total Energy -inf0.0e+00 -infinfnan Temperature Pressure (bar) infinf Step Time Lambda 59843 119.686000.0 Energies (kJ/mol) LJ (SR) Coulomb (SR) PotentialKinetic En. Total Energy nan0.0e+00nannannan Temperature Pressure (bar) nannan and then comes the range error. Clearly something very nasty happens between step 59840 and 59841. Both reducing the timestep and increasing bd_fric delay the crash, but do not prevent it. I have inspected the trajectory visually using VMD and it does not look unusual. Doing an EM minimization step on the original coordinates before running the full BD run does not cure the problem either. I would be very grateful for suggestions as to what is causing this (my main suspect at the moment is the tabulated potential) and what I can do to either fix it or diagnose the true cause. Many thanks in advance, Steve Kirk -- Dr. Steven R. Kirk [EMAIL PROTECTED], [EMAIL PROTECTED] Dept. of Technology, Mathematics Computer Science (P)+46 520 223215 University West (F)+46 520 223299 P.O.
RE: [gmx-users] Coarse-graining and tabulated non-bonded potentials - will write
Berk Hess kindly answered my original query below. I have some more questions that might be of general interest to the list, so I'm conducting another step of the conversation in the list. - Date: Mon, 26 Nov 2007 14:18:36 +0100 From: Berk Hess [EMAIL PROTECTED] Subject: RE: [gmx-users] Coarse-graining and tabulated non-bonded potentials -will write To: gmx-users@gromacs.org Message-ID: [EMAIL PROTECTED] Content-Type: text/plain; format=flowed From: Steven Kirk [EMAIL PROTECTED] Reply-To: Discussion list for GROMACS users gmx-users@gromacs.org To: gmx-users@gromacs.org Subject: [gmx-users] Coarse-graining and tabulated non-bonded potentials - will write up on the wiki Date: Mon, 26 Nov 2007 13:57:37 +0100 Hello all, Firstly, many thanks to everyone who has contributed useful advice to me over a number of years using GROMACS. I have performed a large number of potential of mean force calculations for the forces acting between two approximately spherical amorphous silica particles (various sizes 5 nm diameter) in TIP4P water with PME electrostatics and varying concentrations of background ions. Now I want to 'coarse-grain' the simulation, treating each silica particle as a single point mass, and use the interaction potential between the particles obtained from the PMF results as a tabulated potential in mdrun, to allow longer time- and size-scales to be investigated (I want to examine colloidal aggregation behaviour of these particles). I have tabulated the PMF potential and its derivatives as a function of centre-of-mass separation of the particles as suggested in the manual (but I can only tabulate for COM-COM distances greater than the contact distance [ = 2r in the hard-sphere approximation] out to some cutoff). Will I have to add a short-distance 'hard-sphere' wall to my tabulated potentials? I don't really understand the issue here. You (would) also directly get the repulsive part of the PMF from a constraint simulation. I assume that you mean that you have not done simulations to obtain the PMF at shorter distances? If not, just do so. Yes, I have calculated the PMF at shorter distances, and it becomes repulsive at shorter distances. My confusion was caused by my uncertainty as to what should be in the tabulated potential file for very short distances (tabulated potentials must be tabulated for distances r = 0, even if I only have real tabulated potential data from some r_min (0) upward). A quick check of the example 'table6-9.xvg' file for the directory referenced in the standard GROMACS distribution using the environment variable GMXLIB has shown me that I can safely 'pad out' all the columns for distances 0 = r r_min with zeroes. I have read the appropriate section of the manual on tabulated interactions, and am working on building an appropriate topology. IU am assuming that my coarse-grained particles consist of the silica particles plus a number of surface counterions in such a way that the particles will be electrically neutral, so presumably I can set all the entries in the tabulated potential file for the Coulomb terms to zero. If my understanding of the manual is correct, I can then introduce my PMF potential in one or other of the g() and h() columns, along with its appropriate derivatives. The plan is to randomly place a number of the coarse-grained particles in a simulation box and choose an appropriate time step and thermostat to run aggregation simulations. Some of the literature I have read suggests timesteps of around 10^-6 s and Brownian dynamics - can anyone comment on the advisability of these choices? I anticipate significant aggregation within ~seconds. Another issue is whether or not to include coarse-grained water and explicit ions in the simulation box. Recent postings on this list have suggested that Lagrangian dynamics should not be done in a vacuum, so presumably the same is true for Brownian dynamics? Has anyone on the list used 'coarse-grained' water in their GROMACS simulations (references needed)? If I have extracted different tabulated potentials for each background counterion concentration, this information is presumably 'built in' to my extracted PMF interparticle potentials, and I shouldn't need to include explicit counterions in my coarse-grained simulation, correct? Because of the way you did things, everything is included in the PMF, both water and counterions. If this is an accurate approach is a completely different matter. Now you have the 2-body coarse-grained term correct, but there could of course be multi-body non-additive effects. For such effects you might need coarse-grained water or counterions, but then things get much more complicated. I'm uncertain that if my particles have a persistent (though fluctuating) electric dipole moment throughout the PMF runs, this will also be fully accounted for by the PME electrostatics used during my PMF runs. I have
[gmx-users] Coarse-graining and tabulated non-bonded potentials - will write up on the wiki
Hello all, Firstly, many thanks to everyone who has contributed useful advice to me over a number of years using GROMACS. I have performed a large number of potential of mean force calculations for the forces acting between two approximately spherical amorphous silica particles (various sizes 5 nm diameter) in TIP4P water with PME electrostatics and varying concentrations of background ions. Now I want to 'coarse-grain' the simulation, treating each silica particle as a single point mass, and use the interaction potential between the particles obtained from the PMF results as a tabulated potential in mdrun, to allow longer time- and size-scales to be investigated (I want to examine colloidal aggregation behaviour of these particles). I have tabulated the PMF potential and its derivatives as a function of centre-of-mass separation of the particles as suggested in the manual (but I can only tabulate for COM-COM distances greater than the contact distance [ = 2r in the hard-sphere approximation] out to some cutoff). Will I have to add a short-distance 'hard-sphere' wall to my tabulated potentials? I have read the appropriate section of the manual on tabulated interactions, and am working on building an appropriate topology. IU am assuming that my coarse-grained particles consist of the silica particles plus a number of surface counterions in such a way that the particles will be electrically neutral, so presumably I can set all the entries in the tabulated potential file for the Coulomb terms to zero. If my understanding of the manual is correct, I can then introduce my PMF potential in one or other of the g() and h() columns, along with its appropriate derivatives. The plan is to randomly place a number of the coarse-grained particles in a simulation box and choose an appropriate time step and thermostat to run aggregation simulations. Some of the literature I have read suggests timesteps of around 10^-6 s and Brownian dynamics - can anyone comment on the advisability of these choices? I anticipate significant aggregation within ~seconds. Another issue is whether or not to include coarse-grained water and explicit ions in the simulation box. Recent postings on this list have suggested that Lagrangian dynamics should not be done in a vacuum, so presumably the same is true for Brownian dynamics? Has anyone on the list used 'coarse-grained' water in their GROMACS simulations (references needed)? If I have extracted different tabulated potentials for each background counterion concentration, this information is presumably 'built in' to my extracted PMF interparticle potentials, and I shouldn't need to include explicit counterions in my coarse-grained simulation, correct? Thank you for reading all the way through this posting. As mentioned in the message title, I will use and write up any advice given to me as a draft example tutorial on the wiki, then more qualified people can correct the (probably numerous) mistakes. Many thanks in advance, Steve Kirk -- Dr. Steven R. Kirk [EMAIL PROTECTED], [EMAIL PROTECTED] Dept. of Technology, Mathematics Computer Science (P)+46 520 223215 University West (F)+46 520 223299 P.O. Box 957 Trollhattan 461 29 SWEDEN http://taconet.webhop.org ___ gmx-users mailing listgmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users Please search the archive at http://www.gromacs.org/search before posting! Please don't post (un)subscribe requests to the list. Use the www interface or send it to [EMAIL PROTECTED] Can't post? Read http://www.gromacs.org/mailing_lists/users.php
[gmx-users] Advice needed on choosing a VdW cutoff
Hello, I would be very grateful for advice on the following system: Consider a pair of spherical macromolecules of diameter ~2.5 nm, arranged along an x-axis so that their centres of mass are 4nm apart, then centered in a 14 x 10 x 10 nm periodic box, so that the minimum distance between either macromolecule's centre of mass and its closest simulation cell walls in the x,y, or z directions (standard orthogonal Cartesian axis set) is 5 nm. (the macromolecules are silica particles with negatively charged surface sites and an equal number of Na+ counterions clustered around them, so each macromolecule-counterion 'bundle' is electrically neutral). The box is then filled with TIP4P water, energy minimized, run with 'all-bonds' position restraints at 300K for 100 ps, then the main MD run begins, reusing velocities from the last step of the position restraints to initialize the new run. Here is the mdp file I am using for the main MD run: title = Yo cpp = /lib/cpp constraints = none integrator = md dt = 0.001; ps ! nsteps = 100 ; total 1000 ps. nstcomm = 1 nstxout = 1 nstvout = 1 nstfout = 0 nstlog = 1 nstenergy = 1 nstlist = 10 ns_type = grid coulombtype = pme rlist = 1.0 rcoulomb= 1.0 rvdw= 1.0 ; Berendsen temperature coupling is on in two groups Tcoupl = berendsen tc-grps = SNP SOL NA+ tau_t = 0.1 0.1 0.1 ref_t = 300 300 300 ; Energy monitoring energygrps = SNP SOL NA+ ; Isotropic pressure coupling is not on Pcoupl = no Pcoupltype = isotropic tau_p = 0.5 compressibility = 4.5e-5 ref_p = 1.0 ; Generate velocites is off at 300 K. gen_vel = no gen_temp= 300.0 gen_seed= 173529 The problem is that when I run this simulation, the expected drift of the macromolecules towards each other does not occur. Assuming that I want to force every part of each macromolecule to 'see' every part of the other, this would suggest a value of rvdw of 6.5 nm, but I have several worries about this: 1. I have never seen a recommended rvdw in this forum over 1.4 nm, in any model system 2. Should I use a standard, switched or other type of vdW cutoff? 3. Should I switch on long-range dispersion corrections (DispCorr = Ener) ? My goal is to keep the periodic box for the advantages of PME, but somehow reassure myself that the macromolecules can 'see' each other via the vdW forces, so that they will drift together (the expected behaviour) over the course of the simulation. I have trawled the mailing lists for advice on this topic - the only directly relevant post I could find involved the drifting apart of membranes over the course of a simulation, and if I remember correctly, the value of 'pme-order' was suggested as a culprit. I use the default value of 'pme-order' in my simulations. Can anyone please advise me as to what to do next? I do not want to abandon the investigation before I eliminate all possibilities that I am doing something stupid with the configuration of the vdW treatment. Maybe I do not have a long enough simulation, but the fact that I actually see *repulsion* suggests something more fundamental is wrong. All helpful suggestions very gratefully recieved ! Regards, Steve Kirk ___ gmx-users mailing listgmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users Please don't post (un)subscribe requests to the list. Use the www interface or send it to [EMAIL PROTECTED] Can't post? Read http://www.gromacs.org/mailing_lists/users.php
[gmx-users] Droplet vs. vesicle vs. big periodic box vs. implicit solvent
Hello, I would be grateful if anyone on the list could give me the benefit of their expertise on the following proposed simulation: System: 2 amorphous approximately spherical silica particles of diameter 2.5 nm, each with a net charge of -9e, with about 5 nm spacing between their centers of mass, surrounded by TIP4P water and enough Na+ ions to neutralize the overall system (probably found near the charged silica surfaces). The system temperature would be set at 300K. The idea is to watch the particles move towards each other (e.g monitor particle COM-COM distance as a function of time) due to VdW forces until stopped by the mutual electrostatic repulsion. To my naive thinking, there are four ways to tackle this: 1. An isolated sphere of water molecules of diameter at least 8 nm, centered on the COM of the 2 particles, enclosing both particles, with a Coulomb cutoff, VdW cutoff and rlist set to 8 nm. Disadvantages: slow Coulomb and VdW calculation, possible water ordering effects (see recent DvS paper) with large cutoff, water boiling off the sphere. 2. Like 1. , only with the outermost water molecules in the droplet frozen in place to form a 'vesicle', hopefully keeping the water from boiling off. 3. Particle pair COM centered in a big periodic box (e.g. 12 x 10 x 10 nm), particle pair aligned along the x-axis, box filled with water. Systems arranged so that each particle is closer to every part of its immediate twin that to any periodic copies. rvdw also set to 8 nm. Advantages: gain efficiencies of PME for Coulomb interaction Disadvantages: Painfully large numbers of atoms, might be possible with large parallel computing power 4. Do away with the explicit water molecules completely, using an implicit solvent (but, presumably, keeping explicit Na+ ions and particles). rvdw set to 8 nm. Any advice on how to do it this way - this is unknown territory for me? Any GROMACS tutorials for implicit solvent or GB out there? Advantages: Small number of atoms Disadvantages: Lose details of particle-water and ion-water interactions (but the latter may not be so important, as polarized water models more realistic than TIP4P for ion-water interactions and hydration shell structure anyway?) Insight, suggestions and additional options very gratefully received. -- Dr. Steven R. Kirk [EMAIL PROTECTED], [EMAIL PROTECTED] Dept. of Technology, Mathematics Computer Science (P)+46 520 223215 University West (F)+46 520 223299 P.O. Box 957 Trollhattan 461 29 SWEDEN http://taconet.webhop.org ___ gmx-users mailing listgmx-users@gromacs.org http://www.gromacs.org/mailman/listinfo/gmx-users Please don't post (un)subscribe requests to the list. Use the www interface or send it to [EMAIL PROTECTED] Can't post? Read http://www.gromacs.org/mailing_lists/users.php
