GROMACS development: Issueshttp://redmine.gromacs.org/http://redmine.gromacs.org/favicon.ico?14024030812006-02-18T16:48:54ZGROMACS development
Redmine GROMACS - Bug #49 (Closed): Wrong image convention for evaluating dummy atoms.http://redmine.gromacs.org/issues/492006-02-18T16:48:54ZRamon Garciaramon@juguete.quim.ucm.es
<p>After running mdrun_d with a water system with van der Eerden model, mdrun_d<br />crashed. (Use attached files to reproduce the problem). After debugging, I<br />discovered that the reasons was that the coordinates of a dummy atom were wrong.<br />When they were evaluated, the positions of the atoms have been shifted twice<br />from the original ones. As the coordinates of dummy atoms are a linear<br />combination of other atoms, if the image conventions of them are inconsistent,<br />the dummy can have very wrong values. This caused a crash in the calculation of<br />forces.</p>
<p>Looking at the source code, it looks like the logic of shifting before virtual<br />site evaluation is wrong</p>
<pre><code>if (bVsites) {<br /> if (graph) {<br /> /* Following is necessary because the graph may get out of sync
* with the coordinates if we only have every N'th coordinate set<br /> */<br /> if (bRerunMD || bExchanged)<br /> mk_mshift(log,graph,state->box,state->x);<br /> shift_self(graph,state->box,state->x);<br /> }<br /> construct_vsites(log,state->x,&mynrnb,inputrec->delta_t,state->v,<br /> &top->idef,graph,cr,fr->ePBC,state->box,vsitecomm);</code></pre>
<pre><code>if (graph)<br /> unshift_self(graph,state->box,state->x);<br /> }</code></pre>
<p>This code is written as if the particles were not shifted before, and so they<br />are shifted and then unshifted to restore the previous state. But previously<br />there has been a call to do_pbc_first and therefore the positions have been<br />already shifted. I changed the code above to:</p>
<pre><code>if (bVsites) {<br /> if (graph) {<br /> /* Following is necessary because the graph may get out of sync
* with the coordinates if we only have every N'th coordinate set<br /> */<br /> if (bRerunMD || bExchanged) {<br /> unshift_self(graph, state->box, state->x);<br /> mk_mshift(log,graph,state->box,state->x);<br /> shift_self(graph,state->box,state->x);<br /> }<br /> }<br /> construct_vsites(log,state->x,&mynrnb,inputrec->delta_t,state->v,<br /> &top->idef,graph,cr,fr->ePBC,state->box,vsitecomm);</code></pre>
<pre><code>}</code></pre>
<p>which seems to be correct for me. Please review the proposed change. I have only<br />tested with the system I am working with.</p> GROMACS - Bug #39 (Closed): Generation of initial velocities counts dummy atoms as normal atoms.http://redmine.gromacs.org/issues/392005-12-30T17:09:54ZRamon Garciaramon@juguete.quim.ucm.es
<p>When grompp generates initial velocities, it counts all the atoms for the<br />evaluation of the total kinetic energy. Dummy atoms should not be counted.<br />Because of this, simulations containing water models with dummy atoms such as<br />TIP4P, start with a temperature higher than the expected one. In particular, in<br />a simulation of water with the van der Eerden model (which has the same number<br />of dummy atoms and real atoms), the initial temperature is exactly twice the one<br />the user specified in grompp.mdp .</p> GROMACS - Bug #35 (Closed): trjconv does not read velocities from g96 file.http://redmine.gromacs.org/issues/352005-12-06T18:45:10ZRamon Garciaramon@juguete.quim.ucm.es
<p>When using trjconv to convert a .g96 file to a .trr file, the output file does<br />not contain velocities.</p>
<p>The cause of this problem is a typo in the source src/tools/gmx_trjconv.c</p>
<p>The source reads (line 700 in Gromacs 3.3):</p>
<pre><code>bVels= (ftp==efTRR || ftp==efTRJ || ftp==efGRO || ftp==efG96) <br /> && (ftpin==efTRR || ftpin==efTRJ || ftpin==efGRO || ftp==efG96);</code></pre>
<p>Note the mistake at the end of the second line. The source should read:</p>
<pre><code>bVels= (ftp==efTRR || ftp==efTRJ || ftp==efGRO || ftp==efG96) <br /> && (ftpin==efTRR || ftpin==efTRJ || ftpin==efGRO || ftpin==efG96);</code></pre> GROMACS - Bug #18 (Closed): Implement position rescaling with Rahman Parrinello barostat.http://redmine.gromacs.org/issues/182005-10-10T15:27:48ZRamon Garciaramon@juguete.quim.ucm.es
<p>Althougth there is some ongoing work to enhance the barostat with Gromacs, the<br />current barostat can be made more accurate with a painless change. In<br />do_update_md, modify the updating of the coordinates to include the change of<br />the box. In fact, I believe that not doing so was a misunderstanding of the<br />equations. The differential equation of the coordinates, in the box reference<br />frame is:</p>
<p>m_i d^2s_i/dt^2 = h^(-1) f_i - m_i G^(-1) dG/dt d s_i/dt</p>
<p>s = coordinates in the box frame<br />h = box matrix, transposed of Gromacs box<br />G = h^t t</p>
<p>Looking at Gromacs code in do_update_md and parrinellorahman_pcoupl and<br />comparing to the equation above we see that the velocity is defined as:</p>
<p>v_i = h ds_i/dt</p>
<p>therefore, the updating of the coordinates should be</p>
<p>dx_i/dt = d(h s_i)/dt = dh/dt * s_i + h * d s_i/dt =
= dh/dt * h^(-1) * x_i + v_i =
= M * x_i + v_i</p>
<p>(the last step assumes that dh/dt * h^(-1) is symmetric)</p> GROMACS - Bug #17 (Closed): Double precision gromacs truncates in memory coordinates to single pr...http://redmine.gromacs.org/issues/172005-10-06T11:44:49ZRamon Garciaramon@juguete.quim.ucm.es
<p>When mdrun_d writes an xtc file, in addition to the obvious truncation of the<br />data stored in disk, Gromacs also truncates the coordinates stored in the<br />simulation state.</p>
<p>I noticed because I am modifying Gromacs implementation of Rahman Parrinello. I<br />noticed some unexpected changes in state->x<sup><a href="#fn0">0</a></sup> and so placed a watchpoint in the<br />debugger, which showed that the routine xdr3drcoord in xdrd.c truncates not only<br />the coordinates that will be stored in disk, but also the original copy of the<br />simulation state.</p>
<p>I deleted the code that copies the truncated coordinates to the original copy<br />received by the routine. That is, I removed the lines</p>
<pre><code>for(i=0; (i&lt;isize); i++)<br /> fp[i]=ffp[i];</code></pre>
<p>But I am not sure if this change is correct. Is this routine supposed to return<br />the modified coordinates to report the caller? In that case the caller should<br />have made a copy.</p>
<p>Note that this problem can affect the correctness of simulations in double<br />precision.</p> GROMACS - Bug #16 (Closed): Issue an error if pme_order < fourier_nx/y/zhttp://redmine.gromacs.org/issues/162005-10-03T19:21:50ZRamon Garciaramon@juguete.quim.ucm.es
<p>I found by a random test that if pme_order is less than the number of cells in a<br />dimension of the Fourier grid Gromacs crashes. This is probably an invalid<br />configuration (I know almost nothing about PME), but it is much better to issue<br />an error in grompp than crashing.</p>
<p>Here is the code</p>
<p>in pme.c function make_bspline_moduli</p>
<p>650 for(i=0;i<nmax;i++)<br />651 bsp_data[i]=0;<br />652 for(i=1;i<=order;i++)<br />653 bsp_data[i]=data[i-1];</p>
<p>when order < nmax, the loop below access the array bsp_data beyond its bounds.</p> GROMACS - Bug #14 (Closed): Rahman Parrinello barostat not accurate for the anisotropic case.http://redmine.gromacs.org/issues/142005-09-15T18:03:01ZRamon Garciaramon@juguete.quim.ucm.es
<p>The implementation of the Rahman Parrinello barostat does not match the one<br />published in M. Parrinello, A. Rahaman J. Appl. Phys vol 52 n 12 page 7182</p>
<p>For the case where only isotropic pressure is applied, the implementation is<br />correct. However, when there is an non isotropic stress tensor, Gromacs uses an<br />incorrect generalization that consists in replacing the pressure by the stress<br />tensor in the formula:</p>
<p>W box'' = (pressure_observed - pressure) box^{-1} * volume</p>
<p>This seems a natural generalization. In fact first paper anticipating the Rahman<br />Parrinello thermostat (PRL vol 45, n 14, pg 1196) suggested this extension. But<br />the paper (earlier mentioned) that developed the method gives a different formula.</p>
<p>A posible solution is to support isotropic pressure only, even when allowing box<br />deformations. This is done by DL_POLY. In fact, the review article of Nose and<br />Klein (Molecular Physics, vol 50, n 5 pg 1055) cited by Gromacs does not mention<br />how to implement a barostat to a stress tensor. Furthermore according to the<br />paper by Rahman and Parrinello, if a external stress tensor is applied, both the<br />hydrostatic pressure and the tensor have to be specified. So ref_p for the<br />anisotropic case would need 7 components, one for the isotropic pressure and 6<br />for the external stress tensor.</p> GROMACS - Bug #12 (Closed): Forces between frozen particles contribute to pressure.http://redmine.gromacs.org/issues/122005-09-12T15:49:03ZRamon Garciaramon@juguete.quim.ucm.es
<p>If two particles are frozen, their force should not contribute to pressure. It<br />does currently in Gromacs.</p> GROMACS - Bug #11 (Closed): ngmx exit confirmation is brokenhttp://redmine.gromacs.org/issues/112005-09-12T03:50:12ZRamon Garciaramon@juguete.quim.ucm.es
<p>In ngmx viewer, when one clicks File/Quit, if one chooses "No", the dialog is<br />still there.</p>