C Calculate the components of the gradient in DC and X
C
do l=1,3
- gvdwc(l,i)=gvdwc(l,i)-gg(l)+gg_lipi(k)
- gvdwc(l,j)=gvdwc(l,j)+gg(l)+gg_lipj(k)
+ gvdwc(l,i)=gvdwc(l,i)-gg(l)+gg_lipi(l)
+ gvdwc(l,j)=gvdwc(l,j)+gg(l)+gg_lipj(l)
enddo
return
end
include 'COMMON.VECTORS'
include 'COMMON.FFIELD'
include 'COMMON.TIME1'
+ include 'COMMON.SHIELD'
dimension ggg(3),gggp(3),gggm(3),erij(3),dcosb(3),dcosg(3),
& erder(3,3),uryg(3,3),urzg(3,3),vryg(3,3),vrzg(3,3)
double precision acipa(2,2),agg(3,4),aggi(3,4),aggi1(3,4),
include 'COMMON.VECTORS'
include 'COMMON.FFIELD'
include 'COMMON.TIME1'
+ include 'COMMON.SHIELD'
dimension ggg(3),gggp(3),gggm(3),erij(3),dcosb(3),dcosg(3),
& erder(3,3),uryg(3,3),urzg(3,3),vryg(3,3),vrzg(3,3)
double precision acipa(2,2),agg(3,4),aggi(3,4),aggi1(3,4),
fac3=ael6i*r6ij
fac4=ael3i*r3ij
evdwij=ev1+ev2
+ if (shield_mode.eq.0) then
+ fac_shield(i)=1.0
+ fac_shield(j)=1.0
+ endif
el1=fac3*(4.0D0+fac*fac-3.0D0*(cosb*cosb+cosg*cosg))
el2=fac4*fac
+ el1=el1*fac_shield(i)**2*fac_shield(j)**2
+ el2=el2*fac_shield(i)**2*fac_shield(j)**2
eesij=el1+el2
C 12/26/95 - for the evaluation of multi-body H-bonding interactions
ees0ij=4.0D0+fac*fac-3.0D0*(cosb*cosb+cosg*cosg)
ggg(1)=facel*xj
ggg(2)=facel*yj
ggg(3)=facel*zj
+ if ((fac_shield(i).gt.0).and.(fac_shield(j).gt.0).and.
+ & (shield_mode.gt.0)) then
+C print *,i,j
+ do ilist=1,ishield_list(i)
+ iresshield=shield_list(ilist,i)
+ do k=1,3
+ rlocshield=grad_shield_side(k,ilist,i)*eesij/fac_shield(i)
+ & *2.0
+ gshieldx(k,iresshield)=gshieldx(k,iresshield)+
+ & rlocshield
+ & +grad_shield_loc(k,ilist,i)*eesij/fac_shield(i)*2.0
+ gshieldc(k,iresshield-1)=gshieldc(k,iresshield-1)+rlocshield
+C gshieldc_loc(k,iresshield)=gshieldc_loc(k,iresshield)
+C & +grad_shield_loc(k,ilist,i)*eesij/fac_shield(i)
+C if (iresshield.gt.i) then
+C do ishi=i+1,iresshield-1
+C gshieldc(k,ishi)=gshieldc(k,ishi)+rlocshield
+C & +grad_shield_loc(k,ilist,i)*eesij/fac_shield(i)
+C
+C enddo
+C else
+C do ishi=iresshield,i
+C gshieldc(k,ishi)=gshieldc(k,ishi)-rlocshield
+C & -grad_shield_loc(k,ilist,i)*eesij/fac_shield(i)
+C
+C enddo
+C endif
+ enddo
+ enddo
+ do ilist=1,ishield_list(j)
+ iresshield=shield_list(ilist,j)
+ do k=1,3
+ rlocshield=grad_shield_side(k,ilist,j)*eesij/fac_shield(j)
+ & *2.0
+ gshieldx(k,iresshield)=gshieldx(k,iresshield)+
+ & rlocshield
+ & +grad_shield_loc(k,ilist,j)*eesij/fac_shield(j)*2.0
+ gshieldc(k,iresshield-1)=gshieldc(k,iresshield-1)+rlocshield
+ enddo
+ enddo
+
+ do k=1,3
+ gshieldc(k,i)=gshieldc(k,i)+
+ & grad_shield(k,i)*eesij/fac_shield(i)*2.0
+ gshieldc(k,j)=gshieldc(k,j)+
+ & grad_shield(k,j)*eesij/fac_shield(j)*2.0
+ gshieldc(k,i-1)=gshieldc(k,i-1)+
+ & grad_shield(k,i)*eesij/fac_shield(i)*2.0
+ gshieldc(k,j-1)=gshieldc(k,j-1)+
+ & grad_shield(k,j)*eesij/fac_shield(j)*2.0
+
+ enddo
+ endif
+
c do k=1,3
c ghalf=0.5D0*ggg(k)
c gelc(k,i)=gelc(k,i)+ghalf
cd print '(2i3,2(3(1pd14.5),3x))',i,j,(dcosb(k),k=1,3),
cd & (dcosg(k),k=1,3)
do k=1,3
- ggg(k)=ecosb*dcosb(k)+ecosg*dcosg(k)
+ ggg(k)=(ecosb*dcosb(k)+ecosg*dcosg(k))*
+ & fac_shield(i)**2*fac_shield(j)**2
+
enddo
c do k=1,3
c ghalf=0.5D0*ggg(k)
cgrad enddo
do k=1,3
gelc(k,i)=gelc(k,i)
- & +(ecosa*(dc_norm(k,j)-cosa*dc_norm(k,i))
- & + ecosb*(erij(k)-cosb*dc_norm(k,i)))*vbld_inv(i+1)
+ & +((ecosa*(dc_norm(k,j)-cosa*dc_norm(k,i))
+ & + ecosb*(erij(k)-cosb*dc_norm(k,i)))*vbld_inv(i+1))
+ & *fac_shield(i)**2*fac_shield(j)**2
+
gelc(k,j)=gelc(k,j)
- & +(ecosa*(dc_norm(k,i)-cosa*dc_norm(k,j))
- & + ecosg*(erij(k)-cosg*dc_norm(k,j)))*vbld_inv(j+1)
+ & +((ecosa*(dc_norm(k,i)-cosa*dc_norm(k,j))
+ & + ecosg*(erij(k)-cosg*dc_norm(k,j)))*vbld_inv(j+1))
+ & *fac_shield(i)**2*fac_shield(j)**2
gelc_long(k,j)=gelc_long(k,j)+ggg(k)
gelc_long(k,i)=gelc_long(k,i)-ggg(k)
enddo
if (energy_dec) write (iout,'(a6,2i5,0pf7.3)')
& 'eelloc',i,j,eel_loc_ij
+
+ if (shield_mode.eq.0) then
+ fac_shield(i)=1.0
+ fac_shield(j)=1.0
+C else
+C fac_shield(i)=0.4
+C fac_shield(j)=0.6
+ endif
+ eel_loc_ij=eel_loc_ij
+ & *fac_shield(i)*fac_shield(j)
eel_loc=eel_loc+eel_loc_ij
+
+ if ((fac_shield(i).gt.0).and.(fac_shield(j).gt.0).and.
+ & (shield_mode.gt.0)) then
+C print *,i,j
+
+ do ilist=1,ishield_list(i)
+ iresshield=shield_list(ilist,i)
+ do k=1,3
+ rlocshield=grad_shield_side(k,ilist,i)*eel_loc_ij
+ & /fac_shield(i)
+C & *2.0
+ gshieldx_ll(k,iresshield)=gshieldx_ll(k,iresshield)+
+ & rlocshield
+ & +grad_shield_loc(k,ilist,i)*eel_loc_ij/fac_shield(i)
+ gshieldc_ll(k,iresshield-1)=gshieldc_ll(k,iresshield-1)
+ & +rlocshield
+ enddo
+ enddo
+ do ilist=1,ishield_list(j)
+ iresshield=shield_list(ilist,j)
+ do k=1,3
+ rlocshield=grad_shield_side(k,ilist,j)*eel_loc_ij
+ & /fac_shield(j)
+C & *2.0
+ gshieldx_ll(k,iresshield)=gshieldx_ll(k,iresshield)+
+ & rlocshield
+ & +grad_shield_loc(k,ilist,j)*eel_loc_ij/fac_shield(j)
+ gshieldc_ll(k,iresshield-1)=gshieldc_ll(k,iresshield-1)
+ & +rlocshield
+
+ enddo
+ enddo
+
+ do k=1,3
+ gshieldc_ll(k,i)=gshieldc_ll(k,i)+
+ & grad_shield(k,i)*eel_loc_ij/fac_shield(i)
+ gshieldc_ll(k,j)=gshieldc_ll(k,j)+
+ & grad_shield(k,j)*eel_loc_ij/fac_shield(j)
+ gshieldc_ll(k,i-1)=gshieldc_ll(k,i-1)+
+ & grad_shield(k,i)*eel_loc_ij/fac_shield(i)
+ gshieldc_ll(k,j-1)=gshieldc_ll(k,j-1)+
+ & grad_shield(k,j)*eel_loc_ij/fac_shield(j)
+ enddo
+ endif
+
C Partial derivatives in virtual-bond dihedral angles gamma
if (i.gt.1)
& gel_loc_loc(i-1)=gel_loc_loc(i-1)+
- & a22*muder(1,i)*mu(1,j)+a23*muder(1,i)*mu(2,j)
- & +a32*muder(2,i)*mu(1,j)+a33*muder(2,i)*mu(2,j)
+ & (a22*muder(1,i)*mu(1,j)+a23*muder(1,i)*mu(2,j)
+ & +a32*muder(2,i)*mu(1,j)+a33*muder(2,i)*mu(2,j))
+ & *fac_shield(i)*fac_shield(j)
+
gel_loc_loc(j-1)=gel_loc_loc(j-1)+
- & a22*mu(1,i)*muder(1,j)+a23*mu(1,i)*muder(2,j)
- & +a32*mu(2,i)*muder(1,j)+a33*mu(2,i)*muder(2,j)
+ & (a22*mu(1,i)*muder(1,j)+a23*mu(1,i)*muder(2,j)
+ & +a32*mu(2,i)*muder(1,j)+a33*mu(2,i)*muder(2,j))
+ & *fac_shield(i)*fac_shield(j)
+
C Derivatives of eello in DC(i+1) thru DC(j-1) or DC(nres-2)
do l=1,3
- ggg(l)=agg(l,1)*muij(1)+
- & agg(l,2)*muij(2)+agg(l,3)*muij(3)+agg(l,4)*muij(4)
+ ggg(l)=(agg(l,1)*muij(1)+
+ & agg(l,2)*muij(2)+agg(l,3)*muij(3)+agg(l,4)*muij(4))
+ & *fac_shield(i)*fac_shield(j)
+
gel_loc_long(l,j)=gel_loc_long(l,j)+ggg(l)
gel_loc_long(l,i)=gel_loc_long(l,i)-ggg(l)
cgrad ghalf=0.5d0*ggg(l)
cgrad enddo
C Remaining derivatives of eello
do l=1,3
- gel_loc(l,i)=gel_loc(l,i)+aggi(l,1)*muij(1)+
- & aggi(l,2)*muij(2)+aggi(l,3)*muij(3)+aggi(l,4)*muij(4)
- gel_loc(l,i+1)=gel_loc(l,i+1)+aggi1(l,1)*muij(1)+
- & aggi1(l,2)*muij(2)+aggi1(l,3)*muij(3)+aggi1(l,4)*muij(4)
- gel_loc(l,j)=gel_loc(l,j)+aggj(l,1)*muij(1)+
- & aggj(l,2)*muij(2)+aggj(l,3)*muij(3)+aggj(l,4)*muij(4)
- gel_loc(l,j1)=gel_loc(l,j1)+aggj1(l,1)*muij(1)+
- & aggj1(l,2)*muij(2)+aggj1(l,3)*muij(3)+aggj1(l,4)*muij(4)
+ gel_loc(l,i)=gel_loc(l,i)+(aggi(l,1)*muij(1)+
+ & aggi(l,2)*muij(2)+aggi(l,3)*muij(3)+aggi(l,4)*muij(4))
+ & *fac_shield(i)*fac_shield(j)
+
+ gel_loc(l,i+1)=gel_loc(l,i+1)+(aggi1(l,1)*muij(1)+
+ & aggi1(l,2)*muij(2)+aggi1(l,3)*muij(3)+aggi1(l,4)*muij(4))
+ & *fac_shield(i)*fac_shield(j)
+
+ gel_loc(l,j)=gel_loc(l,j)+(aggj(l,1)*muij(1)+
+ & aggj(l,2)*muij(2)+aggj(l,3)*muij(3)+aggj(l,4)*muij(4))
+ & *fac_shield(i)*fac_shield(j)
+
+ gel_loc(l,j1)=gel_loc(l,j1)+(aggj1(l,1)*muij(1)+
+ & aggj1(l,2)*muij(2)+aggj1(l,3)*muij(3)+aggj1(l,4)*muij(4))
+ & *fac_shield(i)*fac_shield(j)
+
enddo
ENDIF
C Change 12/26/95 to calculate four-body contributions to H-bonding energy
ees0mij=0
endif
c ees0mij=0.0D0
+ if (shield_mode.eq.0) then
+ fac_shield(i)=1.0d0
+ fac_shield(j)=1.0d0
+ else
+ ees0plist(num_conti,i)=j
+C fac_shield(i)=0.4d0
+C fac_shield(j)=0.6d0
+ endif
ees0p(num_conti,i)=0.5D0*fac3*(ees0pij+ees0mij)
+ & *fac_shield(i)*fac_shield(j)
ees0m(num_conti,i)=0.5D0*fac3*(ees0pij-ees0mij)
+ & *fac_shield(i)*fac_shield(j)
+
C Diagnostics. Comment out or remove after debugging!
c ees0p(num_conti,i)=0.5D0*fac3*ees0pij
c ees0m(num_conti,i)=0.5D0*fac3*ees0mij
gacontp_hb1(k,num_conti,i)=!ghalfp
& +(ecosap*(dc_norm(k,j)-cosa*dc_norm(k,i))
& + ecosbp*(erij(k)-cosb*dc_norm(k,i)))*vbld_inv(i+1)
+ & *fac_shield(i)*fac_shield(j)
+
gacontp_hb2(k,num_conti,i)=!ghalfp
& +(ecosap*(dc_norm(k,i)-cosa*dc_norm(k,j))
& + ecosgp*(erij(k)-cosg*dc_norm(k,j)))*vbld_inv(j+1)
+ & *fac_shield(i)*fac_shield(j)
+
gacontp_hb3(k,num_conti,i)=gggp(k)
+ & *fac_shield(i)*fac_shield(j)
+
gacontm_hb1(k,num_conti,i)=!ghalfm
& +(ecosam*(dc_norm(k,j)-cosa*dc_norm(k,i))
& + ecosbm*(erij(k)-cosb*dc_norm(k,i)))*vbld_inv(i+1)
+ & *fac_shield(i)*fac_shield(j)
+
gacontm_hb2(k,num_conti,i)=!ghalfm
& +(ecosam*(dc_norm(k,i)-cosa*dc_norm(k,j))
& + ecosgm*(erij(k)-cosg*dc_norm(k,j)))*vbld_inv(j+1)
+ & *fac_shield(i)*fac_shield(j)
+
gacontm_hb3(k,num_conti,i)=gggm(k)
+ & *fac_shield(i)*fac_shield(j)
+
enddo
ENDIF ! wcorr
endif ! num_conti.le.maxconts
projects / unres.git / blobdiff