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