<div dir="ltr">Dear Vladimir,<div>I just came across the same problem as you.</div><div>First, when you calculate the variance of the totoal moment, you need to use the norm of dipole vector, rather than directly calculate the average number like this: <span style="font-family: monospace; background-color: rgb(250, 250, 250); color: rgb(0, 0, 0);">np</span><span style="font-family: monospace; background-color: rgb(250, 250, 250); color: rgb(102, 102, 0);">.</span><span style="font-family: monospace; background-color: rgb(250, 250, 250); color: rgb(0, 0, 0);">aver<wbr>age</span><span style="font-family: monospace; background-color: rgb(250, 250, 250); color: rgb(102, 102, 0);">(</span><span style="font-family: monospace; background-color: rgb(250, 250, 250); color: rgb(0, 0, 0);">totalMoment</span><span style="font-family: monospace; background-color: rgb(250, 250, 250); color: rgb(102, 102, 0);">).</span></div><div><span style="font-family: monospace; background-color: rgb(250, 250, 250); color: rgb(102, 102, 0);"><br></span></div><div>Second, I recalculate the dipole by myself using sum(coordination X point charge), it gives me a different dipole number compared to the CP2K output dipole moment, which finally gives me a correct dielectric constant, I don't know why there's such a difference as I can not find the way of CP2K calculating dipole moment.<span style="font-family: monospace; background-color: rgb(250, 250, 250); color: rgb(102, 102, 0);"><br></span></div><div><br></div><div>I hope it also works in your case.</div><div><br></div><div>Best,</div><div><br></div><div>Chris</div><div><br>在 2019年11月15日星期五 UTC+8下午4:01:10,Nt写道:<blockquote class="gmail_quote" style="margin: 0;margin-left: 0.8ex;border-left: 1px #ccc solid;padding-left: 1ex;"><div dir="ltr">Dear all:<div><br></div><div>I hope it's OK to repost this question, which didn't get any replies when first posted. I'm hoping that someone who has experience with using CP2K to calculate the dielectric constant will see it because I'm stuck.</div><div><br></div><div>At a temperature of 20 C and atmospheric pressure, the dielectric constant of water has been measured at approximately 80. However, when I calculate it based on a FIST SPC/E simulation and a CP2K tutorial (link provided below), I get a result of <1. I am wondering if maybe there is a problem with my understanding of the units from CP2K or the way I outputted the dipole moment. This is pretty new to me, so I could very well have made a simple mistake somewhere.</div><div><br></div><div>1. I added the following in:</div><div><br><div style="background-color:rgb(250,250,250);border-color:rgb(187,187,187);border-style:solid;border-width:1px"><code><div><span style="color:#660">&</span><span style="color:#606">Force_Eval</span><span style="color:#660">/&</span><span style="color:#000">MM</span><span style="color:#660">/&</span><span style="color:#000">PRINT</span></div></code></div><br></div><div style="background-color:rgb(250,250,250);border-color:rgb(187,187,187);border-style:solid;border-width:1px"><code><div><span style="color:#660">&</span><span style="color:#000">DIPOLE<br> PERIODIC TRUE<br> </span><span style="color:#660">&</span><span style="color:#000">EACH<br> MD </span><span style="color:#066">1</span><span style="color:#000"><br> </span><span style="color:#660">&</span><span style="color:#008">END</span><span style="color:#000"><br> FILENAME </span><span style="color:#660">=</span><span style="color:#000"> moments</span><span style="color:#660">.</span><span style="color:#000">traj<br> ADD_LAST NUMERIC<br></span><span style="color:#660">&</span><span style="color:#008">END</span><span style="color:#000"> DIPOLE</span></div></code></div><div><br></div><div>2. The output file "moments.traj" contains lines like these at each time step with the components of the dipole moment:</div><div><div style="background-color:rgb(250,250,250);border-color:rgb(187,187,187);border-style:solid;border-width:1px"><code><div><span style="color:#000"> MM DIPOLE </span><span style="color:#660">[</span><span style="color:#000">BERRY PHASE</span><span style="color:#660">](</span><span style="color:#000">A</span><span style="color:#660">.</span><span style="color:#000">U</span><span style="color:#660">.)|</span><span style="color:#000"> </span><span style="color:#660">-</span><span style="color:#066">3.333675</span><span style="color:#000"> </span><span style="color:#066">0.611131</span><span style="color:#000"> </span><span style="color:#066">1.550260</span><span style="color:#000"><br> MM DIPOLE </span><span style="color:#660">[</span><span style="color:#000">BERRY PHASE</span><span style="color:#660">](</span><span style="color:#606">Debye</span><span style="color:#660">)|</span><span style="color:#000"> </span><span style="color:#660">-</span><span style="color:#066">8.473356</span><span style="color:#000"> </span><span style="color:#066">1.553341</span><span style="color:#000"> </span><span style="color:#066">3.940367</span><span style="color:#000"><br> MM DIPOLE </span><span style="color:#660">[</span><span style="color:#000">BERRY PHASE</span><span style="color:#660">]</span><span style="color:#000"> DERIVATIVE</span><span style="color:#660">(</span><span style="color:#000">A</span><span style="color:#660">.</span><span style="color:#000">U</span><span style="color:#660">.)|</span><span style="color:#000"> </span><span style="color:#066">0.000767</span><span style="color:#000"> </span><span style="color:#660">-</span><span style="color:#066">0.000211</span><span style="color:#000"> </span><span style="color:#066">0.001268</span></div></code></div><div><br></div></div><div>I extracted the lines in [Debye] units and then turned to a tutorial at cp2k.org (link below), where I followed the procedure on how to calculate the dielectric constant. </div><div><br></div><div><div><span style="color:rgb(0,0,255)"><a href="https://www.cp2k.org/exercises:2014_ethz_mmm:monte_carlo_ice#task_1calculate_the_dielectric_constant" rel="nofollow" target="_blank" onmousedown="this.href='https://www.google.com/url?q\x3dhttps%3A%2F%2Fwww.cp2k.org%2Fexercises%3A2014_ethz_mmm%3Amonte_carlo_ice%23task_1calculate_the_dielectric_constant\x26sa\x3dD\x26sntz\x3d1\x26usg\x3dAFQjCNG9OsRLQiGdf5qrvpsTVJsTk7C19Q';return true;" onclick="this.href='https://www.google.com/url?q\x3dhttps%3A%2F%2Fwww.cp2k.org%2Fexercises%3A2014_ethz_mmm%3Amonte_carlo_ice%23task_1calculate_the_dielectric_constant\x26sa\x3dD\x26sntz\x3d1\x26usg\x3dAFQjCNG9OsRLQiGdf5qrvpsTVJsTk7C19Q';return true;">https://www.cp2k.org/<wbr>exercises:2014_ethz_mmm:monte_<wbr>carlo_ice#task_1calculate_the_<wbr>dielectric_constant</a></span><br></div><div>("Task 1: Calculate the dielectric constant") </div></div><div><br></div><div>- First I used the tutorial to calculate a total dipole moment from the three columns in the output snippet above (sum of squares of each of the three components, followed by taking the square root):</div><div><div><br></div><div><pre style="padding:0.7em 1em;border-width:1px;border-style:solid;border-color:rgb(204,204,204);font-family:Consolas,"Andale Mono WT","Andale Mono","Bitstream Vera Sans Mono","Nimbus Mono L",Monaco,"Courier New",monospace;font-size:14px;direction:ltr;background-color:rgb(251,250,249);color:rgb(51,51,51);border-radius:2px;overflow:auto">M <span style="color:rgb(102,204,102)">=</span> np.<span style="color:rgb(0,102,0)">sqrt</span><span style="color:rgb(102,204,102)">(</span>np.<span style="color:rgb(0,0,0);font-weight:bold">sum</span><span style="color:rgb(102,204,102)">(</span>np.<span style="color:rgb(0,102,0)">square</span><span style="color:rgb(102,204,102)">(</span>M_vec<span style="color:rgb(102,204,102)"><wbr>)</span><span style="color:rgb(102,204,102)">,</span> axis<span style="color:rgb(102,204,102)">=</span><span style="color:rgb(204,102,204)">1</span><span style="color:rgb(102,204,102)">)</span><span style="color:rgb(102,204,102)">)</span> <span style="color:rgb(128,128,128);font-style:italic"># Take norm of dipole moment (from link above)</span>
</pre></div></div><div><br></div><div>- Then I converted the total moment in Debye units to SI units [C*m] by multiplying the total moment in Debye by <span style="color:rgb(39,42,216);font-family:Menlo;font-size:12px">3.335641e-30</span>.<br></div><div><br></div><div>3. </div><div><div>- I executed the following line from the link, where I had converted everything to SI units: </div><div><br></div><div><pre style="padding:0.7em 1em;border-width:1px;border-style:solid;border-color:rgb(204,204,204);font-family:Consolas,"Andale Mono WT","Andale Mono","Bitstream Vera Sans Mono","Nimbus Mono L",Monaco,"Courier New",monospace;font-size:14px;direction:ltr;background-color:rgb(251,250,249);color:rgb(51,51,51);border-radius:2px;overflow:auto">s <span style="color:rgb(102,204,102)">=</span> <span style="color:rgb(102,204,102)">(</span>e*e*<span style="color:rgb(204,102,204)">4</span>*np.<span style="color:rgb(0,102,0)">pi</span><span style="color:rgb(102,204,102)">)</span>/<span style="color:rgb(102,204,102)">(</span><span style="color:rgb(204,102,204)">3</span>*V*kb*T*<wbr>angstrom2meter*epsilon_0<span style="color:rgb(102,204,102)">)</span>
</pre></div><div><br></div><div>When I multiplied the result with the total dipole moment I calculated above, the answer was far from 80. I used another version of the equation that is more common in published literature without the epsilon_0 term in the denominator and without the "e" terms. It gives me an answer much less than 1.<br></div><div><br></div><div>Here is how I calculated the variance of the total dipole moment:</div><div><div style="border-width:1px;border-style:solid;border-color:rgb(187,187,187);background-color:rgb(250,250,250)"><code><span style="color:rgb(0,0,136)">var</span><span style="color:rgb(0,0,0)"> </span><span style="color:rgb(102,102,0)">=</span><span style="color:rgb(0,0,0)"> np</span><span style="color:rgb(102,102,0)">.</span><span style="color:rgb(0,0,0)">average</span><span style="color:rgb(102,102,0)">(</span><span style="color:rgb(0,0,0)">np</span><span style="color:rgb(102,102,0)">.</span><span style="color:rgb(0,0,0)">square</span><span style="color:rgb(102,102,0)">(</span><span style="color:rgb(0,0,0)">tot<wbr>alMoment</span><span style="color:rgb(102,102,0)">))</span><span style="color:rgb(0,0,0)"> </span><span style="color:rgb(102,102,0)">-</span><span style="color:rgb(0,0,0)"> np</span><span style="color:rgb(102,102,0)">.</span><span style="color:rgb(0,0,0)">square</span><span style="color:rgb(102,102,0)">(</span><span style="color:rgb(0,0,0)">np</span><span style="color:rgb(102,102,0)">.</span><span style="color:rgb(0,0,0)">aver<wbr>age</span><span style="color:rgb(102,102,0)">(</span><span style="color:rgb(0,0,0)">totalMoment</span><span style="color:rgb(102,102,0)">))</span><span style="color:rgb(0,0,0)"><br></span></code></div><br>I also tried Python's variance and it seemed to generate a similar result:</div><div><br></div><div><div style="border-width:1px;border-style:solid;border-color:rgb(187,187,187);background-color:rgb(250,250,250)"><code><span style="color:rgb(0,0,136)">var</span><span style="color:rgb(0,0,0)"> </span><span style="color:rgb(102,102,0)">=</span><span style="color:rgb(0,0,0)"> </span><span style="color:rgb(102,102,0)">[</span><span style="color:rgb(0,0,0)">np</span><span style="color:rgb(102,102,0)">.</span><span style="color:rgb(0,0,136)">var</span><span style="color:rgb(102,102,0)">(</span><span style="color:rgb(0,0,0)">totalMoment</span><span style="color:rgb(102,102,0)">[:</span><span style="color:rgb(0,0,0)">i</span><span style="color:rgb(102,102,0)">+</span><span style="color:rgb(0,102,102)"><wbr>1</span><span style="color:rgb(102,102,0)">])</span><span style="color:rgb(0,0,0)"> </span><span style="color:rgb(0,0,136)">for</span><span style="color:rgb(0,0,0)"> i </span><span style="color:rgb(0,0,136)">in</span><span style="color:rgb(0,0,0)"> range</span><span style="color:rgb(102,102,0)">(</span><span style="color:rgb(0,0,0)">len</span><span style="color:rgb(102,102,0)">(</span><span style="color:rgb(0,0,0)">totalMo<wbr>ment</span><span style="color:rgb(102,102,0)">))]</span><span style="color:rgb(0,0,0)"></span></code></div><br>If anyone has any ideas on why I am not having any success, I would be very grateful because I am stuck.</div></div></div></blockquote></div></div>