AMBER Archive (2006)

Subject: AMBER: problems trying to modify Amber8

From: Ed Pate (pate_at_math.wsu.edu)
Date: Thu Dec 28 2006 - 16:41:33 CST

Dear Amber community:

I am trying to modify sander in amber8 to output the coordinates of 4
atoms much more frequently than the .crd output of the entire protein in a
GB simulation. I have modified the source code, using locations of the
parameters ntwx (crd output frequency) and mdcrd (crd output file) in
Amber8 as my guides. Details of the modifications are below. The
modified code runs fine and outputs the desired 4 atom coordinates and
.crd file for the first 1000 or so time steps (tested on shorter runs).
The program then hangs. CPU usage appears to indicate that one processor
drops out as it hangs.

I am running on an Opteron-275 based system, Portland group fortran and
infiniband interconnects.

I seem to be missing something fundamental that I need to take care of in
the code. I would be grateful for any help that anyone could give.

Thanks,

Ed Pate

Details of the modification to the code:

In subroutine dynlib.f

1. The output file for the 4 coordinates is opened by modifying the code
as follows, an open statement instead of a call to amopen as for mdcrd.

    ! subr amopen(lun,fname,fstat,fform,facc) (line 140 of amber8
release)

    if ( master ) then
       if (ioutfm <= 0) then

          ! ----- FORMATTED DUMPING -----
! ADDED CODE FOLLOWS+++++++++++++++++++++++++++++
! patch e.p. 12/28/06 for epr output.
! open unit 22 as epr coordinate output
          write(6,*)" unit 22 opened**********************" ! debug
info
          open(unit=22,file='/home/pate/sl3_myosin.d/eprout.dat', &
              status='replace',form='formatted')
! end patch
!+++++++++++++++++++++++++++++++++++++++++++++++++++
          if (ntwx > 0) then
             if (imin == 5) then
                call amopen(mdcrd_unit,mdcrd,'O','F','R')
                read(mdcrd_unit,1000) title
                write (6,1000) title
             else

2. In dynlib.f, a new subroutine corpac_epr is added to output the
coordinates of the 4 atoms. The new subroutine is added because I need
greater decimal significance than provided by corpac itself. The
arguments are the same as corpac, except the final argument is the length
of the .crd file instead of a logical for formatted output. In the amber8
source code, the subroutine follows immediately after corpac. Using
corpac instead of corpac_epr in test runs does not change the outcome of
the simulation.

!-------------------------------------------------------------------
! patch 12/28/06 e.p. subroutine corpac_epr for epr output
! arguments the same as for corpac, except last argument = length
! of array r
!
subroutine corpac_epr(r,n1,ne1,nf,nlength)
!
!
       implicit none
       integer n1,ne1,nlength,nf
       _REAL_ r(nlength)
       integer j
       write(nf,'(12f12.7)') (r(j),j=n1,ne1)
       return
end subroutine corpac_epr
!-----------------------------------------------------------------------

3. In runmd.f a logical variable iedump is set. It mimicks the logical
variable itdump that controls output of the .crd file. iedump was
declared logical.

! --- determine if velocity or restart writing is imminent or (line
1393 in Amber release)
! if the center of mass motion will be removed;
! these require xdist of velocities or dipoles in parallel runs:

    ivarch = .false.
    ixdump = .false.
    itdump = .false.
! patch 12/28/06 for epr output++++++++++++++++++++++++++++
    iedump = .false.
! end patch++++++++++++++++++++++++++++++++++++++++++

Intervening code follows and iedump is set following protocol for itdump
for coordinate dumping. Parameter ntwepr sets the number of time steps
between output of the desired 4 atoms . equivalent to ntwx for the .crd
file. The value of ntwepr is set as a parameter in runmd.f. Values of 10
and 100 give the same result.

    ! --- determine if trajectory writing is imminent: (line 1418 in
amber release)

    if(ntwx > 0)then
       itdump= mod(nstep+1,ntwx) == 0
    else
       itdump=.false.
    end if
! patch for epr output+++++++++++++++++++++++++++++++++++++++
    if(ntwepr > 0)then
       iedump= mod(nstep+1,ntwepr) == 0
    else
       iedump=.false.
    end if
! end patch++++++++++++++++++++++++++++++++++++++++++++++

4. Testing for coordinate dumping and output of the 4-atom file follows
that for the .crd file in runmd.f. The parameter nepr used in setting up
the call to corpac_epr is the atom number of the first of the 4 sequential
atoms in the pdb file for which coordinates are to be outputted. It is
set as a parameter in runmd.f.

       ! -- Coordinate archive: (line 1821 in Amber release)

       if (itdump) then
          if( iwrap == 0 ) then
             call corpac(x,1,3*natom,12,loutfm)
             if(ntb > 0) call corpac(box,1,3,12,loutfm)
          else
             call get_stack(l_temp,3*natom)
             do m=1,3*natom
                r_stack(l_temp + m - 1) = x(m)
             end do
#ifndef QMMM
             if(ntb > 0) then
                call wrap_molecules(nspm,nsp,r_stack(l_temp))
                if (ifbox == 2) call wrap_to(nspm,nsp,r_stack(l_temp),box)
             end if
#endif
             call corpac(r_stack(l_temp),1,3*natom,12,loutfm)
             if(ntb > 0) call corpac(box,1,3,12,loutfm)
             call free_stack(l_temp)
          end if
       end if ! (itdump)
! MODIFICATION TO AMBER++++++++++++++++++++++++++++++
! patch e.p. 12/28/06 for epr output
       if (iedump) then
            nx1=(nepr-1)*3+1
            nend1=nx1+11
            call corpac_epr(x,nx1,nend1,22,3*natom)
       end if
! end patch
! END MODIFICATION++++++++++++++++++++++++++++++++++

The md.in file follows:

&cntrl
   ntc=2,
   ntf=2,
   cut=16.0,
   igb=2,
   saltcon=0.2,
   gbsa=1,
   ntpr=250,
   nstlim=10000,
   dt=0.002,
   ntt=1,
   tempi=0.0,
   temp0=297.,
   tautp=1.0,
   ntx=1,
   irest=0,
   ntb=0,
   ntwx=500,
   nscm=1000,
   lastist=30000000,
   lastrst=30000000,
&end

Thanks again for the help.

Ed

-----------------------------------------------------------------------
The AMBER Mail Reflector
To post, send mail to amber_at_scripps.edu
To unsubscribe, send "unsubscribe amber" to majordomo_at_scripps.edu