********>Bugfix.33:
Author: Dave Case
Date: 02/10/2005

Programs: sander

Description: When the cap water option is requested, there is no energy 
             computed to go along with the force being applied.  This can
             lead to an apparent lack of energy conservation.

Fix:  apply the following patch to amber8/src/sander/{ene.f,force.f}

------------------------------------------------------------------------------
*** ene.f	2004/03/01 15:17:37	7.71
--- ene.f	2005/02/10 22:10:24
***************
*** 1028,1034 ****
  
  !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  !+ [Enter a one-line description of subroutine capwat here]
! subroutine capwat(nat,x,f)
  
     implicit _REAL_ (a-h,o-z)
  #ifdef MPI
--- 1028,1034 ----
  
  !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
  !+ [Enter a one-line description of subroutine capwat here]
! subroutine capwat(nat,x,f,ecap)
  
     implicit _REAL_ (a-h,o-z)
  #ifdef MPI
***************
*** 1042,1047 ****
--- 1042,1049 ----
     dimension x(3,*),f(3,*)
     data tm34,zero/1.0d-34,0.0d0/
  
+    ecap = 0.d0
+ 
     ! FLOW CONTROL FLAG (debug)
     if ( do_cap == 0 )return
  #ifdef MPI
***************
*** 1054,1060 ****
        ya = ycap-x(2,i)
        za = zcap-x(3,i)
        da = sqrt(xa*xa+ya*ya+za*za+tm34)
!       df = fcap*max(zero,da-cutcap)/da
        f(1,i) = f(1,i)+df*xa
        f(2,i) = f(2,i)+df*ya
        f(3,i) = f(3,i)+df*za
--- 1056,1064 ----
        ya = ycap-x(2,i)
        za = zcap-x(3,i)
        da = sqrt(xa*xa+ya*ya+za*za+tm34)
!       delta = max(zero,da-cutcap)
!       ecap = ecap + 0.5*fcap*delta*delta
!       df = fcap*delta/da
        f(1,i) = f(1,i)+df*xa
        f(2,i) = f(2,i)+df*ya
        f(3,i) = f(3,i)+df*za

------------------------------------------------------------------------------
*** force.f	2004/03/09 00:59:15	7.172
--- force.f	2005/02/10 22:10:34
***************
*** 63,69 ****
     _REAL_  r_stack(*)
     integer i_stack(*)
  
!    _REAL_  enmr(3),devdis(4),devang(4),devtor(4),entr
  
     dimension itrp(3,nt_lim)
     _REAL_  x(*),f(*),ene(30),vir(*)
--- 63,69 ----
     _REAL_  r_stack(*)
     integer i_stack(*)
  
!    _REAL_  enmr(3),devdis(4),devang(4),devtor(4),entr,ecap
  
     dimension itrp(3,nt_lim)
     _REAL_  x(*),f(*),ene(30),vir(*)
***************
*** 524,530 ****
  
     end if
  
!    if(ifcap > 0) call capwat(natom,x,f)
  
  #ifdef MPI
     if(mytaskid == 0)then
--- 524,533 ----
  
     end if
  
!    if(ifcap > 0) then
!       call capwat(natom,x,f,ecap)
!       ene(20) = ene(20) + ecap
!    end if
  
  #ifdef MPI
     if(mytaskid == 0)then
***************
*** 725,731 ****
     !    ene(9):    1-4 electrostatics
     !    ene(10):   constraint energy
     !    ene(11-19):  used a scratch, but not needed further below
!    !    ene(20):   position constraint energy
     !    ene(21):   charging free energy result
     !    ene(22):   noe volume penalty
     !    ene(23):   surface-area dependent energy
--- 728,734 ----
     !    ene(9):    1-4 electrostatics
     !    ene(10):   constraint energy
     !    ene(11-19):  used a scratch, but not needed further below
!    !    ene(20):   position constraint energy + cap energy
     !    ene(21):   charging free energy result
     !    ene(22):   noe volume penalty
     !    ene(23):   surface-area dependent energy
  
------------------------------------------------------------------------------

Workarounds: Don't worry, be happy....Amber has been like this forever.