<html xmlns:o="urn:schemas-microsoft-com:office:office" xmlns:w="urn:schemas-microsoft-com:office:word" xmlns:m="http://schemas.microsoft.com/office/2004/12/omml" xmlns="http://www.w3.org/TR/REC-html40"><head><meta http-equiv=Content-Type content="text/html; charset=utf-8"><meta name=Generator content="Microsoft Word 15 (filtered medium)"><style><!--
/* Font Definitions */
@font-face
        {font-family:"Cambria Math";
        panose-1:2 4 5 3 5 4 6 3 2 4;}
@font-face
        {font-family:DengXian;
        panose-1:2 1 6 0 3 1 1 1 1 1;}
@font-face
        {font-family:Calibri;
        panose-1:2 15 5 2 2 2 4 3 2 4;}
@font-face
        {font-family:"\@DengXian";
        panose-1:2 1 6 0 3 1 1 1 1 1;}
@font-face
        {font-family:"Times New Roman \,serif";}
/* Style Definitions */
p.MsoNormal, li.MsoNormal, div.MsoNormal
        {margin:0in;
        margin-bottom:.0001pt;
        font-size:11.0pt;
        font-family:"Calibri",sans-serif;}
a:link, span.MsoHyperlink
        {mso-style-priority:99;
        color:blue;
        text-decoration:underline;}
a:visited, span.MsoHyperlinkFollowed
        {mso-style-priority:99;
        color:#954F72;
        text-decoration:underline;}
.MsoChpDefault
        {mso-style-type:export-only;}
@page WordSection1
        {size:8.5in 11.0in;
        margin:1.0in 1.25in 1.0in 1.25in;}
div.WordSection1
        {page:WordSection1;}
--></style></head><body lang=EN-US link=blue vlink="#954F72"><div class=WordSection1><p class=MsoNormal><span style='font-size:12.0pt'>I think it should be:<o:p></o:p></span></p><p class=MsoNormal><b><span style='font-size:12.0pt'>            </span></b><span style='font-size:12.0pt'>while (i &lt; top-&gt;excls.nr &amp;&amp; <b>n &lt;= ntarget</b>)<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>in forcerec_set_excl_load() function (in forcerec.c file).<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>It should work then.<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>Please correct me if I am wrong.<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>Thanks!<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>Sheng<o:p></o:p></span></p><p class=MsoNormal><o:p>&nbsp;</o:p></p><p class=MsoNormal>--<br>PhD student<br>State Key Laboratory of Coal Combustion,<br>School of Energy and Power Engineering,<br>Huazhong University of Science and Technology (HUST),<br>Room 302 Power Building, <br>1037 Luoyu Road, Wuhan, Hubei 430074 China<br>Email: chrishengbee@hust.edu.cn<br>Phone: (86)27-87548122 <br>http://itp.energy.hust.edu.cn </p><p class=MsoNormal><span style='font-size:12.0pt;font-family:"Times New Roman",serif'><o:p>&nbsp;</o:p></span></p><div style='mso-element:para-border-div;border:none;border-top:solid #E1E1E1 1.0pt;padding:3.0pt 0in 0in 0in'><p class=MsoNormal style='border:none;padding:0in'><b>From: </b><a href="mailto:chrishengbee@hust.edu.cn">Sheng Bi</a><br><b>Sent: </b>Tuesday, June 7, 2016 11:04 PM<br><b>To: </b><a href="mailto:hess@kth.se">Berk Hess</a>; <a href="mailto:gmx-developers@gromacs.org">gmx-developers@gromacs.org</a><br><b>Cc: </b><a href="mailto:gfeng@hust.edu.cn">gfeng@hust.edu.cn</a><br><b>Subject: </b>Re: [gmx-developers] Question about ewald_LRcorrection() function</p></div><p class=MsoNormal><span style='font-size:12.0pt;font-family:"Times New Roman",serif'><o:p>&nbsp;</o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>Hi Berk,<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; Thanks for your response!<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; However, my question is still there. In fact, I run simulation in one core (with multi-threadings), so domain decomposition could not be the reason for this phenomenon. Meanwhile, we use a quite large box (due to 3DC, Z length usually need a huge vacuum), so I think the water is within the box.<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; Actually, I have printed dipcorrA for both gmx331 and gmx467 and find that there is no difference. So it wouldn’t be the PBC trouble.<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; I think <b>you are right about ewald_LRcorrection().</b> I check this function further and think it works fine. <o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; However, I think there is <b>something wrong in forcerec_set_excl_load() function (in forcerec.c file).</b> We know in ewald_LRcorrection(), for each thread, we loop over (i = fr-&gt;excl_load[t]; i &lt; fr-&gt;excl_load[t+1]) &nbsp;to do dipole correction. The excl_load array is calculated in forcerec_set_excl_load(), as shown in following:<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>/ ****** line 3018 ********* /<o:p></o:p></span></p><div><div><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp; for (t = 1; t &lt;= fr-&gt;nthreads; t++)<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp; {<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; ntarget = (ntot*t)/fr-&gt;nthreads;<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; <b>while (i &lt; top-&gt;excls.nr &amp;&amp; n &lt; ntarget)<o:p></o:p></b></span></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; {<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; for (j = ind[i]; j &lt; ind[i+1]; j++)<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; {<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; if (a[j] &gt; i)<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; {<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; n++;<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; }<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; }<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; i++;<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; }<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; fr-&gt;excl_load[t] = i;<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp; }<o:p></o:p></span></p></div></div><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; <o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; Suppose we have a system consisted by <b>2 types of molecules</b>. The whole system contains 2 <b>atoms</b>, 1 for molecule <b>A(positively charged)</b> and 1 for molecule <b>B(negatively charged)</b>. We perform simulation in one core with one thread. So, for code above, we can deduce parameters like following:<o:p></o:p></span></p><p class=MsoNormal><b><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; top-&gt;excls.nr = 2<o:p></o:p></span></b></p><p class=MsoNormal><b><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; fr-&gt;nthreads = 1 <o:p></o:p></span></b></p><p class=MsoNormal><b><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; ntarget = ntot =1<o:p></o:p></span></b></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; Since it <b>use while() statement here</b>, we can deduce that:<o:p></o:p></span></p><p class=MsoNormal><b><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; fr-&gt;excl_load[0] =0<o:p></o:p></span></b></p><p class=MsoNormal><b><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; fr-&gt;excl_load[1] =1<o:p></o:p></span></b></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; However, <b>fr-&gt;excl_load[1] should be 2</b>. This wrong value lead to <b>missing a dipole correction on the 2<sup>th</sup> atom</b> in ewald_LRcorrection(). Because in ewald_LRcorrection(), we loop from fr-&gt;excl_load[0]=0 to fr-&gt;excl_load[1] =1 using for statement.<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; I think this is the direct reason for wrong additional force on water molecule in my system. This error might result big trouble when someone using umbrella sampling in a system with system dipole.<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; Please correct me if I am wrong.<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>Best regards!<o:p></o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt'>Sheng <o:p></o:p></span></p><p class=MsoNormal>--<br>PhD student<br>State Key Laboratory of Coal Combustion,<br>School of Energy and Power Engineering,<br>Huazhong University of Science and Technology (HUST),<br>Room 302 Power Building, <br>1037 Luoyu Road, Wuhan, Hubei 430074 China<br>Email: chrishengbee@hust.edu.cn<br>Phone: (86)27-87548122 <br>http://itp.energy.hust.edu.cn </p><p class=MsoNormal><span style='font-size:12.0pt;font-family:"Times New Roman",serif'><o:p>&nbsp;</o:p></span></p><div style='border:none;border-top:solid #E1E1E1 1.0pt;padding:3.0pt 0in 0in 0in'><p class=MsoNormal><b>From: </b><a href="mailto:hess@kth.se">Berk Hess</a><br><b>Sent: </b>Tuesday, June 7, 2016 5:49 PM<br><b>To: </b><a href="mailto:gmx-developers@gromacs.org">gmx-developers@gromacs.org</a><br><b>Subject: </b>Re: [gmx-developers] Question about ewald_LRcorrection() function<o:p></o:p></p></div><p class=MsoNormal><span style='font-size:12.0pt;font-family:"Times New Roman",serif'><o:p>&nbsp;</o:p></span></p><p class=MsoNormal><span style='font-size:12.0pt;font-family:"Times New Roman",serif'>Hi,<br><br>I think you misinterpreted the code in ewald_LRcorrection. It does the same things in 3.3.1 and 4.6.7.<br>My guess is that the issue you are observing is due to periodic boundary conditions. With group cut-off scheme charge groups are whole (not split by PBC), so water molecules are usually whole. Then the 3DC correction always works correctly when you have charge groups with zero net charge. In the Verlet cut-off scheme with domain decomposition, or with charged charge groups, (partial) charges moving over PBC change the system dipole. To get your system to run correctly you need to ensure that no water molecule moves over the periodic boundary in along the z-dimension. I think grompp will print a warning for this issue.<br><br>Cheers,<br><br>Berk<br><br>On 2016-06-07 10:12, Sheng Bi wrote:<o:p></o:p></span></p><blockquote style='margin-top:5.0pt;margin-bottom:5.0pt'><p class=MsoNormal>Sorry, some type-offs in previous question.<o:p></o:p></p><p class=MsoNormal>In gmx-4.6.7, ewald_LRcorrection()<b>.</b> We loop over excluded neighbours only in condition <b>calc_excl_corr=1</b>. But when using Verlet cut-off scheme, <b>calc_excl_corr will be 0</b>.</p><p class=MsoNormal>&nbsp;</p><p class=MsoNormal>&nbsp;</p><p class=MsoNormal>--<br>PhD student<br>State Key Laboratory of Coal Combustion,<br>School of Energy and Power Engineering,<br>Huazhong University of Science and Technology (HUST),<br>Room 302 Power Building, <br>1037 Luoyu Road, Wuhan, Hubei 430074 China<br>Email: <a href="mailto:chrishengbee@hust.edu.cn">chrishengbee@hust.edu.cn</a><br>Phone: (86)27-87548122 <br><a href="http://itp.energy.hust.edu.cn">http://itp.energy.hust.edu.cn</a> </p><p class=MsoNormal><span style='font-size:12.0pt;font-family:"Times New Roman \,serif"'>&nbsp;</span></p><div style='border:none;border-top:solid #E1E1E1 1.0pt;padding:3.0pt 0in 0in 0in'><p class=MsoNormal><b>From: </b><a href="mailto:chrishengbee@hust.edu.cn">Sheng Bi</a><br><b>Sent: </b>Tuesday, June 7, 2016 1:28 AM<br><b>To: </b><a href="mailto:gmx-developers@gromacs.org">gmx-developers@gromacs.org</a><br><b>Subject: </b>[gmx-developers] Question about ewald_LRcorrection() function</p></div><p class=MsoNormal><span style='font-size:12.0pt;font-family:"Times New Roman \,serif"'>&nbsp;</span></p><p class=MsoNormal>Dear GMX developers:</p><p class=MsoNormal>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; I find a abnormal phenomenon when using <b>gromacs-4.6.7</b>. I set up a system containing two walls(consist of freeze atoms), one wall has positively charges on atoms, another has negatively charges.&nbsp; I put an water molecule(rigid spc/spce/tip3p model) between these two walls. <b>Whole system seems like below</b>:</p><p class=MsoNormal>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; ****************************************</p><p class=MsoNormal>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; * - &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; H&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; +&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; *&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; </p><p class=MsoNormal>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; * -&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; &nbsp;&nbsp; \&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; +&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; vacuum&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; *</p><p class=MsoNormal>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; * -&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; &nbsp;&nbsp; &nbsp;&nbsp;O&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; +&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; *</p><p class=MsoNormal>&nbsp; &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; * -&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; &nbsp; &nbsp;/&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; +&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; *</p><p class=MsoNormal>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; * -&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; H&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; &nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; +&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; *</p><p class=MsoNormal>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; ****************************************</p><p class=MsoNormal>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; I use <b>Verlet cut-off scheme, PME method, and ewald-geometry=3DC</b>. Temperature of whole system is at <b>0K</b>, controlled by <b>v-rescale method</b>. I run whole simulation <b>only in OPENMP parallel</b>.</p><p class=MsoNormal>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; One can have following expected behavior of water molecule without MD simulation: Since it can be seen as an infinite panel capacitor, <b>there is constant electric field between two walls. So the water will change it orientation, but won’t change it’s position(Since 0K)</b>.</p><p class=MsoNormal>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; This phenomenon can be <b>confirmed by gromcas-3.3.1</b> (although there is position shift, but quite small). Also, if one use higher version of gromacs, such as 4.6.7 or higher <b>with ewald-geometry=3D, either no problem</b>. However, if one use <b>3DC</b> (and it should be), you can see <b>an obvious force on water molecule</b>.</p><p class=MsoNormal>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; I hacked into source code, and find this might due to <b>ewald_LRcorrection()</b> function:</p><p class=MsoNormal>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; Firstly, in gmx-4.6.7, ewald_LRcorrection()<b>.</b> We <b>loop over excluded neighbours only in condition calc_excl_corr=0. But when using verlet cut-off scheme, calc_excl_corr will be 1</b>. In gmx-3.3.1, we do loop over excluded neighbours without additional if() statement.</p><p class=MsoNormal>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; Secondly, in gmx-4.6.7, , ewald_LRcorrection(). When do dipole correction on force, we loop over <b>excl_load</b> (Exclusion load distribution over the threads). In fact, most water model will have nrexel=2. So, as a result, we <b>only do dipole correction on O atom, but do nothing on H atoms</b> and we will see an obvious force on the whole water molecule. I think this is not purpose of 3DC<span lang=ZH-CN style='font-family:DengXian'>。</span></p><p class=MsoNormal>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; Best regards<span lang=ZH-CN style='font-family:DengXian'>!</span></p><p class=MsoNormal>Sheng Bi</p><p class=MsoNormal><b>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; </b></p><p class=MsoNormal>---<br>PhD student<br>State Key Laboratory of Coal Combustion,<br>School of Energy and Power Engineering,<br>Huazhong University of Science and Technology (HUST),<br>Room 302 Power Building, <br>1037 Luoyu Road, Wuhan, Hubei 430074 China<br>Email: <a href="mailto:chrishengbee@hust.edu.cn">chrishengbee@hust.edu.cn</a><br>Phone: (86)27-87548122 <br><a href="http://itp.energy.hust.edu.cn">http://itp.energy.hust.edu.cn</a> </p><p class=MsoNormal><span style='font-size:12.0pt;font-family:"Times New Roman \,serif"'>&nbsp;</span></p><p class=MsoNormal>&nbsp;</p><p class=MsoNormal style='margin-bottom:12.0pt'><span style='font-family:"Times New Roman",serif'><o:p>&nbsp;</o:p></span></p></blockquote><p class=MsoNormal><span style='font-family:"Times New Roman",serif'><o:p>&nbsp;</o:p></span></p><p class=MsoNormal><span style='color:black'><o:p>&nbsp;</o:p></span></p><p class=MsoNormal><o:p>&nbsp;</o:p></p></div></body></html>