<html><head><style type="text/css"><!-- DIV {margin:0px;} --></style></head><body><div style="font-family:times new roman,new york,times,serif;font-size:10pt"><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 <nancy5villa@gmail.com> wrote:</span><br></font><div> <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>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.<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>Just out of curiosity, why are you defining flexible waters at this point?<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>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<br> genbox -cp ethanediol_box.gro -cs spc216.gro -o ethanediol_solv.gro -p ethanediol.top -box 3 3 3<br> 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><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;"><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><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;"><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;"><span style="font-family: Courier New,courier,monaco,monospace,sans-serif;"> nstxout = 1</span><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;"><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>-- <br>Bruce D. Ray,
Ph.D.<br>Associate Scientist<br>IUPUI<br>Physics Dept.<br>402 N. Blackford St.<br>Indianapolis, IN 46202-3273<div><br></div></div><br>
</body></html>