********> bugfix.26, rev.A Correction Author: David A. Pearlman Correction Date: 12/23/91 Revision Date: 1/8/92 Revision notice: The original version of bugfix.26 was flawed. This revised version will effect the proper fix. Replace bugfix.26 with this bugfix report. Programs: minmd/sander/gibbs Severity: moderate Problem: If the user attempts to run a dynamics simulation on a system where 1) the solute consists of < 3 atoms 2) the center of mass motion of the solute is being removed (ntcm >0 or nscm >0) The trajectory will likely become unstable, and/or will be unreliable. This bug will not affect any case except that listed above (solute of one or two atoms _and_ center of mass motion of the solute is being removed). Affects: The reported temperature may become unstable, or unexpected failures in shake may occur. Cause: The value of NDFMIN (the number of degrees of freedom that are to be subtracted from the total when center of mass motion is being removed) was set to 6 when center of mass motion was removed (ntcm =1 or nscm >0). This would result in an inappropriate number of degrees of freedom for one and two atom solutes. These small solutes actually present special cases, which should be dealt with separately. Fix: In routine MDREAD in src/minmd/mdread.f src/sander/mdread.f src/gibbs/giba.f make the following change: --old lines (same in MINMD/SANDER/GIBBS): if (ntcm.ne.0 .or. nscm.gt.0) then ndfmin = 6 else ndfmin = 0 endif --new lines (in SANDER/GIBBS): NSOLUT = NRP IF ((NTT.EQ.2 .OR. NTT.EQ.5)) THEN IF (NSPSTR.GT.0) NSOLUT = NSPSTR IF (ISOLVP.GT.0) NSOLUT = ISOLVP END IF if ((ntcm.eq.1 .or. nscm.gt.0) .and. ntb.eq.0) then ndfmin = 6 if (nsolut.eq.1) ndfmin = 3 if (nsolut.eq.2) ndfmin = 5 else ndfmin = 0 endif c the following 'if..endif' should be added to gibbs (already in sander) if(ibelly .eq. 1) then c --- No COM Motion Removal, ever. --- ntcm = 0 nscm = 9999999 ndfmin = 0 endif --new lines (in MINMD): if ((ntcm.eq.1 .or. nscm.gt.0) .and. ntb.eq.0) then ndfmin = 6 nrpu = nrp if (ibelly.gt.0 .and. ntt.lt.2) nrpu = natom if (nrpu.eq.1) ndfmin = 3 if (nrpu.eq.2) ndfmin = 5 else ndfmin = 0 endif --new lines (denoted by '+'; in SANDER, GIBBS): if (ntcm.ne.0 .and. ntcm.ne.1) then write(6,'(/2x,a,i3,a)') 'NTCM (',ntcm,') must be 0 or 1.' inerr = 1 endif + if (nsolut.eq.0 .and. ndfmin.ne.0) then + write(6,'(/2x,a)') 'NTCM or NSCM set but no solute atoms' + inerr = 1 + endif --new lines (denoted by '+'; in MINMD): if (ntcm.ne.0 .and. ntcm.ne.1) then write(6,'(/2x,a,i3,a)') 'NTCM (', ntcm, ') must be 0 or 1.' inprob = 1 endif + if (ndfmin.ne.0 .and. nrpu.eq.0) then + write(6,'(/2x,a)') 'NTCM or NSCM set but no solute atoms' + inprob = 1 + endif --also, delete the following lines in MINMD mdread.f: if (ndfmin.ne.0 .and. ndfmin.ne.6) then write(6,'(/2x,a,i3,a)') 'NDFMIN (', ndfmin, ') must be 0 or 6.' inprob = 1 endif Temporary workarounds: None, when the above conditions apply to the simulation. However, for a one or two atom system in water, the center of mass motion of the solute should probably NOT be removed periodically for most applications. This would avoid the problem. Routines affected: MINMD : Routine MDREAD in file .../amber4/src/minmd/mdread.f SANDER: Routine MDREAD in file .../amber4/src/sander/mdread.f GIBBS : Routine MDREAD in file .../amber4/src/gibbs/giba.f