1 subroutine sum_gradient
4 include "DIMENSIONS.ZSCOPT"
8 cMS$ATTRIBUTES C :: proc_proc
13 include "COMMON.SETUP"
16 double precision gradbufc(3,-1:maxres),gradbufx(3,-1:maxres),
17 & glocbuf(4*maxres),gradbufc_sum(3,-1:maxres)
18 & ,gloc_scbuf(3,-1:maxres)
19 include 'COMMON.IOUNITS'
20 include 'COMMON.FFIELD'
21 include 'COMMON.DERIV'
22 include 'COMMON.INTERACT'
23 include 'COMMON.SBRIDGE'
24 include 'COMMON.CHAIN'
26 include 'COMMON.CONTROL'
27 include 'COMMON.TIME1'
28 include 'COMMON.SCCOR'
30 double precision time00,time_barrier,time_barrier_g,time_reduce
35 write (iout,*) "sum_gradient gvdwc, gvdwx"
37 write (iout,'(i3,3f10.5,5x,3f10.5,5x,f10.5)')
38 & i,(gvdwx(j,i),j=1,3),(gvdwc(j,i),j=1,3)
43 C FG slaves call the following matching MPI_Bcast in ERGASTULUM
44 if (nfgtasks.gt.1 .and. fg_rank.eq.0)
45 & call MPI_Bcast(1,1,MPI_INTEGER,king,FG_COMM,IERROR)
48 C 9/29/08 AL Transform parts of gradients in site coordinates to the gradient
49 C in virtual-bond-vector coordinates
52 c write (iout,*) "gel_loc gel_loc_long and gel_loc_loc"
54 c write (iout,'(i5,3f10.5,2x,3f10.5,2x,f10.5)')
55 c & i,(gel_loc(j,i),j=1,3),(gel_loc_long(j,i),j=1,3),gel_loc_loc(i)
57 c write (iout,*) "gel_loc_tur3 gel_loc_turn4"
59 c write (iout,'(i5,3f10.5,2x,f10.5)')
60 c & i,(gcorr4_turn(j,i),j=1,3),gel_loc_turn4(i)
62 write (iout,*) "gradcorr5 gradcorr5_long gradcorr5_loc"
64 write (iout,'(i3,3f10.5,5x,3f10.5,5x,f10.5)')
65 & i,(gradcorr5(j,i),j=1,3),(gradcorr5_long(j,i),j=1,3),
73 gradbufc(j,i)=wsc*gvdwc(j,i)+
74 & wscp*(gvdwc_scp(j,i)+gvdwc_scpp(j,i))+
75 & welec*gelc_long(j,i)+wvdwpp*gvdwpp(j,i)+
76 & wel_loc*gel_loc_long(j,i)+
77 & wcorr*gradcorr_long(j,i)+
78 & wcorr5*gradcorr5_long(j,i)+
79 & wcorr6*gradcorr6_long(j,i)+
80 & wturn6*gcorr6_turn_long(j,i)+
82 & +wliptran*gliptranc(j,i)
84 & +welec*gshieldc(j,i)
85 & +wcorr*gshieldc_ec(j,i)
86 & +wturn3*gshieldc_t3(j,i)
87 & +wturn4*gshieldc_t4(j,i)
88 & +wel_loc*gshieldc_ll(j,i)
89 c & +wtube*gg_tube(j,i)
98 gradbufc(j,i)=wsc*gvdwc(j,i)+
99 & wscp*(gvdwc_scp(j,i)+gvdwc_scpp(j,i))+
100 & welec*gelc_long(j,i)+
102 & wel_loc*gel_loc_long(j,i)+
103 & wcorr*gradcorr_long(j,i)+
104 & wcorr5*gradcorr5_long(j,i)+
105 & wcorr6*gradcorr6_long(j,i)+
106 & wturn6*gcorr6_turn_long(j,i)+
108 & +wliptran*gliptranc(j,i)
110 & +welec*gshieldc(j,i)
111 & +wcorr*gshieldc_ec(j,i)
112 & +wturn4*gshieldc_t4(j,i)
113 & +wel_loc*gshieldc_ll(j,i)
114 c & +wtube*gg_tube(j,i)
122 if (nfgtasks.gt.1) then
125 write (iout,*) "gradbufc before allreduce"
127 write (iout,'(i3,3f10.5)') i,(gradbufc(j,i),j=1,3)
133 gradbufc_sum(j,i)=gradbufc(j,i)
136 c call MPI_AllReduce(gradbufc(1,1),gradbufc_sum(1,1),3*nres,
137 c & MPI_DOUBLE_PRECISION,MPI_SUM,FG_COMM,IERR)
138 c time_reduce=time_reduce+MPI_Wtime()-time00
140 c write (iout,*) "gradbufc_sum after allreduce"
142 c write (iout,'(i3,3f10.5)') i,(gradbufc_sum(j,i),j=1,3)
147 c time_allreduce=time_allreduce+MPI_Wtime()-time00
155 write (iout,*) "igrad_start",igrad_start," igrad_end",igrad_end
156 write (iout,*) (i," jgrad_start",jgrad_start(i),
157 & " jgrad_end ",jgrad_end(i),
158 & i=igrad_start,igrad_end)
161 c Obsolete and inefficient code; we can make the effort O(n) and, therefore,
162 c do not parallelize this part.
164 c do i=igrad_start,igrad_end
165 c do j=jgrad_start(i),jgrad_end(i)
167 c gradbufc(k,i)=gradbufc(k,i)+gradbufc_sum(k,j)
172 gradbufc(j,nres-1)=gradbufc_sum(j,nres)
176 gradbufc(j,i)=gradbufc(j,i+1)+gradbufc_sum(j,i+1)
180 write (iout,*) "gradbufc after summing"
182 write (iout,'(i3,3f10.5)') i,(gradbufc(j,i),j=1,3)
189 write (iout,*) "gradbufc"
191 write (iout,'(i3,3f10.5)') i,(gradbufc(j,i),j=1,3)
197 gradbufc_sum(j,i)=gradbufc(j,i)
202 gradbufc(j,nres-1)=gradbufc_sum(j,nres)
206 gradbufc(j,i)=gradbufc(j,i+1)+gradbufc_sum(j,i+1)
211 c gradbufc(k,i)=0.0d0
215 c gradbufc(k,i)=gradbufc(k,i)+gradbufc(k,j)
220 write (iout,*) "gradbufc after summing"
222 write (iout,'(i3,3f10.5)') i,(gradbufc(j,i),j=1,3)
230 gradbufc(k,nres)=0.0d0
235 C print *,gradbufc(1,13)
236 C print *,welec*gelc(1,13)
237 C print *,wel_loc*gel_loc(1,13)
238 C print *,0.5d0*(wscp*gvdwc_scpp(1,13))
239 C print *,welec*gelc_long(1,13)+wvdwpp*gvdwpp(1,13)
240 C print *,wel_loc*gel_loc_long(1,13)
241 C print *,gradafm(1,13),"AFM"
242 gradc(j,i,icg)=gradbufc(j,i)+welec*gelc(j,i)+
243 & wel_loc*gel_loc(j,i)+
244 & 0.5d0*(wscp*gvdwc_scpp(j,i)+
245 & welec*gelc_long(j,i)+wvdwpp*gvdwpp(j,i)+
246 & wel_loc*gel_loc_long(j,i)+
247 & wcorr*gradcorr_long(j,i)+
248 & wcorr5*gradcorr5_long(j,i)+
249 & wcorr6*gradcorr6_long(j,i)+
250 & wturn6*gcorr6_turn_long(j,i))+
252 & wcorr*gradcorr(j,i)+
253 & wturn3*gcorr3_turn(j,i)+
254 & wturn4*gcorr4_turn(j,i)+
255 & wcorr5*gradcorr5(j,i)+
256 & wcorr6*gradcorr6(j,i)+
257 & wturn6*gcorr6_turn(j,i)+
258 & wsccor*gsccorc(j,i)
259 & +wscloc*gscloc(j,i)
260 & +wliptran*gliptranc(j,i)
262 & +welec*gshieldc(j,i)
263 & +welec*gshieldc_loc(j,i)
264 & +wcorr*gshieldc_ec(j,i)
265 & +wcorr*gshieldc_loc_ec(j,i)
266 & +wturn3*gshieldc_t3(j,i)
267 & +wturn3*gshieldc_loc_t3(j,i)
268 & +wturn4*gshieldc_t4(j,i)
269 & +wturn4*gshieldc_loc_t4(j,i)
270 & +wel_loc*gshieldc_ll(j,i)
271 & +wel_loc*gshieldc_loc_ll(j,i)
272 c & +wtube*gg_tube(j,i)
275 gradc(j,i,icg)=gradbufc(j,i)+welec*gelc(j,i)+
276 & wel_loc*gel_loc(j,i)+
277 & 0.5d0*(wscp*gvdwc_scpp(j,i)+
278 & welec*gelc_long(j,i)+
279 & wel_loc*gel_loc_long(j,i)+
280 & wcorr*gcorr_long(j,i)+
281 & wcorr5*gradcorr5_long(j,i)+
282 & wcorr6*gradcorr6_long(j,i)+
283 & wturn6*gcorr6_turn_long(j,i))+
285 & wcorr*gradcorr(j,i)+
286 & wturn3*gcorr3_turn(j,i)+
287 & wturn4*gcorr4_turn(j,i)+
288 & wcorr5*gradcorr5(j,i)+
289 & wcorr6*gradcorr6(j,i)+
290 & wturn6*gcorr6_turn(j,i)+
291 & wsccor*gsccorc(j,i)
292 & +wscloc*gscloc(j,i)
293 & +wliptran*gliptranc(j,i)
295 & +welec*gshieldc(j,i)
296 & +welec*gshieldc_loc(j,i)
297 & +wcorr*gshieldc_ec(j,i)
298 & +wcorr*gshieldc_loc_ec(j,i)
299 & +wturn3*gshieldc_t3(j,i)
300 & +wturn3*gshieldc_loc_t3(j,i)
301 & +wturn4*gshieldc_t4(j,i)
302 & +wturn4*gshieldc_loc_t4(j,i)
303 & +wel_loc*gshieldc_ll(j,i)
304 & +wel_loc*gshieldc_loc_ll(j,i)
305 c & +wtube*gg_tube(j,i)
309 gradx(j,i,icg)=wsc*gvdwx(j,i)+wscp*gradx_scp(j,i)+
311 & wstrain*ghpbx(j,i)+wcorr*gradxorr(j,i)+
312 & wsccor*gsccorx(j,i)
313 & +wscloc*gsclocx(j,i)
314 & +wliptran*gliptranx(j,i)
315 & +welec*gshieldx(j,i)
316 & +wcorr*gshieldx_ec(j,i)
317 & +wturn3*gshieldx_t3(j,i)
318 & +wturn4*gshieldx_t4(j,i)
319 & +wel_loc*gshieldx_ll(j,i)
320 c & +wtube*gg_tube_sc(j,i)
327 write (iout,*) "gloc before adding corr"
329 write (iout,*) i,gloc(i,icg)
333 gloc(i,icg)=gloc(i,icg)+wcorr*gcorr_loc(i)
334 & +wcorr5*g_corr5_loc(i)
335 & +wcorr6*g_corr6_loc(i)
336 & +wturn4*gel_loc_turn4(i)
337 & +wturn3*gel_loc_turn3(i)
338 & +wturn6*gel_loc_turn6(i)
339 & +wel_loc*gel_loc_loc(i)
342 write (iout,*) "gloc after adding corr"
344 write (iout,*) i,gloc(i,icg)
348 if (nfgtasks.gt.1) then
351 gradbufc(j,i)=gradc(j,i,icg)
352 gradbufx(j,i)=gradx(j,i,icg)
356 glocbuf(i)=gloc(i,icg)
360 write (iout,*) "gloc_sc before reduce"
363 write (iout,*) i,j,gloc_sc(j,i,icg)
370 gloc_scbuf(j,i)=gloc_sc(j,i,icg)
374 c write (iout,*) "GRAD: MPI_BARRIER"
376 call MPI_Barrier(FG_COMM,IERR)
377 time_barrier_g=time_barrier_g+MPI_Wtime()-time00
379 call MPI_Reduce(gradbufc(1,1),gradc(1,1,icg),3*nres,
380 & MPI_DOUBLE_PRECISION,MPI_SUM,king,FG_COMM,IERR)
381 call MPI_Reduce(gradbufx(1,1),gradx(1,1,icg),3*nres,
382 & MPI_DOUBLE_PRECISION,MPI_SUM,king,FG_COMM,IERR)
383 call MPI_Reduce(glocbuf(1),gloc(1,icg),4*nres,
384 & MPI_DOUBLE_PRECISION,MPI_SUM,king,FG_COMM,IERR)
385 time_reduce=time_reduce+MPI_Wtime()-time00
386 call MPI_Reduce(gloc_scbuf(1,1),gloc_sc(1,1,icg),3*nres,
387 & MPI_DOUBLE_PRECISION,MPI_SUM,king,FG_COMM,IERR)
388 time_reduce=time_reduce+MPI_Wtime()-time00
391 write (iout,*) "gloc_sc after reduce"
394 write (iout,*) i,j,gloc_sc(j,i,icg)
400 write (iout,*) "gloc after reduce"
402 write (iout,*) i,gloc(i,icg)
408 write (iout,*) "gradc gradx gloc"
410 write (iout,'(i5,3f10.5,5x,3f10.5,5x,f10.5)')
411 & i,(gradc(j,i,icg),j=1,3),(gradx(j,i,icg),j=1,3),gloc(i,icg)
415 time_sumgradient=time_sumgradient+MPI_Wtime()-time01
419 c---------------------------------------------------------------------------
420 subroutine sum_gradient_compon
423 include "DIMENSIONS.ZSCOPT"
427 cMS$ATTRIBUTES C :: proc_proc
432 include "COMMON.SETUP"
435 double precision gradbufc(3,-1:maxres),gradbufc_sum(3,-1:maxres),
436 & glocbuf(max_ene,4*maxres)
437 include 'COMMON.IOUNITS'
438 include 'COMMON.FFIELD'
439 include 'COMMON.DERIV'
440 include 'COMMON.INTERACT'
441 include 'COMMON.SBRIDGE'
442 include 'COMMON.CHAIN'
444 include 'COMMON.CONTROL'
445 include 'COMMON.TIME1'
446 include 'COMMON.SCCOR'
447 include "COMMON.NAMES"
448 include "COMMON.ENERGIES"
449 integer i,j,k,ii,iene
450 double precision time_reduce,time00
454 C Sum up the components of the Cartesian gradient.
461 c Part in coordinates
462 gcompon(1,j,i)=gvdwc(j,i)
463 gcompon(2,j,i)=gvdwc_scp(j,i)+gvdwc_scpp(j,i)
464 gcompon(3,j,i)=gelc_long(j,i)
465 gcompon(16,j,i)=gvdwpp(j,i)
466 gcompon(7,j,i)=gel_loc_long(j,i)
467 gcompon(4,j,i)=gradcorr_long(j,i)
468 gcompon(5,j,i)=gradcorr5_long(j,i)
469 gcompon(6,j,i)=gradcorr6_long(j,i)
470 gcompon(10,j,i)=gcorr6_turn_long(j,i)
471 gcompon(15,j,i)=ghpbc(j,i)
472 gcompon(22,j,i)=gliptranc(j,i)
474 gcompon_loc(2,j,i)=0.5d0*gvdwc_scpp(j,i)
475 gcompon_loc(3,j,i)=0.5d0*gelc_long(j,i)+gelc(j,i)
476 gcompon_loc(16,j,i)=0.5d0*gvdwpp(j,i)
477 gcompon_loc(4,j,i)=0.5d0*gradcorr_long(j,i)+gradcorr(j,i)
478 gcompon_loc(5,j,i)=0.5d0*gradcorr5_long(j,i)+gradcorr5(j,i)
479 gcompon_loc(6,j,i)=0.5d0*gradcorr6_long(j,i)+gradcorr6(j,i)
480 gcompon_loc(7,j,i)=0.5d0*gel_loc_long(j,i)+gel_loc(j,i)
481 gcompon_loc(8,j,i)=gcorr3_turn(j,i)
482 gcompon_loc(9,j,i)=gcorr4_turn(j,i)
483 gcompon_loc(10,j,i)=0.5d0*gcorr6_turn_long(j,i)
485 gcompon_loc(12,j,i)=gscloc(j,i)
486 gcompon_loc(17,j,i)=gradb(j,i)
487 gcompon_loc(19,j,i)=gsccorc(j,i)
488 gcompon_loc(22,j,i)=gliptranc(j,i)
489 c sidechain components
490 gcomponx(1,j,i)=gvdwx(j,i)
491 gcomponx(2,j,i)=gradx_scp(j,i)
492 gcomponx(12,j,i)=gsclocx(j,i)
493 gcomponx(17,j,i)=gradbx(j,i)
494 gcomponx(19,j,i)=gsccorx(j,i)
495 gcomponx(22,j,i)=gliptranx(j,i)
496 if (shield_mode.gt.0) then
497 c Part in coordinates
498 gcompon(3,j,i)=gcompon(3,j,i)+gshieldc(j,i)
499 gcompon(4,j,i)=gcompon(4,j,i)+gshieldc_ec(j,i)
500 gcompon(7,j,i)=gcompon(7,j,i)+gshieldc_ll(j,i)
501 gcompon(8,j,i)=gcompon(8,j,i)+gshieldc_t3(j,i)
502 gcompon(9,j,i)=gcompon(9,j,i)+gshieldc_t4(j,i)
504 gcompon_loc(3,j,i)=gshieldc(j,i)+gshieldc_loc(j,i)
505 gcompon_loc(4,j,i)=gcompon_loc(4,j,i)+gshieldc_ec(j,i)
506 & +gshieldc_loc_ec(j,i)
507 gcompon_loc(7,j,i)=gcompon_loc(7,j,i)+gshieldc_ll(j,i)
508 & +gshieldc_loc_ll(j,i)
509 gcompon_loc(8,j,i)=gcompon_loc(8,j,i)+gshieldc_t3(j,i)
510 & +gshieldc_loc_t3(j,i)
511 gcompon_loc(9,j,i)=gcompon_loc(9,j,i)+gshieldc_t4(j,i)
512 & +gshieldc_loc_t4(j,i)
514 gcomponx(3,j,i)=gshieldx(j,i)
515 gcomponx(4,j,i)=gshieldx_ec(j,i)
516 gcomponx(7,j,i)=gshieldx_ec(j,i)
517 gcomponx(8,j,i)=gshieldx_t3(j,i)
518 gcomponx(9,j,i)=gshieldx_t4(j,i)
521 gcompon(3,j,i)=gcompon(3,j,i)+gcompon(16,j,i)
522 gcompon_loc(3,j,i)=gcompon_loc(3,j,i)+gcompon_loc(16,j,i)
528 c gloc(i,icg)=gloc(i,icg)+wcorr*gcorr_loc(i)
529 c & +wcorr5*g_corr5_loc(i)
530 c & +wcorr6*g_corr6_loc(i)
531 c & +wturn4*gel_loc_turn4(i)
532 c & +wturn3*gel_loc_turn3(i)
533 c & +wturn6*gel_loc_turn6(i)
534 c & +wel_loc*gel_loc_loc(i)
535 gloc_compon(4,i)=gcorr_loc(i)
536 gloc_compon(5,i)=g_corr5_loc(i)
537 gloc_compon(6,i)=g_corr6_loc(i)
538 gloc_compon(7,i)=gel_loc_loc(i)
539 gloc_compon(8,i)=gel_loc_turn3(i)
540 gloc_compon(9,i)=gel_loc_turn4(i)
541 gloc_compon(10,i)=gel_loc_turn6(i)
542 c & +wsccor*gsccor_loc(i)
543 c BYLA ROZNICA Z CLUSTER< OSTATNIA LINIA DODANA
546 write (iout,*) "Gradient components right after assignment"
549 write (iout,'(a,i3,1x,a)') "Component",iene,ename(iene)
552 write (iout,'(a4,i5,3e15.5,5x,e15.5,5x,3e15.5,5x,e15.5)')
553 & restyp(itype(i)),i,(gcompon(iene,j,i),j=1,3),
554 & gloc_compon(iene,i),(gcomponx(iene,j,i),j=1,3),
555 & gloc_compon(iene,nres+i)
561 C FG slaves call the following matching MPI_Bcast in ERGASTULUM
562 if (nfgtasks.gt.1 .and. fg_rank.eq.0)
563 & call MPI_Bcast(1,1,MPI_INTEGER,king,FG_COMM,IERROR)
573 gradbufc_sum(j,i)=gcompon(iene,j,i)
577 gradbufc(j,nres-1)=gradbufc_sum(j,nres)
581 gradbufc(j,i)=gradbufc(j,i+1)+gradbufc_sum(j,i+1)
585 write (iout,*) "gradbufc after summing"
587 write (iout,'(i3,3f10.5)') i,(gradbufc(j,i),j=1,3)
592 gradbufc(k,nres)=0.0d0
596 gcompon(iene,j,i)=gradbufc(j,i)+gcompon_loc(iene,j,i)
602 if (nfgtasks.gt.1) then
603 gcompon_loc = gcompon
605 call MPI_Reduce(gcompon_loc(1,1,1),gcompon(1,1,1),
607 & MPI_DOUBLE_PRECISION,MPI_SUM,king,FG_COMM,IERR)
608 gcompon_loc = gcomponx
609 call MPI_Reduce(gcompon_loc(1,1,1),gcomponx(1,1,1),
611 & MPI_DOUBLE_PRECISION,MPI_SUM,king,FG_COMM,IERR)
612 glocbuf = gloc_compon
613 call MPI_Reduce(glocbuf(1,1),gloc_compon(1,1),4*nres*max_ene,
614 & MPI_DOUBLE_PRECISION,MPI_SUM,king,FG_COMM,IERR)
615 time_reduce=time_reduce+MPI_Wtime()-time00
620 write (iout,*) "Gradient components in sum_gradient_compon"
623 write (iout,'(a,i3,1x,a)') "Component",iene,ename(iene)
626 write (iout,'(a4,i5,3e15.5,5x,e15.5,5x,3e15.5,5x,e15.5)')
627 & restyp(itype(i)),i,(gcompon(iene,j,i),j=1,3),
628 & gloc_compon(iene,i),(gcomponx(iene,j,i),j=1,3),
629 & gloc_compon(iene,nres+i)
635 time_sumgradient=time_sumgradient+MPI_Wtime()-time01