[gmx-users] Instantaneous dipole moment of water molecule
Dear All I want to calculate instantaneous dipole moment of a water molecule(from liquid water) in gromacs.I had gone through the archives of the forum but couldn't found much. I used command g_dipoles with following extension g_dipoles -f 293.trr -s 293.tpr -n atoms.ndx -o Mtot.xvg -mu -1 atoms.ndx contains three atoms of a water molecule. which gave the output as Dipole moment (Debye) - Average = 2.2740 Std. Dev. = 0.0007 Error = 0. The following averages for the complete trajectory have been calculated: Total M_x = -0.00546922 Debye Total M_y = 0.0192573 Debye Total M_z = -0.0446763 Debye Total M_x^2 = 1.68877 Debye^2 Total M_y^2 = 1.71468 Debye^2 Total M_z^2 = 1.76741 Debye^2 Total |M|^2 = 5.17086 Debye^2 Total |M| ^2 = 0.00239673 Debye^2 |M|^2 - |M| ^2 = 5.16846 Debye^2 Finite system Kirkwood g factor G_k = 0.999537 Infinite system Kirkwood g factor g_k = 0.978649 Epsilon = 1.06689 I am convinced with the value in comparison to that reported in literature. But when I plot the Mtot.xvg I found the fluctuations in dipole from -2.5 to 2.5 D. Is it reasonable to have this kind of behaviour or is it due to change in the direction of dipole vector. I need instantaneous fluctuations to calculate dipole correlation function. Further to calculate dipole moment of water I believe I should use Dflexible in mdp file,for without fluctuations the dipole for rigid water molecule must not change .But It has been strongly demotivated in previous discussions.Please help me in this regard. -- DeepaK Ojha School Of Chemistry Selfishness is not living as one wishes to live, it is asking others to live as one wishes to live -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] strange continued jobs
HI Justin: thanks for such kind comments. I restart the job and it works fine now. best Albert On 09/11/2012 12:50 PM, Justin Lemkul wrote: On 9/11/12 3:11 AM, Albert wrote: hello: I've restart my jobs with command: mpiexec -n 160 /opt/gromacs/4.5.5/bin/mdrun -nosum -dlb yes -v -s md.tpr -npme 20 -o md2.trr -g md2.log -e md2.edr -cpi The previous log file is: log---one--- DD step 11346 vol min/aver 0.540! load imb.: force 40.5% pme mesh/force 0.629 Step Time Lambda 11347 226940.00.0 Energies (kJ/mol) AngleU-BProper Dih. Improper Dih. CMAP Dih. 1.70355e+024.99014e+043.43419e+048.47175e+02 -1.29907e+03 LJ-14 Coulomb-14LJ (SR)LJ (LR) Disper. corr. 1.00912e+04 -5.66686e+031.58148e+04 -4.31489e+02 -4.61752e+03 Coulomb (SR) Coul. recip. PotentialKinetic En. Total Energy -4.29387e+05 -9.64431e+04 -4.26678e+051.22078e+05 -3.04600e+05 Temperature Pres. DC (bar) Pressure (bar) Constr. rmsd 3.03352e+02 -1.64273e+026.28299e+011.75083e-05 - After I continued the job, I found my log file is the following: --long--two DD step 8675 vol min/aver 0.450! load imb.: force 38.0% pme mesh/force 0.771 Step Time Lambda 8676 173520.00.0 Energies (kJ/mol) AngleU-BProper Dih. Improper Dih. CMAP Dih. 1.53725e+025.06175e+043.48843e+048.51459e+02 -1.35598e+03 LJ-14 Coulomb-14LJ (SR)LJ (LR) Disper. corr. 1.02049e+04 -6.48180e+031.63099e+04 -4.31662e+02 -4.62554e+03 Coulomb (SR) Coul. recip. PotentialKinetic En. Total Energy -4.29154e+05 -9.60364e+04 -4.25064e+051.21298e+05 -3.03766e+05 Temperature Pres. DC (bar) Pressure (bar) Constr. rmsd 3.01413e+02 -1.64844e+02 -7.07993e+011.74002e-05 - obviously, the continued job didn't wirte log based on previous time. What's present in the .log file is rather irrelevant. What's important is the contents of the .cpt file that was used as input for -cpi. Use gmxcheck to find out what time from which the simulation should have continued. It seems odd indeed that you would have gone back 50 ns, especially since the .cpt file should be written every 15 minutes. I am just wondering is it OK for this continued job when I am trying to merged it with previous after it is finished? Given the information at hand, there's nothing immediately wrong about what you're seeing. Inconvenient, perhaps, to have lost time, but the contents of the .cpt file used will confirm what should have happened. Also I found that the trjconv command doesn't work when I try to extract trajectory for analyis: trjconv -f md2.trr -s md.tpr -pbc mol -o md2.xtc -error- Select a group: 0 Selected 0: 'System' trn version: GMX_trn_file (single precision) Reading frame 0 time 2000.000 Setting output precision to 0.001 (nm) --- Program trjconv, VERSION 4.5.5 Source code file: futil.c, line: 459 File input/output error: md2.xtc For more information and tips for troubleshooting, please check the GROMACS website at http://www.gromacs.org/Documentation/Errors --- Don't Push Me, Cause I'm Close to the Edge (Grandmaster Flash) Do you have sufficient disk space? Read/write permissions? -Justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] dihedraltypes amber forcefield
Hi everybody, I get the error that something is wrong with the dihedraltype of 4 atoms in my protein when I use the grompp command. I looked at the topology file which atoms are involved in this dehedraltype. Those are the four atoms: CA CA C OS Where can I find out the parameter for this dihedraltype so that I can add it to the ffbonded file? Thank you -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] Orientation of protein
Dear users, As a matter of fact I don't know the exact orientation of my protein. To insert a protein in a lipid bilayer, I need to set the orientation of protein and bilayer parallel to z-axis. I'd like to know if it is necessary for protein to be in z-direction exactly? Is it acceptable to set the protein's direction approximately? Thanks in advance. Sincerely, Shima -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] No molecules were defined in the system/ No such molecule type Dio
On 9/12/12 1:58 AM, sarah k wrote: Dear all, I'm trying to simulate several drug-protein complexes on a clustered server. I used PRODRG and editted my .top, .pdb, .gro and posre.itp files. Here is the final lines of my .top file: ; Include Position restraint file ; Include Position restraint file #ifdef POSRES #include posre.itp #endif ; Include water topology #include spc.itp ; Include Dio topology #include dio.top #ifdef POSRES_WATER ; Position restraint for each water oxygen [ position_restraints ] ; i funct fcxfcyfcz 11 1000 1000 1000 #endif ; Include generic topology for ions #include ions.itp [ system ] ; Name Nona in water [ molecules ] ; Compound#mols Protein_X 1 Dio 1 SOL 63598 The statements are closed with #endif. Yet I recieve No molecules were defined in the system. For the same file I sometimes get No such molecule type Dio in the first grompp prompt. What's the solution? Your topology is misformatted somewhere, either with an #ifdef/#endif block not being closed properly or a naming inconsistency somewhere. The error is not apparent from the snippet you've shown, so the best any of us can do is guess at this point. Beware that PRODRG topologies require considerable manual adjustment to be considered sufficiently accurate. -Justin -- Justin A. Lemkul, Ph.D. Research Scientist Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] dihedraltypes amber forcefield
On 9/12/12 3:29 AM, reising...@rostlab.informatik.tu-muenchen.de wrote: Hi everybody, I get the error that something is wrong with the dihedraltype of 4 atoms in my protein when I use the grompp command. I looked at the topology file which atoms are involved in this dehedraltype. Those are the four atoms: CA CA C OS Where can I find out the parameter for this dihedraltype so that I can add it to the ffbonded file? Google search, published literature, or derive the parameters yourself in accordance with prescribed force field methodology. -Justin -- Justin A. Lemkul, Ph.D. Research Scientist Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] Orientation of protein
On 9/12/12 4:28 AM, Shima Arasteh wrote: Dear users, As a matter of fact I don't know the exact orientation of my protein. To insert a protein in a lipid bilayer, I need to set the orientation of protein and bilayer parallel to z-axis. I'd like to know if it is necessary for protein to be in z-direction exactly? Is it acceptable to set the protein's direction approximately? If you simulate for sufficiently long, likely the protein will adopt a reasonable orientation, provided that whatever you start with isn't too terribly far from reality. That depends, of course, on the protein and whether or not forces like tilting and hydrophobic mismatch will play a significant role in the dynamics. -Justin -- Justin A. Lemkul, Ph.D. Research Scientist Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] LINCS warning in md run
Hi everybody, during the mdrun_mpi I get many LINCS warnings like for example: Step 1143, time 1.143 (ps) LINCS WARNING relative constraint deviation after LINCS: rms 0.265373, max 1.416702 (between atoms 976 and 979) bonds that rotated more than 30 degrees: atom 1 atom 2 angle previous, current, constraint length 976979 90.00.1072 0.2320 0.0960 But I don't know why. I already minimized my protein. And also in the grompp run there were no error messages. This was the end of my grompp output: processing index file... Making dummy/rest group for Acceleration containing 80212 elements Making dummy/rest group for Freeze containing 51274 elements Making dummy/rest group for VCM containing 80212 elements Number of degrees of freedom in T-Coupling group System is 102623.00 Making dummy/rest group for User1 containing 80212 elements Making dummy/rest group for User2 containing 80212 elements Making dummy/rest group for XTC containing 80212 elements Making dummy/rest group for Or. Res. Fit containing 80212 elements Making dummy/rest group for QMMM containing 80212 elements T-Coupling has 1 element(s): System Energy Mon. has 2 element(s): Protein non-Protein Acceleration has 1 element(s): rest Freeze has 2 element(s): Protein__!TYP rest User1has 1 element(s): rest User2has 1 element(s): rest VCM has 1 element(s): rest XTC has 1 element(s): rest Or. Res. Fit has 1 element(s): rest QMMM has 1 element(s): rest Checking consistency between energy and charge groups... Estimate for the relative computational load of the PME mesh part: 0.32 writing run input file... Can you please help me? Thank you, Eva -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] Published articles which have used Gromcas for their work
Dear All Actually I am very interested to know the list of works which have used Gromacs for their work. I searched in both google scholar and scopus, but the results are not the same! of course I looked for all articles that have cited principal papers in Gromacs in: http://www.gromacs.org/About_Gromacs/Citations in the result of my search there are some articles which has not used gromacs but because of some reasons it has cited Gromacs. Please let me know how can I have a list of all publications by Gromacs. Thanks in advance for your guidances Best Mohsen Ramezanpour -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] LINCS warning in md run
On 9/12/12 6:57 AM, reising...@rostlab.informatik.tu-muenchen.de wrote: Hi everybody, during the mdrun_mpi I get many LINCS warnings like for example: Step 1143, time 1.143 (ps) LINCS WARNING relative constraint deviation after LINCS: rms 0.265373, max 1.416702 (between atoms 976 and 979) bonds that rotated more than 30 degrees: atom 1 atom 2 angle previous, current, constraint length 976979 90.00.1072 0.2320 0.0960 But I don't know why. I already minimized my protein. And also in the grompp run there were no error messages. This was the end of my grompp output: processing index file... Making dummy/rest group for Acceleration containing 80212 elements Making dummy/rest group for Freeze containing 51274 elements Making dummy/rest group for VCM containing 80212 elements Number of degrees of freedom in T-Coupling group System is 102623.00 Making dummy/rest group for User1 containing 80212 elements Making dummy/rest group for User2 containing 80212 elements Making dummy/rest group for XTC containing 80212 elements Making dummy/rest group for Or. Res. Fit containing 80212 elements Making dummy/rest group for QMMM containing 80212 elements T-Coupling has 1 element(s): System Energy Mon. has 2 element(s): Protein non-Protein Acceleration has 1 element(s): rest Freeze has 2 element(s): Protein__!TYP rest User1has 1 element(s): rest User2has 1 element(s): rest VCM has 1 element(s): rest XTC has 1 element(s): rest Or. Res. Fit has 1 element(s): rest QMMM has 1 element(s): rest Checking consistency between energy and charge groups... Estimate for the relative computational load of the PME mesh part: 0.32 writing run input file... Can you please help me? http://www.gromacs.org/Documentation/Terminology/Blowing_Up Without a complete description of your system (its contents, previous minimization and equilibration and their success/failure) and any relevant .mdp file(s), there's little else that anyone can say. -Justin -- Justin A. Lemkul, Ph.D. Research Scientist Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] Re: freeze value
On 9/12/12 4:52 AM, samira ansari wrote: Dear justin, I have a problem with freezing group. I got no answer from gmx_users. I don't recall seeing the original post. In any case, please keep any Gromacs-related correspondence on the mailing list; I am not a private help service. I am CC'ing this message to the list and ask that anything further be posted there. I made three group Fflex rigid1 rigid2 in my index I changed the fx, fy, fz of atoms belong to rigid1 and rigid2 in .itp file after pdb2gmx. Assuming these are x,y,z force constants in posre.itp, they are irrelevant. Position restraints and freezing are different concepts and different algorithms entirely. also I made this part in em.mdp file: ; freezing group energygrp_excl = rigid1 rigid1 rigid1 SOL rigid2 rigid2 rigid2 SOL ! To remove ;computation of nonbonding interactions between the frozen groups with each other and surroundings (i.e. the solvent, SOL) An exclamation point on this line may cause problems. freezegrps = rigid1 rigid2 freezdim= Y Y Y Y Y Y but during grompp it shows fatal error: Invalid Freezing input: 2 groups and 0 freeze values what should I do to solve this problem? You may have some issue with line endings. Make sure you are always using a plain text editor and make use of dos2unix if necessary. -Justin -- Justin A. Lemkul, Ph.D. Research Scientist Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] Smaller Area Per Lipid for DPPC Bilayer
Hello, I have been basing some DPPC bilayer simulations off of files from Justin Lemkul's tutorial, including the .itp files and .mdp files. Everything has been working fine except that my area/lipid seems to be too low and my diffusion coefficient seems to be too slow compared to experimental values. As a test, I just started with Tieleman's equilibrated 128 DPPC bilayer system, including the waters, and ran a simulation using the mdp file below (note though I selected continuation=yes, this was in fact not continued from a previous equilibration). The simulation has been running for ~75 ns so far, and the area/lipid is on average ~.61-.62 nm^2 . When I do full temperature/pressure equilibrations, even using different thermostats/barostats, I seem to get area/lipid values similar to these. Also, my diffusion coefficients are smaller than those reported in papers invovling DPPC bilayers. I was wondering what the possible reasons for this could be. Any help you could provide would be great. Thanks, David ; Run parameters integrator = md; leap-frog integrator nsteps = 5000 dt = 0.002 ; 2 fs ; Output control nstxout = 5000 ; save coordinates every 2 ps nstvout = 5000 ; save velocities every 2 ps nstxtcout = 5000 ; xtc compressed trajectory output every 2 ps nstenergy = 5000 ; save energies every 2 ps nstlog = 5000 ; update log file every 2 ps ; Bond parameters continuation= yes ; Restarting after NPT constraint_algorithm = lincs; holonomic constraints constraints = all-bonds ; all bonds (even heavy atom-H bonds) constrained lincs_iter = 1 ; accuracy of LINCS lincs_order = 4 ; also related to accuracy ; Neighborsearching ns_type = grid ; search neighboring grid cels nstlist = 5 ; 10 fs rlist = 1.2 ; short-range neighborlist cutoff (in nm) rcoulomb= 1.2 ; short-range electrostatic cutoff (in nm) rvdw= 1.2 ; short-range van der Waals cutoff (in nm) ; Electrostatics coulombtype = PME ; Particle Mesh Ewald for long-range electrostatics pme_order = 4 ; cubic interpolation fourierspacing = 0.12 ; grid spacing for FFT ; Temperature coupling is on tcoupl = Nose-Hoover ; Less accurate thermostat tc-grps = DPPC SOL ; three coupling groups - more accurate tau_t = 0.1 0.1 ; time constant, in ps ref_t = 323 323 ; reference temperature, one for each group, in K ; Pressure coupling is on pcoupl = Parrinello-Rahman ; Pressure coupling on in NPT pcoupltype = semiisotropic ; uniform scaling of x-y box vectors, independent z tau_p = 1.0 ; time constant, in ps ref_p = 0.0 1.0 ; reference pressure, x-y, z (in bar) compressibility = 4.5e-5 4.5e-5 ; isothermal compressibility, bar^-1 ; Periodic boundary conditions pbc = xyz ; 3-D PBC ; Dispersion correction DispCorr= EnerPres ; account for cut-off vdW scheme ; Velocity generation gen_vel = no; Velocity generation is off ;dihre = yes ;dihre_fc = 100 ;dihre_tau = 0.0 ;nstdihreout =50 ; COM motion removal ; These options remove motion of the protein/bilayer relative to the solvent/ions nstcomm = 1 comm-mode = Linear comm-grps = DPPC SOL -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] Smaller Area Per Lipid for DPPC Bilayer
On 9/12/12 10:56 AM, David Ackerman wrote: Hello, I have been basing some DPPC bilayer simulations off of files from Justin Lemkul's tutorial, including the .itp files and .mdp files. Everything has been working fine except that my area/lipid seems to be too low and my diffusion coefficient seems to be too slow compared to experimental values. As a test, I just started with Tieleman's How far off are the diffusion constants? I have never had a lot of luck reproducing experimental values, but this may reflect a limitation of the parameter set, simulation length, or both. equilibrated 128 DPPC bilayer system, including the waters, and ran a simulation using the mdp file below (note though I selected continuation=yes, this was in fact not continued from a previous equilibration). The simulation has been running for ~75 ns so far, and the area/lipid is on average ~.61-.62 nm^2 . When I do full That sounds like the expected outcome for this force field. Why do you say that is too low? temperature/pressure equilibrations, even using different thermostats/barostats, I seem to get area/lipid values similar to these. Also, my diffusion coefficients are smaller than those reported in papers invovling DPPC bilayers. I was wondering what the possible reasons for this could be. Any help you could provide would be great. Curiosities in the .mdp file: tcoupl = Nose-Hoover ; Less accurate thermostat tc-grps = DPPC SOL ; three coupling groups - more accurate tau_t = 0.1 0.1 ; time constant, in ps ref_t = 323 323 ; reference temperature, one for each Why is your tau_t so small? Generally one should use 0.5 - 2.0 with Nose-Hoover. group, in K ; Pressure coupling is on pcoupl = Parrinello-Rahman ; Pressure coupling on in NPT pcoupltype = semiisotropic ; uniform scaling of x-y box vectors, independent z tau_p = 1.0 ; time constant, in ps ref_p = 0.0 1.0 ; reference pressure, x-y, z (in bar) Why are you setting zero pressure along the x-y plane? compressibility = 4.5e-5 4.5e-5 ; isothermal compressibility, bar^-1 ; Periodic boundary conditions pbc = xyz ; 3-D PBC ; Dispersion correction DispCorr= EnerPres ; account for cut-off vdW scheme ; Velocity generation gen_vel = no; Velocity generation is off If you are not continuing from a previous run (as you say above) and you are also not generating velocities, you may be delaying equilibration by allowing the initial forces dictate the velocities. I suppose if the run is stable enough, this is not a huge problem, but in general this combination is not recommended. -Justin -- Justin A. Lemkul, Ph.D. Research Scientist Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] Smaller Area Per Lipid for DPPC Bilayer
Hi, The experimental range of diffusion coefficients are quite large for DPPC, plus the force field and simulation parameters can have a large impact upon the diffusion speeds seen in the simulations. We have just published a study comparing force fields for simulating DPPC and POPC membranes and further details on differences in lipid diffusion are provided in the paper: http://pubs.acs.org/doi/abs/10.1021/ct3003157 We did not test this exact set of cut-offs with this force field. However, from the tests we did perform using these Berger DPPC parameters, I expect that the diffusion speeds should fall within the experimental range using this set of cut-offs. As for the area per lipid, what you are seeing is pretty much as I would expect with the 1.2 nm cut-offs and a dispersion correction. If you want a higher area per lipid, you could try removing the dispersion correction or reducing the cut-off (with the dispersion correction, we saw sensible membrane behaviour with 1.0 nm cut-offs). Do be sure to check the lipid diffusion rate is still sensible if you remove the dispersion correction, as it should substantially increase when doing this (see the paper for some more details). Cheers Tom On 12/09/12 16:29, Justin Lemkul wrote: On 9/12/12 10:56 AM, David Ackerman wrote: Hello, I have been basing some DPPC bilayer simulations off of files from Justin Lemkul's tutorial, including the .itp files and .mdp files. Everything has been working fine except that my area/lipid seems to be too low and my diffusion coefficient seems to be too slow compared to experimental values. As a test, I just started with Tieleman's How far off are the diffusion constants? I have never had a lot of luck reproducing experimental values, but this may reflect a limitation of the parameter set, simulation length, or both. equilibrated 128 DPPC bilayer system, including the waters, and ran a simulation using the mdp file below (note though I selected continuation=yes, this was in fact not continued from a previous equilibration). The simulation has been running for ~75 ns so far, and the area/lipid is on average ~.61-.62 nm^2 . When I do full That sounds like the expected outcome for this force field. Why do you say that is too low? temperature/pressure equilibrations, even using different thermostats/barostats, I seem to get area/lipid values similar to these. Also, my diffusion coefficients are smaller than those reported in papers invovling DPPC bilayers. I was wondering what the possible reasons for this could be. Any help you could provide would be great. Curiosities in the .mdp file: tcoupl = Nose-Hoover ; Less accurate thermostat tc-grps = DPPC SOL ; three coupling groups - more accurate tau_t = 0.1 0.1 ; time constant, in ps ref_t = 323 323 ; reference temperature, one for each Why is your tau_t so small? Generally one should use 0.5 - 2.0 with Nose-Hoover. group, in K ; Pressure coupling is on pcoupl = Parrinello-Rahman ; Pressure coupling on in NPT pcoupltype = semiisotropic ; uniform scaling of x-y box vectors, independent z tau_p = 1.0 ; time constant, in ps ref_p = 0.0 1.0 ; reference pressure, x-y, z (in bar) Why are you setting zero pressure along the x-y plane? compressibility = 4.5e-5 4.5e-5 ; isothermal compressibility, bar^-1 ; Periodic boundary conditions pbc = xyz ; 3-D PBC ; Dispersion correction DispCorr= EnerPres ; account for cut-off vdW scheme ; Velocity generation gen_vel = no; Velocity generation is off If you are not continuing from a previous run (as you say above) and you are also not generating velocities, you may be delaying equilibration by allowing the initial forces dictate the velocities. I suppose if the run is stable enough, this is not a huge problem, but in general this combination is not recommended. -Justin -- Dr Thomas Piggot University of Southampton, UK. -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] Calculating tertiary structure rotation
Dear all, I would like to measure the rotation angle of helix movement during MD simulation course. Can you tell me how can I measure and make that rotation angle plot ? Thanks. Rajiv-- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] Smaller Area Per Lipid for DPPC Bilayer
Hi, Thank you for your response. As to my concern about incorrect areas and diffusion, I am basing it off of other papers that simulate DPPC bilayers. For instance, in this paper: http://www.ncbi.nlm.nih.gov/pmc/articles/PMC3251217/figure/F12/ , they simulate a DPPC bilayer with DiI molecules in it. I did the same simulation, but whereas they get APL of ~.64-.65 nm^2, mine are again ~0.03-0.04 nm^2 smaller. Also, in this paper they show that the lipids diffuse ~1.1 nm^2 over the span of 20 ns, whereas I get a much slower rate of traveling ~1.4-1.6 nm^2 over 90 ns. As mentioned in the other response, if I turn off dispersion correction I get higher APL (~.65-.66 nm^2) and diffusion values that more closely match this and other papers. These APLs and diffusion values are similar for some other papers that simulate DPPC bilayers. Is it ok to have ranges this large compared to these other simulations, and does it make physical sense to turn off the dispersion correction for this force field? Thanks for your time, David On Wed, Sep 12, 2012 at 11:29 AM, Justin Lemkul jalem...@vt.edu wrote: On 9/12/12 10:56 AM, David Ackerman wrote: Hello, I have been basing some DPPC bilayer simulations off of files from Justin Lemkul's tutorial, including the .itp files and .mdp files. Everything has been working fine except that my area/lipid seems to be too low and my diffusion coefficient seems to be too slow compared to experimental values. As a test, I just started with Tieleman's How far off are the diffusion constants? I have never had a lot of luck reproducing experimental values, but this may reflect a limitation of the parameter set, simulation length, or both. equilibrated 128 DPPC bilayer system, including the waters, and ran a simulation using the mdp file below (note though I selected continuation=yes, this was in fact not continued from a previous equilibration). The simulation has been running for ~75 ns so far, and the area/lipid is on average ~.61-.62 nm^2 . When I do full That sounds like the expected outcome for this force field. Why do you say that is too low? temperature/pressure equilibrations, even using different thermostats/barostats, I seem to get area/lipid values similar to these. Also, my diffusion coefficients are smaller than those reported in papers invovling DPPC bilayers. I was wondering what the possible reasons for this could be. Any help you could provide would be great. Curiosities in the .mdp file: tcoupl = Nose-Hoover ; Less accurate thermostat tc-grps = DPPC SOL ; three coupling groups - more accurate tau_t = 0.1 0.1 ; time constant, in ps ref_t = 323 323 ; reference temperature, one for each Why is your tau_t so small? Generally one should use 0.5 - 2.0 with Nose-Hoover. group, in K ; Pressure coupling is on pcoupl = Parrinello-Rahman ; Pressure coupling on in NPT pcoupltype = semiisotropic ; uniform scaling of x-y box vectors, independent z tau_p = 1.0 ; time constant, in ps ref_p = 0.0 1.0 ; reference pressure, x-y, z (in bar) Why are you setting zero pressure along the x-y plane? compressibility = 4.5e-5 4.5e-5 ; isothermal compressibility, bar^-1 ; Periodic boundary conditions pbc = xyz ; 3-D PBC ; Dispersion correction DispCorr= EnerPres ; account for cut-off vdW scheme ; Velocity generation gen_vel = no; Velocity generation is off If you are not continuing from a previous run (as you say above) and you are also not generating velocities, you may be delaying equilibration by allowing the initial forces dictate the velocities. I suppose if the run is stable enough, this is not a huge problem, but in general this combination is not recommended. -Justin -- Justin A. Lemkul, Ph.D. Research Scientist Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read
[gmx-users] About Bond breaking and usage of constraints
Dear justin , Thank you fro your Reply I am doing MD for Cyclic poly Peptide With ARG as N-termina and PRO as C-terminal . I have run pdb2gmx it is oak. . I have solvated and added ions successfully But when i Run the Energy Minimization The Bond between N atom of ARG and C atom of PRO is broken . What Should i Do To keep this Bond Through entire EM and MD ? Can i Use Constraints option in topology file if i need to use Constraints option What is the Syntax ? Also The Bond is not constructed between N atom and C atom of terminal ARG and PRO residues in Topology -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] Smaller Area Per Lipid for DPPC Bilayer
Hi, While comparing to other simulations can be useful, I would argue that the real test for the combination of force field and simulation parameters is to determine if the simulated membrane properties compare well to the experimentally determined values. As long as they do this, you can argue that your choices are sensible. Obviously if you are including proteins (or other molecules) in the system, the parameters need to also be shown to work well for these too. As for not including the dispersion correction, yes it is fine to do this (if it improves the behaviour of your membrane) as a dispersion correction is most appropriately applied to homogeneous systems and not membranes. Regarding the large range of values seen, I would only be concerned if you are exactly reproducing what other people have done, in terms of force field and simulation parameters used and seeing large differences to what they report. As I mentioned before, fairly small changes in some of these parameters can make some pretty substantial impacts upon the membrane properties. You also should be careful to ensure that the properties of the membranes you are analysing are converged (a block analysis is the way to properly check this). Cheers Tom On 12/09/12 18:21, David Ackerman wrote: Hi, Thank you for your response. As to my concern about incorrect areas and diffusion, I am basing it off of other papers that simulate DPPC bilayers. For instance, in this paper: http://www.ncbi.nlm.nih.gov/pmc/articles/PMC3251217/figure/F12/ , they simulate a DPPC bilayer with DiI molecules in it. I did the same simulation, but whereas they get APL of ~.64-.65 nm^2, mine are again ~0.03-0.04 nm^2 smaller. Also, in this paper they show that the lipids diffuse ~1.1 nm^2 over the span of 20 ns, whereas I get a much slower rate of traveling ~1.4-1.6 nm^2 over 90 ns. As mentioned in the other response, if I turn off dispersion correction I get higher APL (~.65-.66 nm^2) and diffusion values that more closely match this and other papers. These APLs and diffusion values are similar for some other papers that simulate DPPC bilayers. Is it ok to have ranges this large compared to these other simulations, and does it make physical sense to turn off the dispersion correction for this force field? Thanks for your time, David On Wed, Sep 12, 2012 at 11:29 AM, Justin Lemkul jalem...@vt.edu wrote: On 9/12/12 10:56 AM, David Ackerman wrote: Hello, I have been basing some DPPC bilayer simulations off of files from Justin Lemkul's tutorial, including the .itp files and .mdp files. Everything has been working fine except that my area/lipid seems to be too low and my diffusion coefficient seems to be too slow compared to experimental values. As a test, I just started with Tieleman's How far off are the diffusion constants? I have never had a lot of luck reproducing experimental values, but this may reflect a limitation of the parameter set, simulation length, or both. equilibrated 128 DPPC bilayer system, including the waters, and ran a simulation using the mdp file below (note though I selected continuation=yes, this was in fact not continued from a previous equilibration). The simulation has been running for ~75 ns so far, and the area/lipid is on average ~.61-.62 nm^2 . When I do full That sounds like the expected outcome for this force field. Why do you say that is too low? temperature/pressure equilibrations, even using different thermostats/barostats, I seem to get area/lipid values similar to these. Also, my diffusion coefficients are smaller than those reported in papers invovling DPPC bilayers. I was wondering what the possible reasons for this could be. Any help you could provide would be great. Curiosities in the .mdp file: tcoupl = Nose-Hoover ; Less accurate thermostat tc-grps = DPPC SOL ; three coupling groups - more accurate tau_t = 0.1 0.1 ; time constant, in ps ref_t = 323 323 ; reference temperature, one for each Why is your tau_t so small? Generally one should use 0.5 - 2.0 with Nose-Hoover. group, in K ; Pressure coupling is on pcoupl = Parrinello-Rahman ; Pressure coupling on in NPT pcoupltype = semiisotropic ; uniform scaling of x-y box vectors, independent z tau_p = 1.0 ; time constant, in ps ref_p = 0.0 1.0 ; reference pressure, x-y, z (in bar) Why are you setting zero pressure along the x-y plane? compressibility = 4.5e-5 4.5e-5 ; isothermal compressibility, bar^-1 ; Periodic boundary conditions pbc = xyz ; 3-D PBC ; Dispersion correction DispCorr= EnerPres ; account for cut-off vdW scheme ; Velocity generation gen_vel = no; Velocity generation is off If you are not continuing from a previous run (as you say above) and you are also not generating velocities, you may be delaying equilibration by
Re: [gmx-users] Smaller Area Per Lipid for DPPC Bilayer
On 9/12/12 1:21 PM, David Ackerman wrote: Hi, Thank you for your response. As to my concern about incorrect areas and diffusion, I am basing it off of other papers that simulate DPPC bilayers. For instance, in this paper: http://www.ncbi.nlm.nih.gov/pmc/articles/PMC3251217/figure/F12/ , they simulate a DPPC bilayer with DiI molecules in it. I did the same simulation, but whereas they get APL of ~.64-.65 nm^2, mine are again ~0.03-0.04 nm^2 smaller. Also, in this paper they show that the lipids diffuse ~1.1 nm^2 over the span of 20 ns, whereas I get a much slower rate of traveling ~1.4-1.6 nm^2 over 90 ns. As mentioned in the other response, if I turn off dispersion correction I get higher APL (~.65-.66 nm^2) and diffusion values that more closely match this and other papers. The paper cited above does not report any error estimates for their values (unless I have missed them) and it appears they produced only one trajectory per condition. Multiple simulations and proper statistical measurements should be made. Experimental ranges for DPPC APL (if I recall them properly) are 0.62-0.65 nm^2 at 323K, so everyone seems to be in the right ballpark. Regarding diffusion, that's harder to compare. You're also seeing less diffusion than in Lindahl and Edholm (2001) J. Chem. Phys. 10: 4938-4950 over the same time period (MSD there was on the order of 4 nm^2 by 90 ns). These APLs and diffusion values are similar for some other papers that simulate DPPC bilayers. Is it ok to have ranges this large compared to these other simulations, and does it make physical sense to turn off the dispersion correction for this force field? I think it all amounts to sampling. One trajectory is not definitive. Averaging over several that have been demonstrated to have converged is much more reliable. Lipids take a long time to be happy; 20-40 ns of equilibration time may have to be neglected before claiming equilibrium properties. -Justin -- Justin A. Lemkul, Ph.D. Research Scientist Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] Re: About Bond breaking and usage of constraints
On 9/12/12 1:34 PM, vidhya sankar wrote: Dear justin , Thank you fro your Reply I am doing MD for Cyclic poly Peptide With ARG as N-termina and PRO as C-terminal . I have run pdb2gmx it is oak. . I have solvated and added ions successfully But when i Run the Energy Minimization The Bond between N atom of ARG and C atom of PRO is broken . What Should i Do To keep this Bond Through entire EM and MD ?Can i Use Constraints option in topology file if i need to use Constraints option What is the Syntax ? Also The Bond is not constructed between N atom and C atom of terminal ARG and PRO residues in Topology Then you already have your answer. Bonds do not break or form during MD. pdb2gmx does not produce cyclic polypeptide topologies unless you have made appropriate modifications. As you are observing, there is no such bond and you need to add it (and any other resulting bonded parameters like angles, dihedrals, and impropers) yourself. -Justin -- Justin A. Lemkul, Ph.D. Research Scientist Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] Smaller Area Per Lipid for DPPC Bilayer
Hello, Right I guess my biggest concern was diffusion. I did in fact do 12 simulations of DPPC bilayers for 100 ns each, and still got the aforementioned APL and diffusion. When I turn off the dispersion, I get more appropriate APL and MSD values, that match other papers, even when only looking at one simulation. To me, it does not seem the mdp file I used is able to get more common APL and diffusion values even when averaging over a large number of simulations. As for the pressure in the x/y direction, is it more appropriate to use 1 atm or 0 atm for bilayer simulations? -David On Wed, Sep 12, 2012 at 1:43 PM, Justin Lemkul jalem...@vt.edu wrote: On 9/12/12 1:21 PM, David Ackerman wrote: Hi, Thank you for your response. As to my concern about incorrect areas and diffusion, I am basing it off of other papers that simulate DPPC bilayers. For instance, in this paper: http://www.ncbi.nlm.nih.gov/pmc/articles/PMC3251217/figure/F12/ , they simulate a DPPC bilayer with DiI molecules in it. I did the same simulation, but whereas they get APL of ~.64-.65 nm^2, mine are again ~0.03-0.04 nm^2 smaller. Also, in this paper they show that the lipids diffuse ~1.1 nm^2 over the span of 20 ns, whereas I get a much slower rate of traveling ~1.4-1.6 nm^2 over 90 ns. As mentioned in the other response, if I turn off dispersion correction I get higher APL (~.65-.66 nm^2) and diffusion values that more closely match this and other papers. The paper cited above does not report any error estimates for their values (unless I have missed them) and it appears they produced only one trajectory per condition. Multiple simulations and proper statistical measurements should be made. Experimental ranges for DPPC APL (if I recall them properly) are 0.62-0.65 nm^2 at 323K, so everyone seems to be in the right ballpark. Regarding diffusion, that's harder to compare. You're also seeing less diffusion than in Lindahl and Edholm (2001) J. Chem. Phys. 10: 4938-4950 over the same time period (MSD there was on the order of 4 nm^2 by 90 ns). These APLs and diffusion values are similar for some other papers that simulate DPPC bilayers. Is it ok to have ranges this large compared to these other simulations, and does it make physical sense to turn off the dispersion correction for this force field? I think it all amounts to sampling. One trajectory is not definitive. Averaging over several that have been demonstrated to have converged is much more reliable. Lipids take a long time to be happy; 20-40 ns of equilibration time may have to be neglected before claiming equilibrium properties. -Justin -- Justin A. Lemkul, Ph.D. Research Scientist Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] Smaller Area Per Lipid for DPPC Bilayer
While dispersion correction is a great idea that helps to reduce the impact of the precise choice of cutoff distance on the results, the Berger parameters (and indeed all other parameters) were not developed with the inclusion of dispersion correction and one could argue that it is thus non-optimal to include dispersion correction here... especially since it leads to poorer results. This is separate from the difference between isotropic dispersion correction and a proper PME-type LJ term. Both of which are expected to lead to smaller APLs. When using anisotropic pressure coupling for lipid bilayers, you should use 1 atm in all dimensions. Chris. -- original message -- Right I guess my biggest concern was diffusion. I did in fact do 12 simulations of DPPC bilayers for 100 ns each, and still got the aforementioned APL and diffusion. When I turn off the dispersion, I get more appropriate APL and MSD values, that match other papers, even when only looking at one simulation. To me, it does not seem the mdp file I used is able to get more common APL and diffusion values even when averaging over a large number of simulations. As for the pressure in the x/y direction, is it more appropriate to use 1 atm or 0 atm for bilayer simulations? -David -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] Error regarding topology and coordinate files
Hi I am trying to run a simulation of the self-assembly of a course-grained lipid bilayer using GROMACS. So far, I have made a box of DSPC lipid molecules (7.5 x 7.5 x 7.5, 128 molecules). I have also added water molecules using the genbox command (768 water molecules added), and named the coordinate file as waterbox.gro. For the next step, when I try to run a minimization, I get an error as follows: Fatal error: number of coordinates in coordinate file (waterbox.gro, 2560) does not match topology (dspc.top, 1792) What I understood from this was that the water molecules (768 water molecules = 2560-1792) added by genbox had not been updated to the topology file. However, my topology file reads as follows: #include martini_v2.1.itp #include martini_v2.0_lipids.itp [ system ] DSPC BILAYER SELF-ASSEMBLY [ molecules ] DSPC 128 ; W 768 which indicates that 768 water molecules have, indeed, been added to the topology file- this appears to be fine. Is there anything I'm missing out on here? Thanks a lot! -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] Error regarding topology and coordinate files
On 9/12/12 2:38 PM, Bharath K. Srikanth wrote: Hi I am trying to run a simulation of the self-assembly of a course-grained lipid bilayer using GROMACS. So far, I have made a box of DSPC lipid molecules (7.5 x 7.5 x 7.5, 128 molecules). I have also added water molecules using the genbox command (768 water molecules added), and named the coordinate file as waterbox.gro. For the next step, when I try to run a minimization, I get an error as follows: Fatal error: number of coordinates in coordinate file (waterbox.gro, 2560) does not match topology (dspc.top, 1792) What I understood from this was that the water molecules (768 water molecules = 2560-1792) added by genbox had not been updated to the topology file. However, my topology file reads as follows: #include martini_v2.1.itp #include martini_v2.0_lipids.itp [ system ] DSPC BILAYER SELF-ASSEMBLY [ molecules ] DSPC 128 ; W 768 which indicates that 768 water molecules have, indeed, been added to the topology file- this appears to be fine. Is there anything I'm missing out on here? The W line is commented out, so grompp doesn't read the fact that there are 768 W particles at all. -Justin -- Justin A. Lemkul, Ph.D. Research Scientist Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] LINCS warning in md run
Hi Justin, my mdp file for the md run looks like this: define = -DPOSRES integrator = md dt = 0.001 nsteps = 5000 nstxout = 100 nstvout = 0 nstfout = 0 nstlog = 1000 nstxtcout = 500 nstenergy = 5 energygrps = Protein Non-Protein nstcalcenergy = 5 nstlist = 10 ns-type = Grid pbc = xyz rlist = 0.9 coulombtype = PME rcoulomb= 0.9 rvdw= 0.9 fourierspacing = 0.12 pme_order = 4 ewald_rtol = 1e-5 gen_vel = yes gen_temp= 200.0 gen_seed= constraints = all-bonds tcoupl = V-rescale tc-grps = Protein System__!Protein tau_t = 0.1 0.1 ref_t = 298 298 pcoupl = no freezegrps = Protein__!TYP freezedim = Y Y Y I added a phosphate (TYP) to my protein and now want to do first a minimization (which I already did), then a short md run (for the case that I am in a local minima) and after that again a minimization. Since I only want to minimize the energy of the phosphate I freezed the rest of the protein (Protein__!TYP). Since on my cluster I am not allowed to run long jobs I had to divide the minimization run in several runs. The output of the last run was: Steepest Descents converged to machine precision in 15 steps, but did not reach the requested Fmax 10. Potential Energy = -7.7264606e+05 Maximum force = 1.6227990e+04 on atom 961 Norm of force = 1.2328220e+02 gcq#192: It's So Fast It's Slow (F. Black) Steepest Descents converged to machine precision in 15 steps, but did not reach the requested Fmax 10. Potential Energy = -7.7264606e+05 Maximum force = 1.6227990e+04 on atom 961 Norm of force = 1.2328220e+02 After this minimization run I wanted to do the md run. My protein is a membrane protein with its surrounding membrane. I already did a minimization and md run with this protein but without the phosphate. The only difference to the time when everything worked fine is the phosphate. Does this help you somehow to see the failure? Thank you!! On 9/12/12 6:57 AM, reising...@rostlab.informatik.tu-muenchen.de wrote: Hi everybody, during the mdrun_mpi I get many LINCS warnings like for example: Step 1143, time 1.143 (ps) LINCS WARNING relative constraint deviation after LINCS: rms 0.265373, max 1.416702 (between atoms 976 and 979) bonds that rotated more than 30 degrees: atom 1 atom 2 angle previous, current, constraint length 976979 90.00.1072 0.2320 0.0960 But I don't know why. I already minimized my protein. And also in the grompp run there were no error messages. This was the end of my grompp output: processing index file... Making dummy/rest group for Acceleration containing 80212 elements Making dummy/rest group for Freeze containing 51274 elements Making dummy/rest group for VCM containing 80212 elements Number of degrees of freedom in T-Coupling group System is 102623.00 Making dummy/rest group for User1 containing 80212 elements Making dummy/rest group for User2 containing 80212 elements Making dummy/rest group for XTC containing 80212 elements Making dummy/rest group for Or. Res. Fit containing 80212 elements Making dummy/rest group for QMMM containing 80212 elements T-Coupling has 1 element(s): System Energy Mon. has 2 element(s): Protein non-Protein Acceleration has 1 element(s): rest Freeze has 2 element(s): Protein__!TYP rest User1has 1 element(s): rest User2has 1 element(s): rest VCM has 1 element(s): rest XTC has 1 element(s): rest Or. Res. Fit has 1 element(s): rest QMMM has 1 element(s): rest Checking consistency between energy and charge groups... Estimate for the relative computational load of the PME mesh part: 0.32 writing run input file... Can you please help me? http://www.gromacs.org/Documentation/Terminology/Blowing_Up Without a complete description of your system (its contents, previous minimization and equilibration and their success/failure) and any relevant .mdp file(s), there's little else that anyone can say. -Justin -- Justin A. Lemkul, Ph.D. Research Scientist Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't
Re: [gmx-users] LINCS warning in md run
On 9/12/12 3:39 PM, reising...@rostlab.informatik.tu-muenchen.de wrote: Hi Justin, my mdp file for the md run looks like this: define = -DPOSRES integrator = md dt = 0.001 nsteps = 5000 nstxout = 100 nstvout = 0 nstfout = 0 nstlog = 1000 nstxtcout = 500 nstenergy = 5 energygrps = Protein Non-Protein nstcalcenergy = 5 nstlist = 10 ns-type = Grid pbc = xyz rlist = 0.9 coulombtype = PME rcoulomb= 0.9 rvdw= 0.9 fourierspacing = 0.12 pme_order = 4 ewald_rtol = 1e-5 gen_vel = yes gen_temp= 200.0 gen_seed= constraints = all-bonds tcoupl = V-rescale tc-grps = Protein System__!Protein tau_t = 0.1 0.1 ref_t = 298 298 pcoupl = no freezegrps = Protein__!TYP freezedim = Y Y Y I added a phosphate (TYP) to my protein and now want to do first a minimization (which I already did), then a short md run (for the case that I am in a local minima) and after that again a minimization. Since I only want to minimize the energy of the phosphate I freezed the rest of the protein (Protein__!TYP). Since on my cluster I am not allowed to run long jobs I had to divide the minimization run in several runs. The output of the last run was: Steepest Descents converged to machine precision in 15 steps, but did not reach the requested Fmax 10. Potential Energy = -7.7264606e+05 Maximum force = 1.6227990e+04 on atom 961 Norm of force = 1.2328220e+02 gcq#192: It's So Fast It's Slow (F. Black) Steepest Descents converged to machine precision in 15 steps, but did not reach the requested Fmax 10. Potential Energy = -7.7264606e+05 Maximum force = 1.6227990e+04 on atom 961 Norm of force = 1.2328220e+02 After this minimization run I wanted to do the md run. My protein is a membrane protein with its surrounding membrane. I already did a minimization and md run with this protein but without the phosphate. The only difference to the time when everything worked fine is the phosphate. Does this help you somehow to see the failure? Your minimization is insufficient. You have a maximum force in excess of 16000 on atom 961. Such a large force explains the crash. You should investigate what this atom is and what is pushing on it so hard. -Justin -- Justin A. Lemkul, Ph.D. Research Scientist Department of Biochemistry Virginia Tech Blacksburg, VA jalemkul[at]vt.edu | (540) 231-9080 http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] Smaller Area Per Lipid for DPPC Bilayer
Hi Chris, I'm not so sure as to say that all of the parameters for simulating lipid membranes were not developed with a dispersion correction. The Berger parameterisation used a dispersion correction for the pentadecane simulations which were used to parameterise the tails (although as far as I can tell from the Berger paper, this correction was not used in the DPPC membrane simulations). Furthermore, the GROMOS 53A6 parameters of Kukol were tested using simulations which applied a dispersion correction (although you could argue that these GROMOS parameters were initially developed without this correction) and if I remember correctly the CHARMM27 (but not CHARMM36) lipid parameters were intended to be used with a dispersion correction applied (although these parameters are not for use with NPT simulations). I would still argue that above all else, you should choose parameters that someone has shown to accurately reproduce the experimental membrane properties, irrespective of whether that is the original parameterisation work or not (it may well just be your own simulation tests). The Berger force field is a good example of where this sort of testing/validation has been important. Several papers have shown that PME should be used with this force field and not the direct 1.8 nm coulombic cut-off used by Berger et al. Furthermore, in our work I mentioned before, we show that with a 1.0 nm cut-off and no dispersion correction (so the van der Waals parameters I believe were used in the Berger DPPC simulations) there are several membrane properties that do not match the experimental range. I do agree with you for the example here though, it seems (from the information provided) the dispersion correction should not be included with 1.2 nm cut-offs (and this also agrees with results from three different cut-offs tested with the Berger force field in our work). Cheers Tom On 12/09/12 19:30, Christopher Neale wrote: While dispersion correction is a great idea that helps to reduce the impact of the precise choice of cutoff distance on the results, the Berger parameters (and indeed all other parameters) were not developed with the inclusion of dispersion correction and one could argue that it is thus non-optimal to include dispersion correction here... especially since it leads to poorer results. This is separate from the difference between isotropic dispersion correction and a proper PME-type LJ term. Both of which are expected to lead to smaller APLs. When using anisotropic pressure coupling for lipid bilayers, you should use 1 atm in all dimensions. Chris. -- original message -- Right I guess my biggest concern was diffusion. I did in fact do 12 simulations of DPPC bilayers for 100 ns each, and still got the aforementioned APL and diffusion. When I turn off the dispersion, I get more appropriate APL and MSD values, that match other papers, even when only looking at one simulation. To me, it does not seem the mdp file I used is able to get more common APL and diffusion values even when averaging over a large number of simulations. As for the pressure in the x/y direction, is it more appropriate to use 1 atm or 0 atm for bilayer simulations? -David -- Dr Thomas Piggot University of Southampton, UK. -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
Re: [gmx-users] How to create two parallel wall?
There's a built-in wall option. Check section 7.3.20 of the manual. Otherwise, you could just include coordinates and topologies for walls in your simulation. For your second question, g_density with appropriate index groups. On Wed, Sep 12, 2012 at 4:27 PM, Ali Alizadeh ali.alizadehmoja...@gmail.com wrote: Dear All users How to create two parallel wall? How the component density profiles be measured between two walls? Sincerely -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists -- Alex Marshall M.Sc. Candidate Department of Applied Mathematics The University of Western Ontario -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] Re: How to determin the simulation time
在 2012年9月11日星期二，xiaojing gong 写道： Dear all, Typically, we use over 100 ns to simulate the transport progress. But *How the sufficient number of iterations and these times has been determined? Are there some convergence tests ?* * * *Many thanks!* *YLK* * * ** -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] Re: time constant in V-rescale algorithm
在 2012年9月11日星期二，xiaojing gong 写道： Dear all, When use the V-rescale algorithm , how to choose the time constant value, I saw some choose 0, some choose 0.1 ps. Are there some standard for choosing the time constant? Best YLK -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] Calculating secondary structure movements
Dear all, I would like to measure the rotation angle of helix movement during MD simulation course. Can you tell me how can I measure and make that rotation angle plot ? Thanks. Regards Rajiv -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
[gmx-users] Resuming of calculation from last *.cpt
Dear Gromacs Users! I'm looking for possible way to resume trajectory calculations after that calculations have been stoped. Typically in such cases I start new task using cpt file from incomplete run. This produce new trajectory which start from the last frame of previous run. After that I have to merge both files. Is there any way to go on sinulation in the existed files from the last frame of incomplete run without creating of new ones? Thanks for help James -- gmx-users mailing listgmx-users@gromacs.org http://lists.gromacs.org/mailman/listinfo/gmx-users * Please search the archive at http://www.gromacs.org/Support/Mailing_Lists/Search before posting! * Please don't post (un)subscribe requests to the list. Use the www interface or send it to gmx-users-requ...@gromacs.org. * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
