update new files
[unres.git] / source / unres / src_MD-M-homology / 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.CONTROL'
19       include 'COMMON.IOUNITS'
20       double precision energia(0:n_ene)
21       include 'COMMON.FFIELD'
22       include 'COMMON.DERIV'
23       include 'COMMON.INTERACT'
24       include 'COMMON.SBRIDGE'
25       include 'COMMON.CHAIN'
26       include 'COMMON.VAR'
27       include 'COMMON.LOCAL'
28       include 'COMMON.MD'
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_(25)=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(25)
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       if (ipot.lt.6) then
139 #ifdef SPLITELE
140          if (welec.gt.0d0.or.wvdwpp.gt.0d0.or.wel_loc.gt.0d0.or.
141      &       wturn3.gt.0d0.or.wturn4.gt.0d0 .or. wcorr.gt.0.0d0
142      &       .or. wcorr4.gt.0.0d0 .or. wcorr5.gt.0.d0
143      &       .or. wcorr6.gt.0.0d0 .or. wturn6.gt.0.0d0 ) then
144 #else
145          if (welec.gt.0d0.or.wel_loc.gt.0d0.or.
146      &       wturn3.gt.0d0.or.wturn4.gt.0d0 .or. wcorr.gt.0.0d0
147      &       .or. wcorr4.gt.0.0d0 .or. wcorr5.gt.0.d0
148      &       .or. wcorr6.gt.0.0d0 .or. wturn6.gt.0.0d0 ) then
149 #endif
150            call eelec_scale(ees,evdw1,eel_loc,eello_turn3,eello_turn4)
151          else
152             ees=0
153             evdw1=0
154             eel_loc=0
155             eello_turn3=0
156             eello_turn4=0
157          endif
158       else
159 c        write (iout,*) "Soft-spheer ELEC potential"
160         call eelec_soft_sphere(ees,evdw1,eel_loc,eello_turn3,
161      &   eello_turn4)
162       endif
163 C
164 C Calculate excluded-volume interaction energy between peptide groups
165 C and side chains.
166 C
167       if (ipot.lt.6) then
168        if(wscp.gt.0d0) then
169         call escp_long(evdw2,evdw2_14)
170        else
171         evdw2=0
172         evdw2_14=0
173        endif
174       else
175         call escp_soft_sphere(evdw2,evdw2_14)
176       endif
177
178 C 12/1/95 Multi-body terms
179 C
180       n_corr=0
181       n_corr1=0
182       if ((wcorr4.gt.0.0d0 .or. wcorr5.gt.0.0d0 .or. wcorr6.gt.0.0d0 
183      &    .or. wturn6.gt.0.0d0) .and. ipot.lt.6) then
184          call multibody_eello(ecorr,ecorr5,ecorr6,eturn6,n_corr,n_corr1)
185 c         write (2,*) 'n_corr=',n_corr,' n_corr1=',n_corr1,
186 c     &" ecorr",ecorr," ecorr5",ecorr5," ecorr6",ecorr6," eturn6",eturn6
187       else
188          ecorr=0.0d0
189          ecorr5=0.0d0
190          ecorr6=0.0d0
191          eturn6=0.0d0
192       endif
193       if ((wcorr4.eq.0.0d0 .and. wcorr.gt.0.0d0) .and. ipot.lt.6) then
194          call multibody_hb(ecorr,ecorr5,ecorr6,n_corr,n_corr1)
195       endif
196
197 C If performing constraint dynamics, call the constraint energy
198 C  after the equilibration time
199       if(usampl.and.totT.gt.eq_time) then
200          call EconstrQ   
201          call Econstr_back
202       else
203          Uconst=0.0d0
204          Uconst_back=0.0d0
205       endif
206
207 C Sum the energies
208 C
209       do i=1,n_ene
210         energia(i)=0.0d0
211       enddo
212       energia(1)=evdw
213 #ifdef SCP14
214       energia(2)=evdw2-evdw2_14
215       energia(18)=evdw2_14
216 #else
217       energia(2)=evdw2
218       energia(18)=0.0d0
219 #endif
220 #ifdef SPLITELE
221       energia(3)=ees
222       energia(16)=evdw1
223 #else
224       energia(3)=ees+evdw1
225       energia(16)=0.0d0
226 #endif
227       energia(4)=ecorr
228       energia(5)=ecorr5
229       energia(6)=ecorr6
230       energia(7)=eel_loc
231       energia(8)=eello_turn3
232       energia(9)=eello_turn4
233       energia(10)=eturn6
234       energia(20)=Uconst+Uconst_back
235       call sum_energy(energia,.true.)
236 c      write (iout,*) "Exit ETOTAL_LONG"
237       call flush(iout)
238       return
239       end
240 c------------------------------------------------------------------------------
241       subroutine etotal_short(energia)
242       implicit real*8 (a-h,o-z)
243       include 'DIMENSIONS'
244 c
245 c Compute the short-range fast-varying contributions to the energy
246 c
247 #ifndef ISNAN
248       external proc_proc
249 #ifdef WINPGI
250 cMS$ATTRIBUTES C ::  proc_proc
251 #endif
252 #endif
253 #ifdef MPI
254       include "mpif.h"
255       double precision weights_(n_ene)
256 #endif
257       include 'COMMON.SETUP'
258       include 'COMMON.CONTROL'
259       include 'COMMON.IOUNITS'
260       double precision energia(0:n_ene)
261       include 'COMMON.FFIELD'
262       include 'COMMON.DERIV'
263       include 'COMMON.INTERACT'
264       include 'COMMON.SBRIDGE'
265       include 'COMMON.CHAIN'
266       include 'COMMON.VAR'
267       include 'COMMON.LOCAL'
268
269 c      write(iout,'(a,i2)')'Calling etotal_short ipot=',ipot
270 c      call flush(iout)
271       if (modecalc.eq.12.or.modecalc.eq.14) then
272 #ifdef MPI
273         if (fg_rank.eq.0) call int_from_cart1(.false.)
274 #else
275         call int_from_cart1(.false.)
276 #endif
277       endif
278 #ifdef MPI      
279 c      write(iout,*) "ETOTAL_SHORT Processor",fg_rank,
280 c     & " absolute rank",myrank," nfgtasks",nfgtasks
281 c      call flush(iout)
282       if (nfgtasks.gt.1) then
283         time00=MPI_Wtime()
284 C FG slaves call the following matching MPI_Bcast in ERGASTULUM
285         if (fg_rank.eq.0) then
286           call MPI_Bcast(2,1,MPI_INTEGER,king,FG_COMM,IERROR)
287 c          write (iout,*) "Processor",myrank," BROADCAST iorder"
288 c          call flush(iout)
289 C FG master sets up the WEIGHTS_ array which will be broadcast to the 
290 C FG slaves as WEIGHTS array.
291           weights_(1)=wsc
292           weights_(2)=wscp
293           weights_(3)=welec
294           weights_(4)=wcorr
295           weights_(5)=wcorr5
296           weights_(6)=wcorr6
297           weights_(7)=wel_loc
298           weights_(8)=wturn3
299           weights_(9)=wturn4
300           weights_(10)=wturn6
301           weights_(11)=wang
302           weights_(12)=wscloc
303           weights_(13)=wtor
304           weights_(14)=wtor_d
305           weights_(15)=wstrain
306           weights_(16)=wvdwpp
307           weights_(17)=wbond
308           weights_(18)=scal14
309           weights_(21)=wsccor
310           weights_(25)=wsaxs
311 C FG Master broadcasts the WEIGHTS_ array
312           call MPI_Bcast(weights_(1),n_ene,
313      &        MPI_DOUBLE_PRECISION,king,FG_COMM,IERROR)
314         else
315 C FG slaves receive the WEIGHTS array
316           call MPI_Bcast(weights(1),n_ene,
317      &        MPI_DOUBLE_PRECISION,king,FG_COMM,IERROR)
318           wsc=weights(1)
319           wscp=weights(2)
320           welec=weights(3)
321           wcorr=weights(4)
322           wcorr5=weights(5)
323           wcorr6=weights(6)
324           wel_loc=weights(7)
325           wturn3=weights(8)
326           wturn4=weights(9)
327           wturn6=weights(10)
328           wang=weights(11)
329           wscloc=weights(12)
330           wtor=weights(13)
331           wtor_d=weights(14)
332           wstrain=weights(15)
333           wvdwpp=weights(16)
334           wbond=weights(17)
335           scal14=weights(18)
336           wsccor=weights(21)
337           wsaxs=weights(25)
338         endif
339 c        write (iout,*),"Processor",myrank," BROADCAST weights"
340         call MPI_Bcast(c(1,1),maxres6,MPI_DOUBLE_PRECISION,
341      &    king,FG_COMM,IERR)
342 c        write (iout,*) "Processor",myrank," BROADCAST c"
343         call MPI_Bcast(dc(1,1),maxres6,MPI_DOUBLE_PRECISION,
344      &    king,FG_COMM,IERR)
345 c        write (iout,*) "Processor",myrank," BROADCAST dc"
346         call MPI_Bcast(dc_norm(1,1),maxres6,MPI_DOUBLE_PRECISION,
347      &    king,FG_COMM,IERR)
348 c        write (iout,*) "Processor",myrank," BROADCAST dc_norm"
349         call MPI_Bcast(theta(1),nres,MPI_DOUBLE_PRECISION,
350      &    king,FG_COMM,IERR)
351 c        write (iout,*) "Processor",myrank," BROADCAST theta"
352         call MPI_Bcast(phi(1),nres,MPI_DOUBLE_PRECISION,
353      &    king,FG_COMM,IERR)
354 c        write (iout,*) "Processor",myrank," BROADCAST phi"
355         call MPI_Bcast(alph(1),nres,MPI_DOUBLE_PRECISION,
356      &    king,FG_COMM,IERR)
357 c        write (iout,*) "Processor",myrank," BROADCAST alph"
358         call MPI_Bcast(omeg(1),nres,MPI_DOUBLE_PRECISION,
359      &    king,FG_COMM,IERR)
360 c        write (iout,*) "Processor",myrank," BROADCAST omeg"
361         call MPI_Bcast(vbld(1),2*nres,MPI_DOUBLE_PRECISION,
362      &    king,FG_COMM,IERR)
363 c        write (iout,*) "Processor",myrank," BROADCAST vbld"
364         call MPI_Bcast(vbld_inv(1),2*nres,MPI_DOUBLE_PRECISION,
365      &    king,FG_COMM,IERR)
366          time_Bcast=time_Bcast+MPI_Wtime()-time00
367 c        write (iout,*) "Processor",myrank," BROADCAST vbld_inv"
368       endif
369 c      write (iout,*) 'Processor',myrank,
370 c     &  ' calling etotal_short ipot=',ipot
371 c      call flush(iout)
372 c      print *,'Processor',myrank,' nnt=',nnt,' nct=',nct
373 #endif     
374 c      call int_from_cart1(.false.)
375 C
376 C Compute the side-chain and electrostatic interaction energy
377 C
378       goto (101,102,103,104,105,106) ipot
379 C Lennard-Jones potential.
380   101 call elj_short(evdw)
381 cd    print '(a)','Exit ELJ'
382       goto 107
383 C Lennard-Jones-Kihara potential (shifted).
384   102 call eljk_short(evdw)
385       goto 107
386 C Berne-Pechukas potential (dilated LJ, angular dependence).
387   103 call ebp_short(evdw)
388       goto 107
389 C Gay-Berne potential (shifted LJ, angular dependence).
390   104 call egb_short(evdw)
391       goto 107
392 C Gay-Berne-Vorobjev potential (shifted LJ, angular dependence).
393   105 call egbv_short(evdw)
394       goto 107
395 C Soft-sphere potential - already dealt with in the long-range part
396   106 evdw=0.0d0
397 c  106 call e_softsphere_short(evdw)
398 C
399 C Calculate electrostatic (H-bonding) energy of the main chain.
400 C
401   107 continue
402 c
403 c Calculate the short-range part of Evdwpp
404 c
405       call evdwpp_short(evdw1)
406 c
407 c Calculate the short-range part of ESCp
408 c
409       if (ipot.lt.6) then
410         call escp_short(evdw2,evdw2_14)
411       endif
412 c
413 c Calculate the bond-stretching energy
414 c
415       call ebond(estr)
416
417 C Calculate the disulfide-bridge and other energy and the contributions
418 C from other distance constraints.
419       call edis(ehpb)
420 C
421 C Calculate the virtual-bond-angle energy.
422 C
423       call ebend(ebe)
424 C
425 C Calculate the SC local energy.
426 C
427       call vec_and_deriv
428       call esc(escloc)
429 C
430 C Calculate the virtual-bond torsional energy.
431 C
432       call etor(etors,edihcnstr)
433 C
434 C 6/23/01 Calculate double-torsional energy
435 C
436       call etor_d(etors_d)
437 C
438 C 21/5/07 Calculate local sicdechain correlation energy
439 C
440       if (wsccor.gt.0.0d0) then
441         call eback_sc_corr(esccor)
442       else
443         esccor=0.0d0
444       endif
445
446 c      write (iout,*) "nsaxs",nsaxs," saxs_mode",saxs_mode
447       if (nsaxs.gt.0 .and. saxs_mode.eq.0) then
448         call e_saxs(Esaxs_constr)
449 c        write (iout,*) "From Esaxs: Esaxs_constr",Esaxs_constr
450       else if (nsaxs.gt.0 .and. saxs_mode.gt.0) then
451         call e_saxsC(Esaxs_constr)
452 c        write (iout,*) "From EsaxsC: Esaxs_constr",Esaxs_constr
453       else
454         Esaxs_constr = 0.0d0
455       endif        
456 C
457 C Put energy components into an array
458 C
459       do i=1,n_ene
460         energia(i)=0.0d0
461       enddo
462       energia(1)=evdw
463 #ifdef SCP14
464       energia(2)=evdw2-evdw2_14
465       energia(18)=evdw2_14
466 #else
467       energia(2)=evdw2
468       energia(18)=0.0d0
469 #endif
470 #ifdef SPLITELE
471       energia(16)=evdw1
472 #else
473       energia(3)=evdw1
474 #endif
475       energia(11)=ebe
476       energia(12)=escloc
477       energia(13)=etors
478       energia(14)=etors_d
479       energia(15)=ehpb
480       energia(17)=estr
481       energia(19)=edihcnstr
482       energia(21)=esccor
483       energia(25)=Esaxs_constr
484 c      write (iout,*) "ETOTAL_SHORT before SUM_ENERGY"
485       call flush(iout)
486       call sum_energy(energia,.true.)
487 c      write (iout,*) "Exit ETOTAL_SHORT"
488       call flush(iout)
489       return
490       end