********>Bugfix 1:
Author: Kim F. Wong
Date: 05/05/2008

Programs: sander

Description: The umbrella forces are not accumulated correctly if the atoms 
             involved in the reaction coordinate do not reside on the master
             node of a sander job.

Fix:  apply the following patch to amber10/src/sander/evb_umb.f,
                                   amber10/src/sander/evb_umb_primitive.f,
                                   amber10/src/sander/runmd.f and
                                   amber10/src/sander/sander.f
------------------------------------------------------------------------------

Index: src/sander/evb_umb.f
===================================================================
RCS file: /thr/loyd/case/cvsroot/amber10/src/sander/evb_umb.f,v
retrieving revision 1.12
diff -c -r1.12 evb_umb.f
*** src/sander/evb_umb.f	7 Apr 2008 16:35:55 -0000	1.12
--- src/sander/evb_umb.f	30 Apr 2008 20:21:50 -0000
***************
*** 5,11 ****
  
  #  include "dprec.h"
  
!    subroutine evb_umb ( f, q, mass, natom )
  
     use constants, only: kB
     use evb_parm,  only: k_umb, r0_umb, evb_dyn, nbias, dbonds_RC, bond_RC &
--- 5,11 ----
  
  #  include "dprec.h"
  
!    subroutine evb_umb ( f, q, mass, natom, istart3, iend3 )
  
     use constants, only: kB
     use evb_parm,  only: k_umb, r0_umb, evb_dyn, nbias, dbonds_RC, bond_RC &
***************
*** 21,27 ****
  
  #include "md.h"
  
!    integer, intent(   in) :: natom
     _REAL_ , intent(   in) :: q(natom*3)
     _REAL_ , intent(   in) :: mass(natom)
     _REAL_ , intent(inout) :: f(natom*3)
--- 21,27 ----
  
  #include "md.h"
  
!    integer, intent(   in) :: natom, istart3, iend3
     _REAL_ , intent(   in) :: q(natom*3)
     _REAL_ , intent(   in) :: mass(natom)
     _REAL_ , intent(inout) :: f(natom*3)
***************
*** 110,116 ****
           enddo
  
           do n = 1, nbias
!             f(:) = f(:) + evb_fbias(:,n)
              evb_frc%evb_nrg  = evb_frc%evb_nrg + evb_bias%nrg_bias(n)
           enddo
  
--- 110,118 ----
           enddo
  
           do n = 1, nbias
!             do m = istart3, iend3
!                f(m) = f(m) + evb_fbias(m,n)
!             enddo
              evb_frc%evb_nrg  = evb_frc%evb_nrg + evb_bias%nrg_bias(n)
           enddo
  
***************
*** 173,179 ****
           enddo
  
           do n = 1, nbias
!             f(:) = f(:) + evb_fbias(:,n)
              evb_frc%evb_nrg  = evb_frc%evb_nrg + evb_bias%nrg_bias(n)
           enddo
  
--- 175,183 ----
           enddo
  
           do n = 1, nbias
!             do m = istart3, iend3
!                f(m) = f(m) + evb_fbias(m,n)
!             enddo
              evb_frc%evb_nrg  = evb_frc%evb_nrg + evb_bias%nrg_bias(n)
           enddo
  
Index: src/sander/evb_umb_primitive.f
===================================================================
RCS file: /thr/loyd/case/cvsroot/amber10/src/sander/evb_umb_primitive.f,v
retrieving revision 1.5
diff -c -r1.5 evb_umb_primitive.f
*** src/sander/evb_umb_primitive.f	7 Apr 2008 16:35:55 -0000	1.5
--- src/sander/evb_umb_primitive.f	30 Apr 2008 20:21:50 -0000
***************
*** 6,12 ****
  #include "assert.h"
  #include "dprec.h"
  
!    subroutine evb_umb_primitive ( f, q, mass, natom )
  
     use constants, only: kB
     use evb_parm,  only: k_umb, r0_umb, evb_dyn, nbias, dbonds_RC, bond_RC &
--- 6,12 ----
  #include "assert.h"
  #include "dprec.h"
  
!    subroutine evb_umb_primitive ( f, q, mass, natom, istart, iend )
  
     use constants, only: kB
     use evb_parm,  only: k_umb, r0_umb, evb_dyn, nbias, dbonds_RC, bond_RC &
***************
*** 26,32 ****
  #  include "files.h"
  #  include "md.h"
  
!    integer, intent(in   ) :: natom
     _REAL_ , intent(in   ) :: q(3,natom)
     _REAL_ , intent(in   ) :: mass(natom)
     _REAL_ , intent(inout) :: f(3,natom)
--- 26,32 ----
  #  include "files.h"
  #  include "md.h"
  
!    integer, intent(in   ) :: natom, istart, iend
     _REAL_ , intent(in   ) :: q(3,natom)
     _REAL_ , intent(in   ) :: mass(natom)
     _REAL_ , intent(inout) :: f(3,natom)
***************
*** 37,43 ****
     integer :: idx, jdx, kdx
  #ifdef LES
     _REAL_  :: q_centroid(3,natomCL), fharm(3,natomCL), evb_fbias(3,natomCL,nbias) &
!             , q_bead(3,natomCL)
  #else
     _REAL_  :: fharm(3,natom), evb_fbias(3,natom,nbias) 
  #endif /* LES */
--- 37,43 ----
     integer :: idx, jdx, kdx
  #ifdef LES
     _REAL_  :: q_centroid(3,natomCL), fharm(3,natomCL), evb_fbias(3,natomCL,nbias) &
!             , q_bead(3,natomCL), f_les(3,natom)
  #else
     _REAL_  :: fharm(3,natom), evb_fbias(3,natom,nbias) 
  #endif /* LES */
***************
*** 179,192 ****
           do n = 1, nbias
              evb_frc%evb_nrg  = evb_frc%evb_nrg + evb_bias%nrg_bias(n)
  #ifdef LES
              do nn = 1, nbead
                 do m = 1, natomCL
                    mm = bead_dcrypt(m,nn)
!                   f(:,mm) = f(:,mm) + evb_fbias(:,m,n)
                 enddo
              enddo
  #else
!             f(:,:) = f(:,:) + evb_fbias(:,:,n)
  #endif /* LES */
           enddo
  
--- 179,198 ----
           do n = 1, nbias
              evb_frc%evb_nrg  = evb_frc%evb_nrg + evb_bias%nrg_bias(n)
  #ifdef LES
+             f_les(:,:) = 0.0d0
              do nn = 1, nbead
                 do m = 1, natomCL
                    mm = bead_dcrypt(m,nn)
!                   f_les(:,mm) = f_les(:,mm) + evb_fbias(:,m,n)
                 enddo
              enddo
+             do m = istart, iend
+                f(:,m) = f(:,m) + f_les(:,m)
+             enddo
  #else
!             do m = istart, iend
!                f(:,m) = f(:,m) + evb_fbias(:,m,n)
!             enddo
  #endif /* LES */
           enddo
  
***************
*** 269,282 ****
           do n = 1, nbias
              evb_frc%evb_nrg  = evb_frc%evb_nrg + evb_bias%nrg_bias(n)
  #ifdef LES
              do nn = 1, nbead
                 do m = 1, natomCL
                    mm = bead_dcrypt(m,nn)
!                   f(:,mm) = f(:,mm) + evb_fbias(:,m,n)
                 enddo
              enddo
  #else
!             f(:,:) = f(:,:) + evb_fbias(:,:,n)
  #endif /* LES */
           enddo
  
--- 275,294 ----
           do n = 1, nbias
              evb_frc%evb_nrg  = evb_frc%evb_nrg + evb_bias%nrg_bias(n)
  #ifdef LES
+             f_les(:,:) = 0.0d0
              do nn = 1, nbead
                 do m = 1, natomCL
                    mm = bead_dcrypt(m,nn)
!                   f_les(:,mm) = f_les(:,mm) + evb_fbias(:,m,n)
                 enddo
              enddo
+             do m = istart, iend
+                f(:,m) = f(:,m) + f_les(:,m)
+             enddo
  #else
!             do m = istart, iend
!                f(:,m) = f(:,m) + evb_fbias(:,m,n)
!             enddo
  #endif /* LES */
           enddo
  
Index: src/sander/runmd.f
===================================================================
RCS file: /thr/loyd/case/cvsroot/amber10/src/sander/runmd.f,v
retrieving revision 9.82
diff -c -r9.82 runmd.f
*** src/sander/runmd.f	7 Apr 2008 16:35:56 -0000	9.82
--- src/sander/runmd.f	30 Apr 2008 20:21:50 -0000
***************
*** 582,588 ****
                    do_list_update)
  #if defined(MPI) && defined(LES)
           if ( ievb == 1 .and. i_qi > 0) then
!             call evb_umb ( f, cartpos, real_mass(1:natom), natom )
              if( i_qi == 2 ) call qi_corrf_les ( cartpos, amass )
              evb_nrg(1) = evb_frc%evb_nrg
              evb_nrg(2) = evb_vel0%evb_nrg
--- 582,588 ----
                    do_list_update)
  #if defined(MPI) && defined(LES)
           if ( ievb == 1 .and. i_qi > 0) then
!             call evb_umb ( f, cartpos, real_mass, natom, istart3, iend3 )
              if( i_qi == 2 ) call qi_corrf_les ( cartpos, amass )
              evb_nrg(1) = evb_frc%evb_nrg
              evb_nrg(2) = evb_vel0%evb_nrg
***************
*** 595,601 ****
  
  #if defined(MPI) && defined(LES)
           if ( ievb /= 0 .and. i_qi == 0 ) then
!             call evb_umb ( f, x, real_mass(1:natom), natom )
              evb_nrg(1) = evb_frc%evb_nrg
              evb_nrg(2) = evb_vel0%evb_nrg
              if( nbias > 0 ) evb_nrg(3) = sum( evb_bias%nrg_bias(:) )
--- 595,601 ----
  
  #if defined(MPI) && defined(LES)
           if ( ievb /= 0 .and. i_qi == 0 ) then
!             call evb_umb ( f, x, real_mass, natom, istart3, iend3 )
              evb_nrg(1) = evb_frc%evb_nrg
              evb_nrg(2) = evb_vel0%evb_nrg
              if( nbias > 0 ) evb_nrg(3) = sum( evb_bias%nrg_bias(:) )
***************
*** 660,668 ****
  #ifdef MPI
           if ( ievb /= 0 ) then
  #ifdef LES
!             call evb_umb_primitive ( f(1:3*natom), x(1:3*natom), real_mass(1:natom), natom )
  #else
!             call evb_umb_primitive ( f(1:3*natom), x(1:3*natom), amass(1:natom), natom )
  #endif
              evb_nrg(1) = evb_frc%evb_nrg
              evb_nrg(2) = evb_vel0%evb_nrg
--- 660,668 ----
  #ifdef MPI
           if ( ievb /= 0 ) then
  #ifdef LES
!             call evb_umb_primitive ( f, x, real_mass, natom, istart, iend )
  #else
!             call evb_umb_primitive ( f, x, amass, natom, istart, iend )
  #endif
              evb_nrg(1) = evb_frc%evb_nrg
              evb_nrg(2) = evb_vel0%evb_nrg
***************
*** 1146,1152 ****
  
  #if defined(MPI) && defined(LES)
           if ( ievb == 1 .and. i_qi > 0) then
!             call evb_umb ( f, cartpos, real_mass(1:natom), natom )
              if( i_qi == 2 ) call qi_corrf_les ( cartpos, amass )
              evb_nrg(1) = evb_frc%evb_nrg
              evb_nrg(2) = evb_vel0%evb_nrg
--- 1146,1152 ----
  
  #if defined(MPI) && defined(LES)
           if ( ievb == 1 .and. i_qi > 0) then
!             call evb_umb ( f, cartpos, real_mass, natom, istart3, iend3 )
              if( i_qi == 2 ) call qi_corrf_les ( cartpos, amass )
              evb_nrg(1) = evb_frc%evb_nrg
              evb_nrg(2) = evb_vel0%evb_nrg
***************
*** 1158,1164 ****
  
  #if defined(MPI) && defined(LES)
           if ( ievb /= 0 .and. i_qi == 0 ) then
!             call evb_umb ( f, x, real_mass(1:natom), natom )
              evb_nrg(1) = evb_frc%evb_nrg
              evb_nrg(2) = evb_vel0%evb_nrg
              if( nbias > 0 ) evb_nrg(3) = sum( evb_bias%nrg_bias(:) )
--- 1158,1164 ----
  
  #if defined(MPI) && defined(LES)
           if ( ievb /= 0 .and. i_qi == 0 ) then
!             call evb_umb ( f, x, real_mass, natom, istart3, iend3 )
              evb_nrg(1) = evb_frc%evb_nrg
              evb_nrg(2) = evb_vel0%evb_nrg
              if( nbias > 0 ) evb_nrg(3) = sum( evb_bias%nrg_bias(:) )
***************
*** 1178,1186 ****
  #if defined(MPI)
           if ( ievb /= 0 ) then
  #ifdef LES
!             call evb_umb_primitive ( f(1:3*natom), x(1:3*natom), real_mass(1:natom), natom )
  #else
!             call evb_umb_primitive ( f(1:3*natom), x(1:3*natom), amass(1:natom), natom )
  #endif
              evb_nrg(1) = evb_frc%evb_nrg
              evb_nrg(2) = evb_vel0%evb_nrg
--- 1178,1186 ----
  #if defined(MPI)
           if ( ievb /= 0 ) then
  #ifdef LES
!             call evb_umb_primitive ( f, x, real_mass, natom, istart, iend )
  #else
!             call evb_umb_primitive ( f, x, amass, natom, istart, iend )
  #endif
              evb_nrg(1) = evb_frc%evb_nrg
              evb_nrg(2) = evb_vel0%evb_nrg
Index: src/sander/sander.f
===================================================================
RCS file: /thr/loyd/case/cvsroot/amber10/src/sander/sander.f,v
retrieving revision 9.61
diff -c -r9.61 sander.f
*** src/sander/sander.f	27 Mar 2008 15:47:22 -0000	9.61
--- src/sander/sander.f	30 Apr 2008 20:21:50 -0000
***************
*** 697,702 ****
--- 697,711 ----
  
     call mpi_bcast ( ievb, 1, MPI_INTEGER, 0, commworld, ierr )
     call mpi_bcast ( ipimd, 1, MPI_INTEGER, 0, commworld, ierr )
+ !KFW
+    if( ievb /= 0 ) then
+       call mpi_bcast ( evbin, 256, MPI_CHARACTER, 0, commworld, ierr )
+       if( .not. master ) then
+          call evb_input
+          call evb_init
+       endif
+    endif
+ !KFW
  #if defined(LES)
     call mpi_bcast ( ncopy, 1, MPI_INTEGER, 0, commworld, ierr )
     call mpi_bcast ( cnum(1:natom), natom, MPI_INTEGER, 0, commworld, ierr )

------------------------------------------------------------------------------

Temporary workarounds: none