subroutine etotal_long(energia)
      implicit real*8 (a-h,o-z)
      include 'DIMENSIONS'
c
c Compute the long-range slow-varying contributions to the energy
c
#ifndef ISNAN
      external proc_proc
#ifdef WINPGI
cMS$ATTRIBUTES C ::  proc_proc
#endif
#endif
#ifdef MPI
      include "mpif.h"
      double precision weights_(n_ene)
#endif
      include 'COMMON.SETUP'
      include 'COMMON.IOUNITS'
      double precision energia(0:n_ene)
      include 'COMMON.FFIELD'
      include 'COMMON.DERIV'
      include 'COMMON.INTERACT'
      include 'COMMON.SBRIDGE'
      include 'COMMON.CHAIN'
      include 'COMMON.VAR'
      include 'COMMON.LOCAL'
      include 'COMMON.MD'
c      write(iout,'(a,i2)')'Calling etotal_long ipot=',ipot
      if (modecalc.eq.12.or.modecalc.eq.14) then
#ifdef MPI
c        if (fg_rank.eq.0) call int_from_cart1(.false.)
#else
        call int_from_cart1(.false.)
#endif
      endif
#ifdef MPI      
c      write(iout,*) "ETOTAL_LONG Processor",fg_rank,
c     & " absolute rank",myrank," nfgtasks",nfgtasks
      call flush(iout)
      if (nfgtasks.gt.1) then
        time00=MPI_Wtime()
C FG slaves call the following matching MPI_Bcast in ERGASTULUM
        if (fg_rank.eq.0) then
          call MPI_Bcast(3,1,MPI_INTEGER,king,FG_COMM,IERROR)
c          write (iout,*) "Processor",myrank," BROADCAST iorder"
c          call flush(iout)
C FG master sets up the WEIGHTS_ array which will be broadcast to the 
C FG slaves as WEIGHTS array.
          weights_(1)=wsc
          weights_(2)=wscp
          weights_(3)=welec
          weights_(4)=wcorr
          weights_(5)=wcorr5
          weights_(6)=wcorr6
          weights_(7)=wel_loc
          weights_(8)=wturn3
          weights_(9)=wturn4
          weights_(10)=wturn6
          weights_(11)=wang
          weights_(12)=wscloc
          weights_(13)=wtor
          weights_(14)=wtor_d
          weights_(15)=wstrain
          weights_(16)=wvdwpp
          weights_(17)=wbond
          weights_(18)=scal14
          weights_(21)=wsccor
C FG Master broadcasts the WEIGHTS_ array
          call MPI_Bcast(weights_(1),n_ene,
     &        MPI_DOUBLE_PRECISION,king,FG_COMM,IERROR)
        else
C FG slaves receive the WEIGHTS array
          call MPI_Bcast(weights(1),n_ene,
     &        MPI_DOUBLE_PRECISION,king,FG_COMM,IERROR)
          wsc=weights(1)
          wscp=weights(2)
          welec=weights(3)
          wcorr=weights(4)
          wcorr5=weights(5)
          wcorr6=weights(6)
          wel_loc=weights(7)
          wturn3=weights(8)
          wturn4=weights(9)
          wturn6=weights(10)
          wang=weights(11)
          wscloc=weights(12)
          wtor=weights(13)
          wtor_d=weights(14)
          wstrain=weights(15)
          wvdwpp=weights(16)
          wbond=weights(17)
          scal14=weights(18)
          wsccor=weights(21)
        endif
        call MPI_Bcast(dc(1,1),6*nres,MPI_DOUBLE_PRECISION,
     &    king,FG_COMM,IERR)
         time_Bcast=time_Bcast+MPI_Wtime()-time00
         time_Bcastw=time_Bcastw+MPI_Wtime()-time00
c        call chainbuild_cart
c        call int_from_cart1(.false.)
      endif
c      write (iout,*) 'Processor',myrank,
c     &  ' calling etotal_short ipot=',ipot
c      call flush(iout)
c      print *,'Processor',myrank,' nnt=',nnt,' nct=',nct
#endif     
cd    print *,'nnt=',nnt,' nct=',nct
C
C Compute the side-chain and electrostatic interaction energy
C
      goto (101,102,103,104,105,106) ipot
C Lennard-Jones potential.
  101 call elj_long(evdw)
cd    print '(a)','Exit ELJ'
      goto 107
C Lennard-Jones-Kihara potential (shifted).
  102 call eljk_long(evdw)
      goto 107
C Berne-Pechukas potential (dilated LJ, angular dependence).
  103 call ebp_long(evdw)
      goto 107
C Gay-Berne potential (shifted LJ, angular dependence).
  104 call egb_long(evdw,evdw_p,evdw_m)
      goto 107
C Gay-Berne-Vorobjev potential (shifted LJ, angular dependence).
  105 call egbv_long(evdw)
      goto 107
C Soft-sphere potential
  106 call e_softsphere(evdw)
C
C Calculate electrostatic (H-bonding) energy of the main chain.
C
  107 continue
      call vec_and_deriv
      if (ipot.lt.6) then
#ifdef SPLITELE
         if (welec.gt.0d0.or.wvdwpp.gt.0d0.or.wel_loc.gt.0d0.or.
     &       wturn3.gt.0d0.or.wturn4.gt.0d0 .or. wcorr.gt.0.0d0
     &       .or. wcorr4.gt.0.0d0 .or. wcorr5.gt.0.d0
     &       .or. wcorr6.gt.0.0d0 .or. wturn6.gt.0.0d0 ) then
#else
         if (welec.gt.0d0.or.wel_loc.gt.0d0.or.
     &       wturn3.gt.0d0.or.wturn4.gt.0d0 .or. wcorr.gt.0.0d0
     &       .or. wcorr4.gt.0.0d0 .or. wcorr5.gt.0.d0
     &       .or. wcorr6.gt.0.0d0 .or. wturn6.gt.0.0d0 ) then
#endif
           call eelec_scale(ees,evdw1,eel_loc,eello_turn3,eello_turn4)
         else
            ees=0
            evdw1=0
            eel_loc=0
            eello_turn3=0
            eello_turn4=0
         endif
      else
c        write (iout,*) "Soft-spheer ELEC potential"
        call eelec_soft_sphere(ees,evdw1,eel_loc,eello_turn3,
     &   eello_turn4)
      endif
C
C Calculate excluded-volume interaction energy between peptide groups
C and side chains.
C
      if (ipot.lt.6) then
       if(wscp.gt.0d0) then
        call escp_long(evdw2,evdw2_14)
       else
        evdw2=0
        evdw2_14=0
       endif
      else
        call escp_soft_sphere(evdw2,evdw2_14)
      endif
C 
C 12/1/95 Multi-body terms
C
      n_corr=0
      n_corr1=0
      if ((wcorr4.gt.0.0d0 .or. wcorr5.gt.0.0d0 .or. wcorr6.gt.0.0d0 
     &    .or. wturn6.gt.0.0d0) .and. ipot.lt.6) then
         call multibody_eello(ecorr,ecorr5,ecorr6,eturn6,n_corr,n_corr1)
c         write (2,*) 'n_corr=',n_corr,' n_corr1=',n_corr1,
c     &" ecorr",ecorr," ecorr5",ecorr5," ecorr6",ecorr6," eturn6",eturn6
      else
         ecorr=0.0d0
         ecorr5=0.0d0
         ecorr6=0.0d0
         eturn6=0.0d0
      endif
      if ((wcorr4.eq.0.0d0 .and. wcorr.gt.0.0d0) .and. ipot.lt.6) then
         call multibody_hb(ecorr,ecorr5,ecorr6,n_corr,n_corr1)
      endif
C 
C If performing constraint dynamics, call the constraint energy
C  after the equilibration time
      if(usampl.and.totT.gt.eq_time) then
         call EconstrQ   
         call Econstr_back
      else
         Uconst=0.0d0
         Uconst_back=0.0d0
      endif
C 
C Sum the energies
C
      do i=1,n_ene
        energia(i)=0.0d0
      enddo
      energia(1)=evdw
#ifdef SCP14
      energia(2)=evdw2-evdw2_14
      energia(18)=evdw2_14
#else
      energia(2)=evdw2
      energia(18)=0.0d0
#endif
#ifdef SPLITELE
      energia(3)=ees
      energia(16)=evdw1
#else
      energia(3)=ees+evdw1
      energia(16)=0.0d0
#endif
      energia(4)=ecorr
      energia(5)=ecorr5
      energia(6)=ecorr6
      energia(7)=eel_loc
      energia(8)=eello_turn3
      energia(9)=eello_turn4
      energia(10)=eturn6
      energia(20)=Uconst+Uconst_back
      energia(22)=evdw_p
      energia(23)=evdw_m
      call sum_energy(energia,.true.)
c      write (iout,*) "Exit ETOTAL_LONG"
      call flush(iout)
      return
      end
c------------------------------------------------------------------------------
      subroutine etotal_short(energia)
      implicit real*8 (a-h,o-z)
      include 'DIMENSIONS'
c
c Compute the short-range fast-varying contributions to the energy
c
#ifndef ISNAN
      external proc_proc
#ifdef WINPGI
cMS$ATTRIBUTES C ::  proc_proc
#endif
#endif
#ifdef MPI
      include "mpif.h"
      double precision weights_(n_ene)
#endif
      include 'COMMON.SETUP'
      include 'COMMON.IOUNITS'
      double precision energia(0:n_ene)
      include 'COMMON.FFIELD'
      include 'COMMON.DERIV'
      include 'COMMON.INTERACT'
      include 'COMMON.SBRIDGE'
      include 'COMMON.CHAIN'
      include 'COMMON.VAR'
      include 'COMMON.LOCAL'
      include 'COMMON.CONTROL'
#ifndef DFA
      edfadis=0.0d0
      edfator=0.0d0
      edfanei=0.0d0
      edfabet=0.0d0
#endif
c      write(iout,'(a,i2)')'Calling etotal_short ipot=',ipot
c      call flush(iout)
      if (modecalc.eq.12.or.modecalc.eq.14) then
#ifdef MPI
        if (fg_rank.eq.0) call int_from_cart1(.false.)
#else
        call int_from_cart1(.false.)
#endif
      endif
#ifdef MPI      
c      write(iout,*) "ETOTAL_SHORT Processor",fg_rank,
c     & " absolute rank",myrank," nfgtasks",nfgtasks
c      call flush(iout)
      if (nfgtasks.gt.1) then
        time00=MPI_Wtime()
C FG slaves call the following matching MPI_Bcast in ERGASTULUM
        if (fg_rank.eq.0) then
          call MPI_Bcast(2,1,MPI_INTEGER,king,FG_COMM,IERROR)
c          write (iout,*) "Processor",myrank," BROADCAST iorder"
c          call flush(iout)
C FG master sets up the WEIGHTS_ array which will be broadcast to the 
C FG slaves as WEIGHTS array.
          weights_(1)=wsc
          weights_(2)=wscp
          weights_(3)=welec
          weights_(4)=wcorr
          weights_(5)=wcorr5
          weights_(6)=wcorr6
          weights_(7)=wel_loc
          weights_(8)=wturn3
          weights_(9)=wturn4
          weights_(10)=wturn6
          weights_(11)=wang
          weights_(12)=wscloc
          weights_(13)=wtor
          weights_(14)=wtor_d
          weights_(15)=wstrain
          weights_(16)=wvdwpp
          weights_(17)=wbond
          weights_(18)=scal14
          weights_(21)=wsccor
C FG Master broadcasts the WEIGHTS_ array
          call MPI_Bcast(weights_(1),n_ene,
     &        MPI_DOUBLE_PRECISION,king,FG_COMM,IERROR)
        else
C FG slaves receive the WEIGHTS array
          call MPI_Bcast(weights(1),n_ene,
     &        MPI_DOUBLE_PRECISION,king,FG_COMM,IERROR)
          wsc=weights(1)
          wscp=weights(2)
          welec=weights(3)
          wcorr=weights(4)
          wcorr5=weights(5)
          wcorr6=weights(6)
          wel_loc=weights(7)
          wturn3=weights(8)
          wturn4=weights(9)
          wturn6=weights(10)
          wang=weights(11)
          wscloc=weights(12)
          wtor=weights(13)
          wtor_d=weights(14)
          wstrain=weights(15)
          wvdwpp=weights(16)
          wbond=weights(17)
          scal14=weights(18)
          wsccor=weights(21)
        endif
c        write (iout,*),"Processor",myrank," BROADCAST weights"
        call MPI_Bcast(c(1,1),maxres6,MPI_DOUBLE_PRECISION,
     &    king,FG_COMM,IERR)
c        write (iout,*) "Processor",myrank," BROADCAST c"
        call MPI_Bcast(dc(1,1),maxres6,MPI_DOUBLE_PRECISION,
     &    king,FG_COMM,IERR)
c        write (iout,*) "Processor",myrank," BROADCAST dc"
        call MPI_Bcast(dc_norm(1,1),maxres6,MPI_DOUBLE_PRECISION,
     &    king,FG_COMM,IERR)
c        write (iout,*) "Processor",myrank," BROADCAST dc_norm"
        call MPI_Bcast(theta(1),nres,MPI_DOUBLE_PRECISION,
     &    king,FG_COMM,IERR)
c        write (iout,*) "Processor",myrank," BROADCAST theta"
        call MPI_Bcast(phi(1),nres,MPI_DOUBLE_PRECISION,
     &    king,FG_COMM,IERR)
c        write (iout,*) "Processor",myrank," BROADCAST phi"
        call MPI_Bcast(alph(1),nres,MPI_DOUBLE_PRECISION,
     &    king,FG_COMM,IERR)
c        write (iout,*) "Processor",myrank," BROADCAST alph"
        call MPI_Bcast(omeg(1),nres,MPI_DOUBLE_PRECISION,
     &    king,FG_COMM,IERR)
c        write (iout,*) "Processor",myrank," BROADCAST omeg"
        call MPI_Bcast(vbld(1),2*nres,MPI_DOUBLE_PRECISION,
     &    king,FG_COMM,IERR)
c        write (iout,*) "Processor",myrank," BROADCAST vbld"
        call MPI_Bcast(vbld_inv(1),2*nres,MPI_DOUBLE_PRECISION,
     &    king,FG_COMM,IERR)
         time_Bcast=time_Bcast+MPI_Wtime()-time00
c        write (iout,*) "Processor",myrank," BROADCAST vbld_inv"
      endif
c      write (iout,*) 'Processor',myrank,
c     &  ' calling etotal_short ipot=',ipot
c      call flush(iout)
c      print *,'Processor',myrank,' nnt=',nnt,' nct=',nct
#endif     
c      call int_from_cart1(.false.)
C
C Compute the side-chain and electrostatic interaction energy
C
      goto (101,102,103,104,105,106) ipot
C Lennard-Jones potential.
  101 call elj_short(evdw)
cd    print '(a)','Exit ELJ'
      goto 107
C Lennard-Jones-Kihara potential (shifted).
  102 call eljk_short(evdw)
      goto 107
C Berne-Pechukas potential (dilated LJ, angular dependence).
  103 call ebp_short(evdw)
      goto 107
C Gay-Berne potential (shifted LJ, angular dependence).
  104 call egb_short(evdw,evdw_p,evdw_m)
      goto 107
C Gay-Berne-Vorobjev potential (shifted LJ, angular dependence).
  105 call egbv_short(evdw)
      goto 107
C Soft-sphere potential - already dealt with in the long-range part
  106 evdw=0.0d0
c  106 call e_softsphere_short(evdw)
C
C Calculate electrostatic (H-bonding) energy of the main chain.
C
  107 continue
#ifdef DFA
C     BARTEK for dfa test!
      if (wdfa_dist.gt.0) then 
         call edfad(edfadis)
      else
         edfadis=0.0
      endif
c      print*, 'edfad is finished!', edfadis
      if (wdfa_tor.gt.0) then 
         call edfat(edfator)
      else
         edfator=0.0
      endif
c      print*, 'edfat is finished!', edfator
      if (wdfa_nei.gt.0) then 
         call edfan(edfanei)
      else
         edfanei=0.0
      endif
c      print*, 'edfan is finished!', edfanei
      if (wdfa_beta.gt.0) then 
         call edfab(edfabet)
      else
         edfabet=0.0
      endif
c      print*, 'edfab is finished!', edfabet
#endif
c
c Calculate the short-range part of Evdwpp
c
      call evdwpp_short(evdw1)
c
c Calculate the short-range part of ESCp
c
      if (ipot.lt.6) then
        call escp_short(evdw2,evdw2_14)
      endif
c
c Calculate the bond-stretching energy
c
      call ebond(estr)
C 
C Calculate the disulfide-bridge and other energy and the contributions
C from other distance constraints.
      call edis(ehpb)
C
C Calculate the virtual-bond-angle energy.
C
      call ebend(ebe)
C
C Calculate the SC local energy.
C
      call vec_and_deriv
      call esc(escloc)
C
C Calculate the virtual-bond torsional energy.
C
      call etor(etors,edihcnstr)
c
c Homology restraints
c
      if (constr_homology.ge.1) then
        call e_modeller(ehomology_constr)
      else
        ehomology_constr=0.0d0
      endif
C
C 6/23/01 Calculate double-torsional energy
C
      call etor_d(etors_d)
C
C 21/5/07 Calculate local sicdechain correlation energy
C
      if (wsccor.gt.0.0d0) then
        call eback_sc_corr(esccor)
      else
        esccor=0.0d0
      endif
C
C Put energy components into an array
C
      do i=1,n_ene
        energia(i)=0.0d0
      enddo
      energia(1)=evdw
#ifdef SCP14
      energia(2)=evdw2-evdw2_14
      energia(18)=evdw2_14
#else
      energia(2)=evdw2
      energia(18)=0.0d0
#endif
#ifdef SPLITELE
      energia(16)=evdw1
#else
      energia(3)=evdw1
#endif
      energia(11)=ebe
      energia(12)=escloc
      energia(13)=etors
      energia(14)=etors_d
      energia(15)=ehpb
      energia(17)=estr
      energia(19)=edihcnstr
      energia(21)=esccor
      energia(22)=evdw_p
      energia(23)=evdw_m
      energia(24)=ehomology_constr
      energia(25)=edfadis
      energia(26)=edfator
      energia(27)=edfanei
      energia(28)=edfabet
c      write (iout,*) "ETOTAL_SHORT before SUM_ENERGY"
      call flush(iout)
      call sum_energy(energia,.true.)
c      write (iout,*) "Exit ETOTAL_SHORT"
      call flush(iout)
      return
      end