<html><head></head><body style="word-wrap: break-word; -webkit-nbsp-mode: space; -webkit-line-break: after-white-space; "><div>I want to add a term E to the MD potential energy function minimized by gromacs. The function E operates on the set of input atoms in a complicated way and returns a real value. I also have a closed-form representation of the first derivative of E with respect to each atom position. <div><br></div></div><div>I've managed to dig up a few questions related to mine in the mailing list archives. Specifically, </div><div><br></div><div><a href="http://lists.gromacs.org/pipermail/gmx-developers/2004-February/000767.html">http://lists.gromacs.org/pipermail/gmx-developers/2004-February/000767.html</a></div><div><br></div><div>says that this is as easy as modifying <font class="Apple-style-span" face="Courier">calc_f_el</font> in <font class="Apple-style-span" face="Courier">sim_util.c</font>. However, the very next message</div><div><br></div><div><a href="http://lists.gromacs.org/pipermail/gmx-developers/2004-February/000768.html">http://lists.gromacs.org/pipermail/gmx-developers/2004-February/000768.html</a></div><div><br></div><div>seems to contradict this, saying that <font class="Apple-style-span" face="Courier">force.c</font> has to be modified. Further, these messages are from 2004, and I haven't been able to locate anything more recent that is also more concrete. For instance, this message asks roughly the same question</div><div><br></div><div><a href="http://lists.gromacs.org/pipermail/gmx-developers/2009-March/003118.html">http://lists.gromacs.org/pipermail/gmx-developers/2009-March/003118.html</a></div><div><br></div><div><br></div><div>but the thread peters out after a few messages with no obvious conclusion/solution. </div><div><br></div><div>Anyway, I've explored the code a little bit, and I believe the following changes are necessary if I'm to add a term to the energy function. </div><div><br></div><div>1) The function most relevant seems to be <font class="Apple-style-span" face="Courier">do_force_lowlevel (force.c).</font> This function is called by <font class="Apple-style-span" face="Courier">do_force</font> in <font class="Apple-style-span" face="Courier">sim_util.c</font>, and appears to compute both the potential energy, stored in <font class="Apple-style-span" face="Courier">gmx_enerdata_t* enerd</font>, and the force, stored in <font class="Apple-style-span" face="Courier">rvec f[].</font> It seems that the right way to proceed would be to do the following: </div><div><br></div><div>a) Edit the <font class="Apple-style-span" face="Courier">gmx_enerdata_t</font> data structure so that it contains an entry for the new term E. This entails editing <font class="Apple-style-span" face="Courier">F_NRE</font>, the size of the array <font class="Apple-style-span" face="Courier">term</font> in that data structure (<font class="Apple-style-span" face="Courier">include/types/forcerec.h</font>). Presumably, since <font class="Apple-style-span" face="Courier">F_NRE </font>is the last entry in the enumeration data structure in <span class="Apple-style-span" style="font-family: Courier; ">include/types/idef.h, </span>adding an entry like<span class="Apple-style-span" style="font-family: Courier; "> F_NEW </span>to the enum data type just before<span class="Apple-style-span" style="font-family: Courier; "> F_NRE </span>will<font class="Apple-style-span" face="Courier"> </font>automagically make (int)F_NRE return the right number of interactions.<span class="Apple-style-span" style="font-family: Courier; "> </span>Correct?</div><div><font class="Apple-style-span" face="Courier"><br></font></div><div>b) Edit <font class="Apple-style-span" face="Courier">do_force_lowlevel</font> such that it adds the gradient due to my new energy function E to <font class="Apple-style-span" face="Courier">f[i]</font>, where <font class="Apple-style-span" face="Courier">i</font> is each atom. Also add the energy contribution due to E thus: <font class="Apple-style-span" face="Courier">E[F_NEW] = //funky new addition to energy function. </font></div><div><font class="Apple-style-span" face="Courier"><br></font></div><div>2) One worry is that the force <font class="Apple-style-span" face="Courier">f[i]</font>/position <font class="Apple-style-span" face="Courier">x[i]</font> doesn't really point to the force on/position of the i-th atom. I suspect this because of lines such as the following: </div><div><br></div><div><font class="Apple-style-span" face="Courier">//src/gmxlib/nonbonded.c</font></div><div><div><font class="Apple-style-span" face="Courier"> nb_kernel_allvsallgb_sse2_double(fr,mdatoms,excl,x[0],f[0],egcoul,egnb,egpol,</font></div><div><font class="Apple-style-span" face="Courier"> &outeriter,&inneriter,&fr->AllvsAll_work);</font></div></div><div><font class="Apple-style-span" face="Courier"><br></font></div><div>Why are <font class="Apple-style-span" face="Courier">f[0]</font> and <font class="Apple-style-span" face="Courier">x[0]</font> special in the above function call? </div><div><br></div><div>3) Finally, I need to edit the input processing subroutines in such a way that I can turn my interaction on and off at will. I'm not sure how exactly to do this, but I suspect that the function <font class="Apple-style-span" face="Courier">parse_common_args</font> in <font class="Apple-style-span" face="Courier">src/gmxlib/statutil.c</font> is relevant. I'd appreciate any "input" on this (pun unintended). </div><div><br></div><div>4) I'd like to know if there's something I haven't covered in the above comments/ questions. </div><div><br></div><div>Also, it'd be great if someone reading this were a) at the University of Texas at Austin, where I'm a graduate student, and b) willing to spare some of their time, not more than a half-hour, to go over the code with me. I promise I'll make it worth their while; we could get a coffee at JP's Java! </div><div><br></div><div><br></div><div>thanks for your time,</div><br><div apple-content-edited="true">
<div>Radhakrishna Bettadapura</div><div>PhD candidate,</div><div>Dept. of Mechanical Engineering/Computational Visualization Center,</div><div>University of Texas at Austin</div><div><a href="http://users.ices.utexas.edu/~brrk">http://users.ices.utexas.edu/~brrk</a></div></div></body></html>