2 subroutine gradient(n,x,nf,g,uiparm,urparm,ufparm)
5 include 'COMMON.CONTROL'
9 include 'COMMON.INTERACT'
10 include 'COMMON.FFIELD'
12 include 'COMMON.QRESTR'
13 include 'COMMON.IOUNITS'
15 double precision ufparm
18 double precision urparm(1)
19 double precision x(n),g(n)
20 integer i,j,k,ind,ind1
21 double precision f,gthetai,gphii,galphai,gomegai
23 c This subroutine calculates total internal coordinate gradient.
24 c Depending on the number of function evaluations, either whole energy
25 c is evaluated beforehand, Cartesian coordinates and their derivatives in
26 c internal coordinates are reevaluated or only the cartesian-in-internal
27 c coordinate derivatives are evaluated. The subroutine was designed to work
33 cd print *,'grad',nf,icg
34 if (nf-nfl+1) 20,30,40
35 20 call func(n,x,nf,f,uiparm,urparm,ufparm)
36 c write (iout,*) 'grad 20'
39 30 call var_to_geom(n,x)
40 call chainbuild_extconf
41 c write (iout,*) 'grad 30'
43 C Transform the gradient to the gradient in angles.
45 40 call cart2intgrad(n,g)
47 C Add the components corresponding to local energy terms.
50 c Add the usampl contributions
53 gloc(i,icg)=gloc(i,icg)+dugamma(i)
56 gloc(nphi+i,icg)=gloc(nphi+i,icg)+dutheta(i)
60 cd write (iout,*) 'i=',i,'g=',g(i),' gloc=',gloc(i,icg)
63 C Uncomment following three lines for diagnostics.
65 cd call briefout(0,0.0d0)
66 cd write (iout,'(i3,1pe15.5)') (k,g(k),k=1,n)
69 C-------------------------------------------------------------------------
70 subroutine grad_restr(n,x,nf,g,uiparm,urparm,ufparm)
73 include 'COMMON.CHAIN'
74 include 'COMMON.DERIV'
76 include 'COMMON.INTERACT'
77 include 'COMMON.FFIELD'
78 include 'COMMON.IOUNITS'
80 double precision ufparm
83 double precision urparm(1)
84 double precision x(maxvar),g(maxvar),gg(maxvar)
85 integer i,j,k,ig,ind,ij,igall
86 double precision f,gthetai,gphii,galphai,gomegai
89 if (nf-nfl+1) 20,30,40
90 20 call func_restr(n,x,nf,f,uiparm,urparm,ufparm)
91 c write (iout,*) 'grad 20'
96 c Intercept NaNs in the coordinates
97 c write(iout,*) (var(i),i=1,nvar)
102 if (x_sum.ne.x_sum) then
103 write(iout,*)" *** grad_restr : Found NaN in coordinates"
105 print *," *** grad_restr : Found NaN in coordinates"
109 call var_to_geom_restr(n,x)
112 C Evaluate the derivatives of virtual bond lengths and SC vectors in variables.
114 40 call cart2intgrad(n,gg)
116 C Convert the Cartesian gradient into internal-coordinate gradient.
122 IF (mask_phi(i+2).eq.1) THEN
130 IF (mask_theta(i+2).eq.1) THEN
137 if (itype(i).ne.10) then
138 IF (mask_side(i).eq.1) THEN
147 if (itype(i).ne.10) then
148 IF (mask_side(i).eq.1) THEN
150 g(ig)=gg(ialph(i,1)+nside)
156 C Add the components corresponding to local energy terms.
163 if (mask_phi(i).eq.1) then
165 g(ig)=g(ig)+gloc(igall,icg)
171 if (mask_theta(i).eq.1) then
173 g(ig)=g(ig)+gloc(igall,icg)
179 if (itype(i).ne.10) then
181 if (mask_side(i).eq.1) then
183 g(ig)=g(ig)+gloc(igall,icg)
190 cd write (iout,'(a2,i5,a3,f25.8)') 'i=',i,' g=',g(i)
195 C-------------------------------------------------------------------------
202 include 'COMMON.CONTROL'
203 include 'COMMON.CHAIN'
204 include 'COMMON.DERIV'
206 include 'COMMON.INTERACT'
207 include 'COMMON.FFIELD'
209 include 'COMMON.QRESTR'
210 include 'COMMON.IOUNITS'
211 include 'COMMON.TIME1'
214 c This subrouting calculates total Cartesian coordinate gradient.
215 c The subroutine chainbuild_cart and energy MUST be called beforehand.
222 write (iout,*) "Before sum_gradient"
224 write (iout,*) i," gradc ",(gradc(j,i,icg),j=1,3)
225 write (iout,*) i," gradx ",(gradx(j,i,icg),j=1,3)
227 write (iout,*) "gsaxsc, gsaxcx"
229 write (iout,*) i," gsaxsc ",(gsaxsc(j,i),j=1,3)
230 write (iout,*) i," gsaxsx ",(gsaxsx(j,i),j=1,3)
237 write (iout,*) "After sum_gradient"
239 write (iout,*) i," gradc ",(gradc(j,i,icg),j=1,3)
240 write (iout,*) i," gradx ",(gradx(j,i,icg),j=1,3)
243 c If performing constraint dynamics, add the gradients of the constraint energy
244 if(usampl.and.totT.gt.eq_time) then
246 write (iout,*) "dudconst, duscdiff, dugamma,dutheta"
247 write (iout,*) "wumb",wumb
249 write (iout,'(i5,3f10.5,5x,3f10.5,5x,2f10.5)')
250 & i,(dudconst(j,i),j=1,3),(duscdiff(j,i),j=1,3),
251 & dugamma(i),dutheta(i)
256 gradc(j,i,icg)=gradc(j,i,icg)+
257 & wumb*(dudconst(j,i)+duscdiff(j,i))
258 gradx(j,i,icg)=gradx(j,i,icg)+
259 & wumb*(dudxconst(j,i)+duscdiffx(j,i))
263 gloc(i,icg)=gloc(i,icg)+wumb*dugamma(i)
266 gloc(nphi+i,icg)=gloc(nphi+i,icg)+wumb*dutheta(i)
274 time_intcartderiv=time_intcartderiv+MPI_Wtime()-time01
276 cd call checkintcartgrad
277 cd write(iout,*) 'calling int_to_cart'
279 write (iout,*) "gcart, gxcart, gloc before int_to_cart"
283 gcart(j,i)=gradc(j,i,icg)
284 gxcart(j,i)=gradx(j,i,icg)
287 if((itype(i).ne.10).and.(itype(i).ne.ntyp1)) then
288 write (iout,'(i5,2(3f10.5,5x),4f10.5)') i,(gcart(j,i),j=1,3),
289 & (gxcart(j,i),j=1,3),gloc(i,icg),gloc(i+nphi,icg),
290 & gloc(ialph(i,1),icg),gloc(ialph(i,1)+nside,icg)
292 write (iout,'(i5,2(3f10.5,5x),4f10.5)') i,(gcart(j,i),j=1,3),
293 & (gxcart(j,i),j=1,3),gloc(i,icg),gloc(i+nphi,icg)
303 time_inttocart=time_inttocart+MPI_Wtime()-time01
306 write (iout,*) "gcart and gxcart after int_to_cart"
308 write (iout,'(i5,3f10.5,5x,3f10.5)') i,(gcart(j,i),j=1,3),
309 & (gxcart(j,i),j=1,3)
313 time_cartgrad=time_cartgrad+MPI_Wtime()-time00
317 c---------------------------------------------------------------------------
319 subroutine grad_transform
325 include 'COMMON.CONTROL'
326 include 'COMMON.CHAIN'
327 include 'COMMON.DERIV'
329 include 'COMMON.INTERACT'
330 include 'COMMON.FFIELD'
332 include 'COMMON.QRESTR'
333 include 'COMMON.IOUNITS'
334 include 'COMMON.TIME1'
337 write (iout,*)"Converting virtual-bond gradient to CA/SC gradient"
338 write (iout,*) "dC/dX gradient"
340 write (iout,'(i5,3f10.5,5x,3f10.5)') i,(gcart(j,i),j=1,3),
341 & (gxcart(j,i),j=1,3)
346 gcart(j,i)=-gcart(j,i)+gcart(j,i-1)-gxcart(j,i)
347 ! gcart_new(j,i)=-gcart(j,i)+gcart(j,i-1)-gxcart(j,i)
349 ! write (iout,'(i5,3f10.5,5x,3f10.5,5x,3f10.5)') i,(gcart(j,i),j=1,3), &
350 ! (gcart_new(j,i),j=1,3),(gxcart(j,i),j=1,3)
352 ! Correction: dummy residues
354 if (itype(i-1).eq.ntyp1 .and. itype(i).ne.ntyp1) then
355 gcart(:,i)=gcart(:,i)+gcart(:,i-1)
356 else if (itype(i-1).ne.ntyp1 .and. itype(i).eq.ntyp1) then
357 gcart(:,i-1)=gcart(:,i-1)+gcart(:,i)
362 c gcart(j,nnt)=gcart(j,nnt)+gcart(j,1)
365 c if (nct.lt.nres) then
367 c! gcart_new(j,nct)=gcart_new(j,nct)+gcart_new(j,nres)
368 c gcart(j,nct)=gcart(j,nct)+gcart(j,nres)
372 write (iout,*) "CA/SC gradient"
374 write (iout,'(i5,3f10.5,5x,3f10.5)') i,(gcart(j,i),j=1,3),
375 & (gxcart(j,i),j=1,3)
381 C-------------------------------------------------------------------------
385 include 'COMMON.DERIV'
386 include 'COMMON.CHAIN'
389 include 'COMMON.SCCOR'
390 include 'COMMON.SHIELD'
391 integer i,j,kk,intertyp,maxshieldlist
394 C Initialize Cartesian-coordinate gradient
402 gvdwc_scpp(j,i)=0.0d0
404 C below is zero grad for shielding in order: ees (p-p)
405 C ecorr4, eturn3, eturn4, eel_loc, c denotes calfa,x is side-chain
408 gshieldc_loc(j,i)=0.0d0
409 gshieldx_ec(j,i)=0.0d0
410 gshieldc_ec(j,i)=0.0d0
411 gshieldc_loc_ec(j,i)=0.0d0
412 gshieldx_t3(j,i)=0.0d0
413 gshieldc_t3(j,i)=0.0d0
414 gshieldc_loc_t3(j,i)=0.0d0
415 gshieldx_t4(j,i)=0.0d0
416 gshieldc_t4(j,i)=0.0d0
417 gshieldc_loc_t4(j,i)=0.0d0
418 gshieldx_ll(j,i)=0.0d0
419 gshieldc_ll(j,i)=0.0d0
420 gshieldc_loc_ll(j,i)=0.0d0
421 C end of zero grad for shielding
427 gel_loc_long(j,i)=0.0d0
432 gcorr3_turn(j,i)=0.0d0
433 gcorr4_turn(j,i)=0.0d0
435 gradcorr_long(j,i)=0.0d0
436 gradcorr5_long(j,i)=0.0d0
437 gradcorr6_long(j,i)=0.0d0
438 gcorr6_turn_long(j,i)=0.0d0
441 gcorr6_turn(j,i)=0.0d0
451 grad_shield(j,i)=0.0d0
453 gg_tube_sc(j,i)=0.0d0
454 C grad_shield_side is Cbeta sidechain gradient
455 do kk=1,maxshieldlist
456 grad_shield_side(j,kk,i)=0.0d0
457 grad_shield_loc(j,kk,i)=0.0d0
459 C grad_shield_side_ca is Calfa sidechain gradient
462 C grad_shield_side_ca(j,kk,i)=0.0d0
465 gloc_sc(intertyp,i,icg)=0.0d0
480 C Initialize the gradient of local energy terms.
490 gel_loc_turn3(i)=0.0d0
491 gel_loc_turn4(i)=0.0d0
492 gel_loc_turn6(i)=0.0d0
495 c initialize gcart and gxcart
504 c-------------------------------------------------------------------------
505 double precision function fdum()