I ran the minimisation, and mdrun gave the following last few lines of output:<br><br>==================================<br>Step=19990, Dmax= 7.8e-06 nm, Epot= -7.14922e+04 Fmax= 3.19585e+04, atom= 2395<br>Step=19992, Dmax= 4.7e-06 nm, Epot= -7.14924e+04 Fmax= 8.72456e+03, atom= 2395<br>
Step=19993, Dmax= 5.6e-06 nm, Epot= -7.14927e+04 Fmax= 3.70693e+04, atom= 2395<br>Step=19994, Dmax= 6.8e-06 nm, Epot= -7.14929e+04 Fmax= 2.10513e+04, atom= 2395<br>Step=19995, Dmax= 8.1e-06 nm, Epot= -7.14929e+04 Fmax= 4.71234e+04, atom= 2395<br>
Step=19996, Dmax= 9.7e-06 nm, Epot= -7.14931e+04 Fmax= 3.68827e+04, atom= 2395<br>Step=19998, Dmax= 5.8e-06 nm, Epot= -7.14933e+04 Fmax= 1.35461e+04, atom= 2395<br>Step=19999, Dmax= 7.0e-06 nm, Epot= -7.14933e+04 Fmax= 4.77130e+04, atom= 2395<br>
Step=20000, Dmax= 8.4e-06 nm, Epot= -7.14936e+04 Fmax= 2.41860e+04, atom= 2395<br><br>writing lowest energy coordinates.<br><br>Steepest Descents did not converge to Fmax < 10 in 20001 steps.<br>Potential Energy = -7.1493609e+04<br>
Maximum force = 2.4185994e+04 on atom 2395<br>Norm of force = 8.1511212e+02<br>==================================<br><br>As you can be seen, the forces still do not converge to Fmax < 10, even after 20,000 steps.<br>
<br>Does anyone know what the problem might be?<br><br>Thanks,<br><br>Nancy<br><br><br><div class="gmail_quote">On Thu, Aug 6, 2009 at 6:48 PM, Bruce D. Ray <span dir="ltr"><<a href="mailto:brucedray@yahoo.com">brucedray@yahoo.com</a>></span> wrote:<br>
<blockquote class="gmail_quote" style="border-left: 1px solid rgb(204, 204, 204); margin: 0pt 0pt 0pt 0.8ex; padding-left: 1ex;"><div><div style="font-family: times new roman,new york,times,serif; font-size: 10pt;"><div><div>
</div><div class="h5"><span style="font-family: times new roman,new york,times,serif;">On </span><font style="font-family: times new roman,new york,times,serif;" face="Tahoma" size="2">Thursday, August 6, 2009 at 4:08:18 PM, </font><font face="Tahoma" size="2"><span style="font-family: times new roman,new york,times,serif;">Nancy <<a href="mailto:nancy5villa@gmail.com" target="_blank">nancy5villa@gmail.com</a>> wrote:</span><br>
</font></div></div><div><div><div></div><div class="h5"> <br>> I have attempted to perform energy minimisation from scratch again, and these are the<br>> commands that I am using to do so:<br>> <br>> $ .../topolbuild1_2_1/src/topolbuild -n ethanediol -dir .../topolbuild1_2_1/dat/gromacs -ff gmx53a6<br>
> <br>> which generates the files:<br>> <br>> ethanediol.gro<br>> ethanediol.log<br>> ethanediolMOL.mol2<br>> ethanediol.top<br>> ffethanediol.itp<br>> posreethanediol.itp<br>> <br>> in the "ethanediol.log" file I noticed the following lines with asterisks:<br>
> <br>> ========================================<br>> Angles Force Field Results<br>> Angle Atoms force angle method measured<br>> 1 H12- C1- H11 ****** ****** ****** 109.403<br>
> 2 C2- C1- H11 ****** ****** ****** 109.456<br>>
3 O1- C1- H11 ****** ****** ****** 109.438<br>> 4 C2- C1- H12 ****** ****** ****** 109.537<br>> 5 O1- C1- H12 ****** ****** ****** 109.468<br>> 6 O1- C1- C2 320.0 109.50 1 109.526<br>
>
7 C1- C2- O2 320.0 109.50 1 109.526<br>> 8 C1- C2- H22 ****** ****** ****** 109.453<br>> 9 C1- C2- H21 ****** ****** ****** 109.537<br>> 10 H22- C2- O2 ****** ****** ****** 109.439<br>
>
11 H21- C2- O2 ****** ****** ****** 109.468<br>> 12 C2- O2- HO2 450.0 109.50 1 106.864<br>> 13 H21- C2- H22 ****** ****** ****** 109.404<br>> 14 C1- O1- HO1 450.0 109.50 1 106.864<br>
>
========================================<br><br></div></div>These are fine. Only the nonpolar hydrogens lack parameters and they are removed from<br>consideration in the final topology because this is a united atoms model force field.<div class="im">
<br><br>> I then enlarge the box and solvate the molecule:<br>> <br>> $ editconf -f ethanediol.gro -o ethanediol_box.gro -box 3 3 3 -c<br>> <br>> $ genbox -cp ethanediol_box.gro -cs spc216.gro -o ethanediol_solv.gro -p ethanediol.top -box 3 3 3<br>
> <br>> There are 2679 atoms (891 waters and 1 ethanediol). I then configure the minimisation with grompp:<br>> <br>> $ grompp -f minim.mdp -c ethanediol_solv.gro -p ethanediol.top -o em.tpr<br>> <br>> using the following .mdp file:<br>
> <br>> ========================================<br>> define = -DFLEXIBLE<br><br></div>Just out of curiosity, why are you defining flexible waters at this point?<div><div></div><div class="h5"><br><br>> integrator = steep<br>
> emtol = 10.0<br>> emstep = 0.1<br>> nsteps = 15000<br>> <br>> nstlist = 1<br>> ns_type = grid<br>>
rlist = 1.0<br>> coulombtype = PME<br>> rcoulomb_switch = 1.0<br>> rvdw_switch = 1.3<br>> pbc = xyz<br>> <br>> pme_order = 4<br>> constraints = none<br>> <br>> nstxout = 1<br>
> nstvout = 1<br>> nstenergy = 1<br>>
nstlog = 1<br>> nstcomm = 1<br>> Tcoupl = no<br>> Pcoupl = no<br>> gen_vel = no<br>> ========================================<br>> <br>> I then run the simulation for 15,000 steps:<br>
> <br>> $ mdrun -v -deffnm em<br>
> <br>> the last few lines of output are as follows:<br>> <br>> ========================================<br>> Step=14991, Dmax= 5.3e-06 nm, Epot= -6.24538e+04 Fmax= 2.50258e+04, atom= 2397<br>> Step=14992, Dmax= 6.4e-06 nm, Epot= -6.24539e+04 Fmax= 3.52982e+04, atom= 1171<br>
>
Step=14993, Dmax= 7.6e-06 nm, Epot= -6.24540e+04 Fmax= 3.71731e+04, atom= 2397<br>> Step=14994, Dmax= 9.2e-06 nm, Epot= -6.24540e+04 Fmax= 4.93409e+04, atom= 1171<br>> Step=14995, Dmax= 1.1e-05 nm, Epot= -6.24541e+04 Fmax= 5.53109e+04, atom= 2397<br>
>
Step=14997, Dmax= 6.6e-06 nm, Epot= -6.24546e+04 Fmax= 8.50832e+03, atom= 1171<br>> Step=14998, Dmax= 7.9e-06 nm, Epot= -6.24547e+04 Fmax= 6.37691e+04, atom= 2397<br>> Step=14999, Dmax= 9.5e-06 nm, Epot= -6.24553e+04 Fmax= 2.67488e+04, atom= 1171<br>
>
Step=15000, Dmax= 1.1e-05 nm, Epot= -6.24548e+04 Fmax= 8.09451e+04, atom= 2397<br>> writing lowest energy coordinates.<br>> <br>> Steepest Descents did not converge to Fmax < 10 in 15001 steps.<br>> Potential Energy = -6.2455254e+04<br>
>
Maximum force = 2.6748828e+04 on atom 1171<br>> Norm of force = 9.5851715e+02<br>> ========================================<br>> <br>> After
viewing the trajectory with ngmx, I noticed that most of the motion of
the water<br>> molecules occurs within the first 50 picoseconds.<br>
> <br>> Please advise on how to further minimise this system.<br><br><br></div></div>I downloaded NSC93876 ethylene glycol as a Sybyl (.sy2) file from NCI, processed it with<br>dos2unix, added the missing@<TRIPOS> last line, and named it ethylene_glycol.mol2<br>
Charge information is absent on this file as downloaded, and I did not do anything to<br>account for charges. This may make my results less than useful. I then ran the following commands:<br><br> /usr/bdray/bin/topolbuild1_3 -dir /usr/bdray/tables -ff gmx53a6 -n glycol/ethylene_glycol<br>
editconf -f ethylene_glycol.gro -o ethanediol_box.gro -box 3 3 3 -c<div class="im"><br> genbox -cp ethanediol_box.gro -cs spc216.gro -o ethanediol_solv.gro -p ethanediol.top -box 3 3 3<br></div> grompp -f em.mdp -c ethanediol_solv.gro -p ethanediol.top -o em.tpr<br>
<br>File em.mdp includes:<br><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"></span><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> constraints = none</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;">
<span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> integrator = steep</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;"><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> nsteps = 20000</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;">
<span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> ;</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;"><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> ; Energy minimizing stuff</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;">
<span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> ;</span><div class="im"><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;"><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> emtol = 10.0</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;">
<span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> emstep =
0.01</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;"><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> </span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;">
</div><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> nstlist = 1</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;"><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> coulombtype = PME</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;">
<span style="font-family: Courier New,courier,monaco,monospace,sans-serif;">
nstcomm = 1</span><div class="im"><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;"><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> ns_type = grid</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;">
<span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> rlist = 1</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;"></div><div class="im"><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> rcoulomb = 1.0</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;">
<span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> rvdw = 1.3</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;"></div><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> nstxout = 1</span><div class="im">
<br style="font-family: Courier New,courier,monaco,monospace,sans-serif;"><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> pbc = xyz</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;">
<span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> pme_order = 4</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;"></div><span style="font-family: times new roman,new york,times,serif;"><br>
I deliberately did not use the flexible waters model with this energy minimization<br>I then ran<br><br></span> mdrun -v -deffnm em<br><br>I got the result:<br><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> Stepsize too small, or no change in energy.</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;">
<span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> </span><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;">Converged to machine precision,</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;">
<span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> </span><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;">but not to the requested precision Fmax < 10</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;">
<span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> </span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;"><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> </span><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;">Double precision normally gives you higher accuracy.</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;">
<span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> </span><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;">You might need to increase your constraint accuracy, or turn</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;">
<span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> </span><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;">off constraints alltogether (set constraints = none in mdp file)</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;">
<span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> </span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;"><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> </span><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;">Steepest Descents converged to machine precision in 6291
steps,</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;"><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> </span><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;">but did not reach the requested Fmax < 10.</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;">
<span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> </span><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;">Potential Energy = -4.9487391e+04</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;">
<span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> </span><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;">Maximum force = 9.0675652e+01 on atom 5</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;">
<span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> </span><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;">Norm of force = 2.5599347e+02</span><br style="font-family: Courier New,courier,monaco,monospace,sans-serif;">
<span style="font-family: times new roman,new york,times,serif;"><br>I note that atom 5 is part of the ethylene glycol molecule. Furthermore, with very few<br>exceptions, most of the activity in the energy minimization was devoted to the ethylene<br>
glycol molecule. I presume that I would reach the desired Fmax if I were to use double<br>precision. However, an Fmax < 90.7 would seem to me to be usable.<br><br>Perhaps the problem stems from the use of flexible waters in the energy minimization?<br>
</span><br><br></div><div class="im">-- <br>Bruce D. Ray,
Ph.D.<br>Associate Scientist<br></div><div class="im">IUPUI<br>Physics Dept.<br>402 N. Blackford St.<br>Indianapolis, IN 46202-3273<div><br></div></div></div><br>
</div><br>_______________________________________________<br>
gmx-users mailing list <a href="mailto:gmx-users@gromacs.org">gmx-users@gromacs.org</a><br>
<a href="http://lists.gromacs.org/mailman/listinfo/gmx-users" target="_blank">http://lists.gromacs.org/mailman/listinfo/gmx-users</a><br>
Please search the archive at <a href="http://www.gromacs.org/search" target="_blank">http://www.gromacs.org/search</a> before posting!<br>
Please don't post (un)subscribe requests to the list. Use the<br>
www interface or send it to <a href="mailto:gmx-users-request@gromacs.org">gmx-users-request@gromacs.org</a>.<br>
Can't post? Read <a href="http://www.gromacs.org/mailing_lists/users.php" target="_blank">http://www.gromacs.org/mailing_lists/users.php</a><br></blockquote></div><br>