source/unres/src_MD/energy_split-sep.F

   1       subroutine etotal_long(energia)
   2       implicit real*8 (a-h,o-z)
   3       include 'DIMENSIONS'
   4 c
   5 c Compute the long-range slow-varying contributions to the energy
   6 c
   7 #ifndef ISNAN
   8       external proc_proc
   9 #ifdef WINPGI
  10 cMS$ATTRIBUTES C ::  proc_proc
  11 #endif
  12 #endif
  13 #ifdef MPI
  14       include "mpif.h"
  15       double precision weights_(n_ene)
  16 #endif
  17       include 'COMMON.SETUP'
  18       include 'COMMON.IOUNITS'
  19       double precision energia(0:n_ene)
  20       include 'COMMON.FFIELD'
  21       include 'COMMON.DERIV'
  22       include 'COMMON.INTERACT'
  23       include 'COMMON.SBRIDGE'
  24       include 'COMMON.CHAIN'
  25       include 'COMMON.VAR'
  26       include 'COMMON.LOCAL'
  27       include 'COMMON.MD'
  28 c      write(iout,'(a,i2)')'Calling etotal_long ipot=',ipot
  29       if (modecalc.eq.12.or.modecalc.eq.14) then
  30 #ifdef MPI
  31 c        if (fg_rank.eq.0) call int_from_cart1(.false.)
  32 #else
  33         call int_from_cart1(.false.)
  34 #endif
  35       endif
  36 #ifdef MPI
  37 c      write(iout,*) "ETOTAL_LONG Processor",fg_rank,
  38 c     & " absolute rank",myrank," nfgtasks",nfgtasks
  39       call flush(iout)
  40       if (nfgtasks.gt.1) then
  41         time00=MPI_Wtime()
  42 C FG slaves call the following matching MPI_Bcast in ERGASTULUM
  43         if (fg_rank.eq.0) then
  44           call MPI_Bcast(3,1,MPI_INTEGER,king,FG_COMM,IERROR)
  45 c          write (iout,*) "Processor",myrank," BROADCAST iorder"
  46 c          call flush(iout)
  47 C FG master sets up the WEIGHTS_ array which will be broadcast to the
  48 C FG slaves as WEIGHTS array.
  49           weights_(1)=wsc
  50           weights_(2)=wscp
  51           weights_(3)=welec
  52           weights_(4)=wcorr
  53           weights_(5)=wcorr5
  54           weights_(6)=wcorr6
  55           weights_(7)=wel_loc
  56           weights_(8)=wturn3
  57           weights_(9)=wturn4
  58           weights_(10)=wturn6
  59           weights_(11)=wang
  60           weights_(12)=wscloc
  61           weights_(13)=wtor
  62           weights_(14)=wtor_d
  63           weights_(15)=wstrain
  64           weights_(16)=wvdwpp
  65           weights_(17)=wbond
  66           weights_(18)=scal14
  67           weights_(21)=wsccor
  68 C FG Master broadcasts the WEIGHTS_ array
  69           call MPI_Bcast(weights_(1),n_ene,
  70      &        MPI_DOUBLE_PRECISION,king,FG_COMM,IERROR)
  71         else
  72 C FG slaves receive the WEIGHTS array
  73           call MPI_Bcast(weights(1),n_ene,
  74      &        MPI_DOUBLE_PRECISION,king,FG_COMM,IERROR)
  75           wsc=weights(1)
  76           wscp=weights(2)
  77           welec=weights(3)
  78           wcorr=weights(4)
  79           wcorr5=weights(5)
  80           wcorr6=weights(6)
  81           wel_loc=weights(7)
  82           wturn3=weights(8)
  83           wturn4=weights(9)
  84           wturn6=weights(10)
  85           wang=weights(11)
  86           wscloc=weights(12)
  87           wtor=weights(13)
  88           wtor_d=weights(14)
  89           wstrain=weights(15)
  90           wvdwpp=weights(16)
  91           wbond=weights(17)
  92           scal14=weights(18)
  93           wsccor=weights(21)
  94         endif
  95         call MPI_Bcast(dc(1,1),6*nres,MPI_DOUBLE_PRECISION,
  96      &    king,FG_COMM,IERR)
  97          time_Bcast=time_Bcast+MPI_Wtime()-time00
  98          time_Bcastw=time_Bcastw+MPI_Wtime()-time00
  99 c        call chainbuild_cart
 100 c        call int_from_cart1(.false.)
 101       endif
 102 c      write (iout,*) 'Processor',myrank,
 103 c     &  ' calling etotal_short ipot=',ipot
 104 c      call flush(iout)
 105 c      print *,'Processor',myrank,' nnt=',nnt,' nct=',nct
 106 #endif
 107 cd    print *,'nnt=',nnt,' nct=',nct
 108 C
 109 C Compute the side-chain and electrostatic interaction energy
 110 C
 111       goto (101,102,103,104,105,106) ipot
 112 C Lennard-Jones potential.
 113   101 call elj_long(evdw)
 114 cd    print '(a)','Exit ELJ'
 115       goto 107
 116 C Lennard-Jones-Kihara potential (shifted).
 117   102 call eljk_long(evdw)
 118       goto 107
 119 C Berne-Pechukas potential (dilated LJ, angular dependence).
 120   103 call ebp_long(evdw)
 121       goto 107
 122 C Gay-Berne potential (shifted LJ, angular dependence).
 123   104 call egb_long(evdw,evdw_p,evdw_m)
 124       goto 107
 125 C Gay-Berne-Vorobjev potential (shifted LJ, angular dependence).
 126   105 call egbv_long(evdw)
 127       goto 107
 128 C Soft-sphere potential
 129   106 call e_softsphere(evdw)
 130 C
 131 C Calculate electrostatic (H-bonding) energy of the main chain.
 132 C
 133   107 continue
 134       call vec_and_deriv
 135       if (ipot.lt.6) then
 136 #ifdef SPLITELE
 137          if (welec.gt.0d0.or.wvdwpp.gt.0d0.or.wel_loc.gt.0d0.or.
 138      &       wturn3.gt.0d0.or.wturn4.gt.0d0 .or. wcorr.gt.0.0d0
 139      &       .or. wcorr4.gt.0.0d0 .or. wcorr5.gt.0.d0
 140      &       .or. wcorr6.gt.0.0d0 .or. wturn6.gt.0.0d0 ) then
 141 #else
 142          if (welec.gt.0d0.or.wel_loc.gt.0d0.or.
 143      &       wturn3.gt.0d0.or.wturn4.gt.0d0 .or. wcorr.gt.0.0d0
 144      &       .or. wcorr4.gt.0.0d0 .or. wcorr5.gt.0.d0
 145      &       .or. wcorr6.gt.0.0d0 .or. wturn6.gt.0.0d0 ) then
 146 #endif
 147            call eelec_scale(ees,evdw1,eel_loc,eello_turn3,eello_turn4)
 148          else
 149             ees=0
 150             evdw1=0
 151             eel_loc=0
 152             eello_turn3=0
 153             eello_turn4=0
 154          endif
 155       else
 156 c        write (iout,*) "Soft-spheer ELEC potential"
 157         call eelec_soft_sphere(ees,evdw1,eel_loc,eello_turn3,
 158      &   eello_turn4)
 159       endif
 160 C
 161 C Calculate excluded-volume interaction energy between peptide groups
 162 C and side chains.
 163 C
 164       if (ipot.lt.6) then
 165        if(wscp.gt.0d0) then
 166         call escp_long(evdw2,evdw2_14)
 167        else
 168         evdw2=0
 169         evdw2_14=0
 170        endif
 171       else
 172         call escp_soft_sphere(evdw2,evdw2_14)
 173       endif
 174 C
 175 C 12/1/95 Multi-body terms
 176 C
 177       n_corr=0
 178       n_corr1=0
 179       if ((wcorr4.gt.0.0d0 .or. wcorr5.gt.0.0d0 .or. wcorr6.gt.0.0d0
 180      &    .or. wturn6.gt.0.0d0) .and. ipot.lt.6) then
 181          call multibody_eello(ecorr,ecorr5,ecorr6,eturn6,n_corr,n_corr1)
 182 c         write (2,*) 'n_corr=',n_corr,' n_corr1=',n_corr1,
 183 c     &" ecorr",ecorr," ecorr5",ecorr5," ecorr6",ecorr6," eturn6",eturn6
 184       else
 185          ecorr=0.0d0
 186          ecorr5=0.0d0
 187          ecorr6=0.0d0
 188          eturn6=0.0d0
 189       endif
 190       if ((wcorr4.eq.0.0d0 .and. wcorr.gt.0.0d0) .and. ipot.lt.6) then
 191          call multibody_hb(ecorr,ecorr5,ecorr6,n_corr,n_corr1)
 192       endif
 193 C
 194 C If performing constraint dynamics, call the constraint energy
 195 C  after the equilibration time
 196       if(usampl.and.totT.gt.eq_time) then
 197          call EconstrQ
 198          call Econstr_back
 199       else
 200          Uconst=0.0d0
 201          Uconst_back=0.0d0
 202       endif
 203 C
 204 C Sum the energies
 205 C
 206       do i=1,n_ene
 207         energia(i)=0.0d0
 208       enddo
 209       energia(1)=evdw
 210 #ifdef SCP14
 211       energia(2)=evdw2-evdw2_14
 212       energia(18)=evdw2_14
 213 #else
 214       energia(2)=evdw2
 215       energia(18)=0.0d0
 216 #endif
 217 #ifdef SPLITELE
 218       energia(3)=ees
 219       energia(16)=evdw1
 220 #else
 221       energia(3)=ees+evdw1
 222       energia(16)=0.0d0
 223 #endif
 224       energia(4)=ecorr
 225       energia(5)=ecorr5
 226       energia(6)=ecorr6
 227       energia(7)=eel_loc
 228       energia(8)=eello_turn3
 229       energia(9)=eello_turn4
 230       energia(10)=eturn6
 231       energia(20)=Uconst+Uconst_back
 232       energia(22)=evdw_p
 233       energia(23)=evdw_m
 234       call sum_energy(energia,.true.)
 235 c      write (iout,*) "Exit ETOTAL_LONG"
 236       call flush(iout)
 237       return
 238       end
 239 c------------------------------------------------------------------------------
 240       subroutine etotal_short(energia)
 241       implicit real*8 (a-h,o-z)
 242       include 'DIMENSIONS'
 243 c
 244 c Compute the short-range fast-varying contributions to the energy
 245 c
 246 #ifndef ISNAN
 247       external proc_proc
 248 #ifdef WINPGI
 249 cMS$ATTRIBUTES C ::  proc_proc
 250 #endif
 251 #endif
 252 #ifdef MPI
 253       include "mpif.h"
 254       double precision weights_(n_ene)
 255 #endif
 256       include 'COMMON.SETUP'
 257       include 'COMMON.IOUNITS'
 258       double precision energia(0:n_ene)
 259       include 'COMMON.FFIELD'
 260       include 'COMMON.DERIV'
 261       include 'COMMON.INTERACT'
 262       include 'COMMON.SBRIDGE'
 263       include 'COMMON.CHAIN'
 264       include 'COMMON.VAR'
 265       include 'COMMON.LOCAL'
 266       include 'COMMON.CONTROL'
 267 #ifndef DFA
 268       edfadis=0.0d0
 269       edfator=0.0d0
 270       edfanei=0.0d0
 271       edfabet=0.0d0
 272 #endif
 273 c      write(iout,'(a,i2)')'Calling etotal_short ipot=',ipot
 274 c      call flush(iout)
 275       if (modecalc.eq.12.or.modecalc.eq.14) then
 276 #ifdef MPI
 277         if (fg_rank.eq.0) call int_from_cart1(.false.)
 278 #else
 279         call int_from_cart1(.false.)
 280 #endif
 281       endif
 282 #ifdef MPI
 283 c      write(iout,*) "ETOTAL_SHORT Processor",fg_rank,
 284 c     & " absolute rank",myrank," nfgtasks",nfgtasks
 285 c      call flush(iout)
 286       if (nfgtasks.gt.1) then
 287         time00=MPI_Wtime()
 288 C FG slaves call the following matching MPI_Bcast in ERGASTULUM
 289         if (fg_rank.eq.0) then
 290           call MPI_Bcast(2,1,MPI_INTEGER,king,FG_COMM,IERROR)
 291 c          write (iout,*) "Processor",myrank," BROADCAST iorder"
 292 c          call flush(iout)
 293 C FG master sets up the WEIGHTS_ array which will be broadcast to the
 294 C FG slaves as WEIGHTS array.
 295           weights_(1)=wsc
 296           weights_(2)=wscp
 297           weights_(3)=welec
 298           weights_(4)=wcorr
 299           weights_(5)=wcorr5
 300           weights_(6)=wcorr6
 301           weights_(7)=wel_loc
 302           weights_(8)=wturn3
 303           weights_(9)=wturn4
 304           weights_(10)=wturn6
 305           weights_(11)=wang
 306           weights_(12)=wscloc
 307           weights_(13)=wtor
 308           weights_(14)=wtor_d
 309           weights_(15)=wstrain
 310           weights_(16)=wvdwpp
 311           weights_(17)=wbond
 312           weights_(18)=scal14
 313           weights_(21)=wsccor
 314 C FG Master broadcasts the WEIGHTS_ array
 315           call MPI_Bcast(weights_(1),n_ene,
 316      &        MPI_DOUBLE_PRECISION,king,FG_COMM,IERROR)
 317         else
 318 C FG slaves receive the WEIGHTS array
 319           call MPI_Bcast(weights(1),n_ene,
 320      &        MPI_DOUBLE_PRECISION,king,FG_COMM,IERROR)
 321           wsc=weights(1)
 322           wscp=weights(2)
 323           welec=weights(3)
 324           wcorr=weights(4)
 325           wcorr5=weights(5)
 326           wcorr6=weights(6)
 327           wel_loc=weights(7)
 328           wturn3=weights(8)
 329           wturn4=weights(9)
 330           wturn6=weights(10)
 331           wang=weights(11)
 332           wscloc=weights(12)
 333           wtor=weights(13)
 334           wtor_d=weights(14)
 335           wstrain=weights(15)
 336           wvdwpp=weights(16)
 337           wbond=weights(17)
 338           scal14=weights(18)
 339           wsccor=weights(21)
 340         endif
 341 c        write (iout,*),"Processor",myrank," BROADCAST weights"
 342         call MPI_Bcast(c(1,1),maxres6,MPI_DOUBLE_PRECISION,
 343      &    king,FG_COMM,IERR)
 344 c        write (iout,*) "Processor",myrank," BROADCAST c"
 345         call MPI_Bcast(dc(1,1),maxres6,MPI_DOUBLE_PRECISION,
 346      &    king,FG_COMM,IERR)
 347 c        write (iout,*) "Processor",myrank," BROADCAST dc"
 348         call MPI_Bcast(dc_norm(1,1),maxres6,MPI_DOUBLE_PRECISION,
 349      &    king,FG_COMM,IERR)
 350 c        write (iout,*) "Processor",myrank," BROADCAST dc_norm"
 351         call MPI_Bcast(theta(1),nres,MPI_DOUBLE_PRECISION,
 352      &    king,FG_COMM,IERR)
 353 c        write (iout,*) "Processor",myrank," BROADCAST theta"
 354         call MPI_Bcast(phi(1),nres,MPI_DOUBLE_PRECISION,
 355      &    king,FG_COMM,IERR)
 356 c        write (iout,*) "Processor",myrank," BROADCAST phi"
 357         call MPI_Bcast(alph(1),nres,MPI_DOUBLE_PRECISION,
 358      &    king,FG_COMM,IERR)
 359 c        write (iout,*) "Processor",myrank," BROADCAST alph"
 360         call MPI_Bcast(omeg(1),nres,MPI_DOUBLE_PRECISION,
 361      &    king,FG_COMM,IERR)
 362 c        write (iout,*) "Processor",myrank," BROADCAST omeg"
 363         call MPI_Bcast(vbld(1),2*nres,MPI_DOUBLE_PRECISION,
 364      &    king,FG_COMM,IERR)
 365 c        write (iout,*) "Processor",myrank," BROADCAST vbld"
 366         call MPI_Bcast(vbld_inv(1),2*nres,MPI_DOUBLE_PRECISION,
 367      &    king,FG_COMM,IERR)
 368          time_Bcast=time_Bcast+MPI_Wtime()-time00
 369 c        write (iout,*) "Processor",myrank," BROADCAST vbld_inv"
 370       endif
 371 c      write (iout,*) 'Processor',myrank,
 372 c     &  ' calling etotal_short ipot=',ipot
 373 c      call flush(iout)
 374 c      print *,'Processor',myrank,' nnt=',nnt,' nct=',nct
 375 #endif
 376 c      call int_from_cart1(.false.)
 377 C
 378 C Compute the side-chain and electrostatic interaction energy
 379 C
 380       goto (101,102,103,104,105,106) ipot
 381 C Lennard-Jones potential.
 382   101 call elj_short(evdw)
 383 cd    print '(a)','Exit ELJ'
 384       goto 107
 385 C Lennard-Jones-Kihara potential (shifted).
 386   102 call eljk_short(evdw)
 387       goto 107
 388 C Berne-Pechukas potential (dilated LJ, angular dependence).
 389   103 call ebp_short(evdw)
 390       goto 107
 391 C Gay-Berne potential (shifted LJ, angular dependence).
 392   104 call egb_short(evdw,evdw_p,evdw_m)
 393       goto 107
 394 C Gay-Berne-Vorobjev potential (shifted LJ, angular dependence).
 395   105 call egbv_short(evdw)
 396       goto 107
 397 C Soft-sphere potential - already dealt with in the long-range part
 398   106 evdw=0.0d0
 399 c  106 call e_softsphere_short(evdw)
 400 C
 401 C Calculate electrostatic (H-bonding) energy of the main chain.
 402 C
 403   107 continue
 404 #ifdef DFA
 405 C     BARTEK for dfa test!
 406       if (wdfa_dist.gt.0) then
 407          call edfad(edfadis)
 408       else
 409          edfadis=0.0
 410       endif
 411 c      print*, 'edfad is finished!', edfadis
 412       if (wdfa_tor.gt.0) then
 413          call edfat(edfator)
 414       else
 415          edfator=0.0
 416       endif
 417 c      print*, 'edfat is finished!', edfator
 418       if (wdfa_nei.gt.0) then
 419          call edfan(edfanei)
 420       else
 421          edfanei=0.0
 422       endif
 423 c      print*, 'edfan is finished!', edfanei
 424       if (wdfa_beta.gt.0) then
 425          call edfab(edfabet)
 426       else
 427          edfabet=0.0
 428       endif
 429 c      print*, 'edfab is finished!', edfabet
 430 #endif
 431 c
 432 c Calculate the short-range part of Evdwpp
 433 c
 434       call evdwpp_short(evdw1)
 435 c
 436 c Calculate the short-range part of ESCp
 437 c
 438       if (ipot.lt.6) then
 439         call escp_short(evdw2,evdw2_14)
 440       endif
 441 c
 442 c Calculate the bond-stretching energy
 443 c
 444       call ebond(estr)
 445 C
 446 C Calculate the disulfide-bridge and other energy and the contributions
 447 C from other distance constraints.
 448       call edis(ehpb)
 449 C
 450 C Calculate the virtual-bond-angle energy.
 451 C
 452       call ebend(ebe)
 453 C
 454 C Calculate the SC local energy.
 455 C
 456       call vec_and_deriv
 457       call esc(escloc)
 458 C
 459 C Calculate the virtual-bond torsional energy.
 460 C
 461       call etor(etors,edihcnstr)
 462 c
 463 c Homology restraints
 464 c
 465       if (constr_homology.ge.1) then
 466         call e_modeller(ehomology_constr)
 467       else
 468         ehomology_constr=0.0d0
 469       endif
 470 C
 471 C 6/23/01 Calculate double-torsional energy
 472 C
 473       call etor_d(etors_d)
 474 C
 475 C 21/5/07 Calculate local sicdechain correlation energy
 476 C
 477       if (wsccor.gt.0.0d0) then
 478         call eback_sc_corr(esccor)
 479       else
 480         esccor=0.0d0
 481       endif
 482 C
 483 C Put energy components into an array
 484 C
 485       do i=1,n_ene
 486         energia(i)=0.0d0
 487       enddo
 488       energia(1)=evdw
 489 #ifdef SCP14
 490       energia(2)=evdw2-evdw2_14
 491       energia(18)=evdw2_14
 492 #else
 493       energia(2)=evdw2
 494       energia(18)=0.0d0
 495 #endif
 496 #ifdef SPLITELE
 497       energia(16)=evdw1
 498 #else
 499       energia(3)=evdw1
 500 #endif
 501       energia(11)=ebe
 502       energia(12)=escloc
 503       energia(13)=etors
 504       energia(14)=etors_d
 505       energia(15)=ehpb
 506       energia(17)=estr
 507       energia(19)=edihcnstr
 508       energia(21)=esccor
 509       energia(22)=evdw_p
 510       energia(23)=evdw_m
 511       energia(24)=ehomology_constr
 512       energia(25)=edfadis
 513       energia(26)=edfator
 514       energia(27)=edfanei
 515       energia(28)=edfabet
 516 c      write (iout,*) "ETOTAL_SHORT before SUM_ENERGY"
 517       call flush(iout)
 518       call sum_energy(energia,.true.)
 519 c      write (iout,*) "Exit ETOTAL_SHORT"
 520       call flush(iout)
 521       return
 522       end