1 subroutine etotal_long(energia)
5 c Compute the long-range slow-varying contributions to the energy
10 cMS$ATTRIBUTES C :: proc_proc
15 double precision weights_(n_ene)
16 double precision time00
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'
28 include 'COMMON.LOCAL'
29 include 'COMMON.QRESTR'
31 include 'COMMON.CONTROL'
32 include 'COMMON.TIME1'
33 double precision evdw,evdw1,evdw2,evdw2_14,ees,eel_loc,
34 & eello_turn3,eello_turn4,edfadis,estr,ehpb,ebe,ethetacnstr,
35 & escloc,etors,edihcnstr,etors_d,esccor,ecorr,ecorr5,ecorr6,eturn6,
36 & eliptran,Eafmforce,Etube,
37 & esaxs_constr,ehomology_constr,edfator,edfanei,edfabet
38 integer i,n_corr,n_corr1
39 c write(iout,'(a,i2)')'Calling etotal_long ipot=',ipot
41 double precision time01
43 if (modecalc.eq.12.or.modecalc.eq.14) then
45 c if (fg_rank.eq.0) call int_from_cart1(.false.)
47 call int_from_cart1(.false.)
55 ehomology_constr=0.0d0
58 c write(iout,*) "ETOTAL_LONG Processor",fg_rank,
59 c & " absolute rank",myrank," nfgtasks",nfgtasks
61 if (nfgtasks.gt.1) then
63 C FG slaves call the following matching MPI_Bcast in ERGASTULUM
64 if (fg_rank.eq.0) then
65 call MPI_Bcast(3,1,MPI_INTEGER,king,FG_COMM,IERROR)
66 c write (iout,*) "Processor",myrank," BROADCAST iorder"
68 C FG master sets up the WEIGHTS_ array which will be broadcast to the
69 C FG slaves as WEIGHTS array.
92 weights_(28)=wdfa_dist
95 weights_(31)=wdfa_beta
96 C FG Master broadcasts the WEIGHTS_ array
97 call MPI_Bcast(weights_(1),n_ene,
98 & MPI_DOUBLE_PRECISION,king,FG_COMM,IERROR)
100 C FG slaves receive the WEIGHTS array
101 call MPI_Bcast(weights(1),n_ene,
102 & MPI_DOUBLE_PRECISION,king,FG_COMM,IERROR)
125 wdfa_dist=weights(28)
128 wdfa_beta=weights(31)
130 call MPI_Bcast(dc(1,1),6*nres,MPI_DOUBLE_PRECISION,
132 time_Bcast=time_Bcast+MPI_Wtime()-time00
133 time_Bcastw=time_Bcastw+MPI_Wtime()-time00
134 c call chainbuild_cart
135 c call int_from_cart1(.false.)
137 c write (iout,*) 'Processor',myrank,
138 c & ' calling etotal_short ipot=',ipot
140 c print *,'Processor',myrank,' nnt=',nnt,' nct=',nct
141 if (nfgtasks.gt.1) then
142 call MPI_Bcast(itime_mat,1,MPI_INT,king,FG_COMM,IERROR)
144 if (mod(itime_mat,imatupdate).eq.0) then
148 call make_SCp_inter_list_RESPA
149 call make_SCSC_inter_list_RESPA
150 call make_pp_inter_list
151 call make_pp_vdw_inter_list_RESPA
153 time_list=time_list+MPI_Wtime()-time01
157 cd print *,'nnt=',nnt,' nct=',nct
159 C Compute the side-chain and electrostatic interaction energy
164 goto (101,102,103,104,105,106) ipot
165 C Lennard-Jones potential.
166 101 call elj_long(evdw)
167 cd print '(a)','Exit ELJ'
169 C Lennard-Jones-Kihara potential (shifted).
170 102 call eljk_long(evdw)
172 C Berne-Pechukas potential (dilated LJ, angular dependence).
173 103 call ebp_long(evdw)
175 C Gay-Berne potential (shifted LJ, angular dependence).
176 104 call egb_long(evdw)
178 C Gay-Berne-Vorobjev potential (shifted LJ, angular dependence).
179 105 call egbv_long(evdw)
181 C Soft-sphere potential
182 106 call e_softsphere(evdw)
184 C Calculate electrostatic (H-bonding) energy of the main chain.
188 time_evdw_long=time_evdw_long+MPI_Wtime()-time01
195 time_vec=time_vec+MPI_Wtime()-time01
197 c write (iout,*) "etotal_long: shield_mode",shield_mode
201 if (shield_mode.eq.1) then
203 else if (shield_mode.eq.2) then
209 if (welec.gt.0d0.or.wvdwpp.gt.0d0.or.wel_loc.gt.0d0.or.
210 & wturn3.gt.0d0.or.wturn4.gt.0d0 .or. wcorr.gt.0.0d0
211 & .or. wcorr4.gt.0.0d0 .or. wcorr5.gt.0.d0
212 & .or. wcorr6.gt.0.0d0 .or. wturn6.gt.0.0d0 ) then
214 if (welec.gt.0d0.or.wel_loc.gt.0d0.or.
215 & wturn3.gt.0d0.or.wturn4.gt.0d0 .or. wcorr.gt.0.0d0
216 & .or. wcorr4.gt.0.0d0 .or. wcorr5.gt.0.d0
217 & .or. wcorr6.gt.0.0d0 .or. wturn6.gt.0.0d0 ) then
219 call eelec_scale(ees,evdw1,eel_loc,eello_turn3,eello_turn4)
228 c write (iout,*) "Soft-spheer ELEC potential"
229 call eelec_soft_sphere(ees,evdw1,eel_loc,eello_turn3,
233 time_eelec_long=time_eelec_long+MPI_Wtime()-time01
236 C Calculate excluded-volume interaction energy between peptide groups
244 call escp_long(evdw2,evdw2_14)
250 call escp_soft_sphere(evdw2,evdw2_14)
253 time_escp_long=time_escp_long+MPI_Wtime()-time01
257 C 12/1/95 Multi-body terms
261 if ((wcorr4.gt.0.0d0 .or. wcorr5.gt.0.0d0 .or. wcorr6.gt.0.0d0
262 & .or. wturn6.gt.0.0d0) .and. ipot.lt.6) then
263 call multibody_eello(ecorr,ecorr5,ecorr6,eturn6,n_corr,n_corr1)
264 c write (2,*) 'n_corr=',n_corr,' n_corr1=',n_corr1,
265 c &" ecorr",ecorr," ecorr5",ecorr5," ecorr6",ecorr6," eturn6",eturn6
272 if ((wcorr4.eq.0.0d0 .and. wcorr.gt.0.0d0) .and. ipot.lt.6) then
273 call multibody_hb(ecorr,ecorr5,ecorr6,n_corr,n_corr1)
287 C If performing constraint dynamics, call the constraint energy
288 C after the equilibration time
289 if(usampl.and.totT.gt.eq_time) then
292 call Econstr_back_qlike
308 energia(2)=evdw2-evdw2_14
325 energia(8)=eello_turn3
326 energia(9)=eello_turn4
328 energia(20)=Uconst+Uconst_back
329 energia(27)=ehomology_constr
334 call sum_energy(energia,.true.)
335 c write (iout,*) "Exit ETOTAL_LONG"
339 c------------------------------------------------------------------------------
340 subroutine etotal_short(energia)
341 implicit real*8 (a-h,o-z)
344 c Compute the short-range fast-varying contributions to the energy
349 cMS$ATTRIBUTES C :: proc_proc
354 double precision weights_(n_ene)
355 double precision time00
358 include 'COMMON.SETUP'
359 include 'COMMON.IOUNITS'
360 double precision energia(0:n_ene)
361 include 'COMMON.FFIELD'
362 include 'COMMON.DERIV'
363 include 'COMMON.INTERACT'
364 include 'COMMON.SBRIDGE'
365 include 'COMMON.CHAIN'
367 include 'COMMON.LOCAL'
368 include 'COMMON.CONTROL'
369 include 'COMMON.SAXS'
370 include 'COMMON.TORCNSTR'
371 include 'COMMON.TIME1'
372 double precision evdw,evdw1,evdw2,evdw2_14,ees,eel_loc,
373 & eello_turn3,eello_turn4,edfadis,estr,ehpb,ebe,ethetacnstr,
374 & escloc,etors,edihcnstr,etors_d,esccor,ecorr,ecorr5,ecorr6,eturn6,
375 & eliptran,Eafmforce,Etube,
376 & esaxs_constr,ehomology_constr,edfator,edfanei,edfabet
377 integer i,n_corr,n_corr1
379 double precision time01
381 c write(iout,'(a,i2)')'Calling etotal_short ipot=',ipot
383 if (modecalc.eq.12.or.modecalc.eq.14) then
385 if (fg_rank.eq.0) call int_from_cart1(.false.)
387 call int_from_cart1(.false.)
409 c write(iout,*) "ETOTAL_SHORT Processor",fg_rank,
410 c & " absolute rank",myrank," nfgtasks",nfgtasks
412 if (nfgtasks.gt.1) then
414 C FG slaves call the following matching MPI_Bcast in ERGASTULUM
415 if (fg_rank.eq.0) then
416 call MPI_Bcast(2,1,MPI_INTEGER,king,FG_COMM,IERROR)
417 c write (iout,*) "Processor",myrank," BROADCAST iorder"
419 C FG master sets up the WEIGHTS_ array which will be broadcast to the
420 C FG slaves as WEIGHTS array.
441 weights_(29)=wdfa_tor
442 weights_(30)=wdfa_nei
443 weights_(31)=wdfa_beta
444 C FG Master broadcasts the WEIGHTS_ array
445 call MPI_Bcast(weights_(1),n_ene,
446 & MPI_DOUBLE_PRECISION,king,FG_COMM,IERROR)
448 C FG slaves receive the WEIGHTS array
449 call MPI_Bcast(weights(1),n_ene,
450 & MPI_DOUBLE_PRECISION,king,FG_COMM,IERROR)
472 c write (iout,*),"Processor",myrank," BROADCAST weights"
473 call MPI_Bcast(c(1,1),6*nres,MPI_DOUBLE_PRECISION,
475 c write (iout,*) "Processor",myrank," BROADCAST c"
476 call MPI_Bcast(dc(1,1),6*nres,MPI_DOUBLE_PRECISION,
478 c write (iout,*) "Processor",myrank," BROADCAST dc"
479 call MPI_Bcast(dc_norm(1,1),6*nres,MPI_DOUBLE_PRECISION,
481 c write (iout,*) "Processor",myrank," BROADCAST dc_norm"
482 call MPI_Bcast(theta(1),nres,MPI_DOUBLE_PRECISION,
484 c write (iout,*) "Processor",myrank," BROADCAST theta"
485 call MPI_Bcast(phi(1),nres,MPI_DOUBLE_PRECISION,
487 c write (iout,*) "Processor",myrank," BROADCAST phi"
488 call MPI_Bcast(alph(1),nres,MPI_DOUBLE_PRECISION,
490 c write (iout,*) "Processor",myrank," BROADCAST alph"
491 call MPI_Bcast(omeg(1),nres,MPI_DOUBLE_PRECISION,
493 c write (iout,*) "Processor",myrank," BROADCAST omeg"
494 call MPI_Bcast(vbld(1),2*nres,MPI_DOUBLE_PRECISION,
496 c write (iout,*) "Processor",myrank," BROADCAST vbld"
497 call MPI_Bcast(vbld_inv(1),2*nres,MPI_DOUBLE_PRECISION,
499 time_Bcast=time_Bcast+MPI_Wtime()-time00
500 c write (iout,*) "Processor",myrank," BROADCAST vbld_inv"
502 c write (iout,*) 'Processor',myrank,
503 c & ' calling etotal_short ipot=',ipot
505 c print *,'Processor',myrank,' nnt=',nnt,' nct=',nct
507 c call int_from_cart1(.false.)
509 C Compute the side-chain and electrostatic interaction energy
514 goto (101,102,103,104,105,106) ipot
515 C Lennard-Jones potential.
516 101 call elj_short(evdw)
517 cd print '(a)','Exit ELJ'
519 C Lennard-Jones-Kihara potential (shifted).
520 102 call eljk_short(evdw)
522 C Berne-Pechukas potential (dilated LJ, angular dependence).
523 103 call ebp_short(evdw)
525 C Gay-Berne potential (shifted LJ, angular dependence).
526 104 call egb_short(evdw)
528 C Gay-Berne-Vorobjev potential (shifted LJ, angular dependence).
529 105 call egbv_short(evdw)
531 C Soft-sphere potential - already dealt with in the long-range part
533 c 106 call e_softsphere_short(evdw)
535 C Calculate electrostatic (H-bonding) energy of the main chain.
539 time_evdw_short=time_evdw_short+MPI_Wtime()-time01
542 c Calculate the short-range part of Evdwpp
547 call evdwpp_short(evdw1)
549 time_eelec_short=time_eelec_short+MPI_Wtime()-time01
552 c Calculate the short-range part of ESCp
558 call escp_short(evdw2,evdw2_14)
561 time_escp_short=time_escp_short+MPI_Wtime()-time01
564 c Calculate the bond-stretching energy
568 C Calculate the disulfide-bridge and other energy and the contributions
569 C from other distance constraints.
572 C Calculate the virtual-bond-angle energy.
574 if (wang.gt.0d0) then
575 if (tor_mode.eq.0) then
578 C ebend kcc is Kubo cumulant clustered rigorous attemp to derive the
586 if (with_theta_constr) call etheta_constr(ethetacnstr)
588 C Calculate the SC local energy.
595 time_vec=time_vec+MPI_Wtime()-time01
599 C Calculate the virtual-bond torsional energy.
601 if (wtor.gt.0.0d0) then
602 if (tor_mode.eq.0) then
605 C etor kcc is Kubo cumulant clustered rigorous attemp to derive the
614 if (wliptran.gt.0) then
615 call Eliptransfer(eliptran)
619 if (AFMlog.gt.0) then
620 call AFMforce(Eafmforce)
621 else if (selfguide.gt.0) then
622 call AFMvel(Eafmforce)
626 if (TUBElog.eq.1) then
627 C print *,"just before call"
629 elseif (TUBElog.eq.2) then
630 call calctube2(Etube)
635 if (ndih_constr.gt.0) call etor_constr(edihcnstr)
636 c print *,"Processor",myrank," computed Utor"
638 C 6/23/01 Calculate double-torsional energy
640 if ((wtor_d.gt.0.0d0).and.(tor_mode.eq.0)) then
646 c Homology restraints
648 if (constr_homology.ge.1) then
649 call e_modeller(ehomology_constr)
651 ehomology_constr=0.0d0
654 C BARTEK for dfa test!
655 if (wdfa_dist.gt.0) then
660 c print*, 'edfad is finished!', edfadis
661 if (wdfa_tor.gt.0) then
666 c print*, 'edfat is finished!', edfator
667 if (wdfa_nei.gt.0) then
672 c print*, 'edfan is finished!', edfanei
673 if (wdfa_beta.gt.0) then
678 c print*, 'edfab is finished!', edfabet
681 C 21/5/07 Calculate local sicdechain correlation energy
683 if (wsccor.gt.0.0d0) then
684 call eback_sc_corr(esccor)
688 c write (iout,*) "nsaxs",nsaxs," saxs_mode",saxs_mode
689 if (nsaxs.gt.0 .and. saxs_mode.eq.0) then
690 call e_saxs(Esaxs_constr)
691 c write (iout,*) "From Esaxs: Esaxs_constr",Esaxs_constr
692 else if (nsaxs.gt.0 .and. saxs_mode.gt.0) then
693 call e_saxsC(Esaxs_constr)
694 c write (iout,*) "From EsaxsC: Esaxs_constr",Esaxs_constr
699 C Put energy components into an array
706 energia(2)=evdw2-evdw2_14
723 energia(8)=eello_turn3
724 energia(9)=eello_turn4
732 energia(19)=edihcnstr
735 energia(24)=ethetacnstr
737 energia(26)=Esaxs_constr
738 energia(27)=ehomology_constr
743 c write (iout,*) "ETOTAL_SHORT before SUM_ENERGY"
745 call sum_energy(energia,.true.)
746 c write (iout,*) "Exit ETOTAL_SHORT"