source/unres/src-HCD-5D/energy_split-sep.F

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