! include 'COMMON.VECTORS'
! include 'COMMON.FFIELD'
! include 'COMMON.TIME1'
- real(kind=8),dimension(3) :: ggg,gggp,gggm,erij,dcosb,dcosg
+ real(kind=8),dimension(3) :: ggg,gggp,gggm,erij,dcosb,dcosg,xtemp
real(kind=8),dimension(3,3) :: erder,uryg,urzg,vryg,vrzg
real(kind=8),dimension(2,2) :: acipa !el,a_temp
!el real(kind=8),dimension(3,4) :: agg,aggi,aggi1,aggj,aggj1
rij=xj*xj+yj*yj+zj*zj
rrmij=1.0D0/rij
rij=dsqrt(rij)
+!C print *,xmedi,ymedi,zmedi,xj,yj,zj,boxxsize,rij
sss_ele_cut=sscale_ele(rij)
sss_ele_grad=sscagrad_ele(rij)
+! sss_ele_cut=1.0d0
+! sss_ele_grad=0.0d0
! print *,sss_ele_cut,sss_ele_grad,&
! (rij),r_cut_ele,rlamb_ele
- if (sss_ele_cut.le.0.0) go to 128
+! if (sss_ele_cut.le.0.0) go to 128
rmij=1.0D0/rij
r3ij=rrmij*rmij
eel_loc_ij=a22*muij(1)+a23*muij(2)+a32*muij(3) &
+a33*muij(4)
! write (iout,*) 'i',i,' j',j,' eel_loc_ij',eel_loc_ij
-
+! eel_loc_ij=0.0
if (energy_dec) write (iout,'(a6,2i5,0pf7.3)') &
'eelloc',i,j,eel_loc_ij
! if (energy_dec) write (iout,*) "a22",a22," a23",a23," a32",a32," a33",a33
! if (energy_dec) write (iout,*) "muij",muij
! write (iout,*) a22,muij(1),a23,muij(2),a32,muij(3)
-
+
eel_loc=eel_loc+eel_loc_ij*sss_ele_cut
! Partial derivatives in virtual-bond dihedral angles gamma
if (i.gt.1) &
*sss_ele_cut
! Derivatives of eello in DC(i+1) thru DC(j-1) or DC(nres-2)
! do l=1,3
- ggg(1)=(agg(1,1)*muij(1)+ &
- agg(1,2)*muij(2)+agg(1,3)*muij(3)+agg(1,4)*muij(4)) &
- *sss_ele_cut+eel_loc_ij*sss_ele_grad*rmij*xj
- ggg(2)=(agg(2,1)*muij(1)+ &
- agg(2,2)*muij(2)+agg(2,3)*muij(3)+agg(2,4)*muij(4)) &
- *sss_ele_cut+eel_loc_ij*sss_ele_grad*rmij*yj
- ggg(3)=(agg(3,1)*muij(1)+ &
- agg(3,2)*muij(2)+agg(3,3)*muij(3)+agg(3,4)*muij(4)) &
- *sss_ele_cut+eel_loc_ij*sss_ele_grad*rmij*zj
+! ggg(1)=(agg(1,1)*muij(1)+ &
+! agg(1,2)*muij(2)+agg(1,3)*muij(3)+agg(1,4)*muij(4)) &
+! *sss_ele_cut &
+! +eel_loc_ij*sss_ele_grad*rmij*xj
+! ggg(2)=(agg(2,1)*muij(1)+ &
+! agg(2,2)*muij(2)+agg(2,3)*muij(3)+agg(2,4)*muij(4)) &
+! *sss_ele_cut &
+! +eel_loc_ij*sss_ele_grad*rmij*yj
+! ggg(3)=(agg(3,1)*muij(1)+ &
+! agg(3,2)*muij(2)+agg(3,3)*muij(3)+agg(3,4)*muij(4)) &
+! *sss_ele_cut &
+! +eel_loc_ij*sss_ele_grad*rmij*zj
+ xtemp(1)=xj
+ xtemp(2)=yj
+ xtemp(3)=zj
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))&
+ *sss_ele_cut &
+ +eel_loc_ij*sss_ele_grad*rmij*xtemp(l)
+
gel_loc_long(l,j)=gel_loc_long(l,j)+ggg(l)
gel_loc_long(l,i)=gel_loc_long(l,i)-ggg(l)
!grad ghalf=0.5d0*ggg(l)
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))&
*sss_ele_cut
+!+eel_loc_ij*sss_ele_grad*rmij*xtemp(l)
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))&
+ aggi1(l,2)*muij(2)+aggi1(l,3)*muij(3) &
+ +aggi1(l,4)*muij(4))&
*sss_ele_cut
+!+eel_loc_ij*sss_ele_grad*rmij*xtemp(l)
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))&
*sss_ele_cut
+!+eel_loc_ij*sss_ele_grad*rmij*xtemp(l)
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))&
+ aggj1(l,2)*muij(2)+aggj1(l,3)*muij(3) &
+ +aggj1(l,4)*muij(4))&
*sss_ele_cut
+!+eel_loc_ij*sss_ele_grad*rmij*xtemp(l)
enddo
ENDIF
! Change 12/26/95 to calculate four-body contributions to H-bonding energy
! include 'COMMON.CONTROL'
real(kind=8),dimension(3) :: ggg
!el local variables
- integer :: i,iint,j,k,iteli,itypj
+ integer :: i,iint,j,k,iteli,itypj,subchap
real(kind=8) :: evdw2,evdw2_14,xi,yi,zi,xj,yj,zj,rrij,fac,&
- e1,e2,evdwij
+ e1,e2,evdwij,rij
+ real(kind=8) :: xj_safe,yj_safe,zj_safe,xj_temp,yj_temp,zj_temp,&
+ dist_temp, dist_init
+ integer xshift,yshift,zshift
evdw2=0.0D0
evdw2_14=0.0d0
xi=0.5D0*(c(1,i)+c(1,i+1))
yi=0.5D0*(c(2,i)+c(2,i+1))
zi=0.5D0*(c(3,i)+c(3,i+1))
+ xi=mod(xi,boxxsize)
+ if (xi.lt.0) xi=xi+boxxsize
+ yi=mod(yi,boxysize)
+ if (yi.lt.0) yi=yi+boxysize
+ zi=mod(zi,boxzsize)
+ if (zi.lt.0) zi=zi+boxzsize
do iint=1,nscp_gr(i)
! yj=c(2,nres+j)-yi
! zj=c(3,nres+j)-zi
! Uncomment following three lines for Ca-p interactions
- xj=c(1,j)-xi
- yj=c(2,j)-yi
- zj=c(3,j)-zi
+! xj=c(1,j)-xi
+! yj=c(2,j)-yi
+! zj=c(3,j)-zi
+ xj=c(1,j)
+ yj=c(2,j)
+ zj=c(3,j)
+ xj=mod(xj,boxxsize)
+ if (xj.lt.0) xj=xj+boxxsize
+ yj=mod(yj,boxysize)
+ if (yj.lt.0) yj=yj+boxysize
+ zj=mod(zj,boxzsize)
+ if (zj.lt.0) zj=zj+boxzsize
+ dist_init=(xj-xi)**2+(yj-yi)**2+(zj-zi)**2
+ xj_safe=xj
+ yj_safe=yj
+ zj_safe=zj
+ subchap=0
+ do xshift=-1,1
+ do yshift=-1,1
+ do zshift=-1,1
+ xj=xj_safe+xshift*boxxsize
+ yj=yj_safe+yshift*boxysize
+ zj=zj_safe+zshift*boxzsize
+ dist_temp=(xj-xi)**2+(yj-yi)**2+(zj-zi)**2
+ if(dist_temp.lt.dist_init) then
+ dist_init=dist_temp
+ xj_temp=xj
+ yj_temp=yj
+ zj_temp=zj
+ subchap=1
+ endif
+ enddo
+ enddo
+ enddo
+ if (subchap.eq.1) then
+ xj=xj_temp-xi
+ yj=yj_temp-yi
+ zj=zj_temp-zi
+ else
+ xj=xj_safe-xi
+ yj=yj_safe-yi
+ zj=zj_safe-zi
+ endif
+
rrij=1.0D0/(xj*xj+yj*yj+zj*zj)
+ rij=dsqrt(1.0d0/rrij)
+ sss_ele_cut=sscale_ele(rij)
+ sss_ele_grad=sscagrad_ele(rij)
+! print *,sss_ele_cut,sss_ele_grad,&
+! (rij),r_cut_ele,rlamb_ele
+ if (sss_ele_cut.le.0.0) cycle
fac=rrij**expon2
e1=fac*fac*aad(itypj,iteli)
e2=fac*bad(itypj,iteli)
if (iabs(j-i) .le. 2) then
e1=scal14*e1
e2=scal14*e2
- evdw2_14=evdw2_14+e1+e2
+ evdw2_14=evdw2_14+(e1+e2)*sss_ele_cut
endif
evdwij=e1+e2
- evdw2=evdw2+evdwij
+ evdw2=evdw2+evdwij*sss_ele_cut
! if (energy_dec) write (iout,'(a6,2i5,0pf7.3,2i3,3e11.3)') &
! 'evdw2',i,j,evdwij,iteli,itypj,fac,aad(itypj,iteli),&
if (energy_dec) write (iout,'(a6,2i5,0pf7.3)') &
!
! Calculate contributions to the gradient in the virtual-bond and SC vectors.
!
- fac=-(evdwij+e1)*rrij
+ fac=-(evdwij+e1)*rrij*sss_ele_cut
+ fac=fac+evdwij*sss_ele_grad/rij/expon
ggg(1)=xj*fac
ggg(2)=yj*fac
ggg(3)=zj*fac
! Check the gradient of the virtual-bond and SC vectors in the internal
! coordinates.
!
- aincr=1.0d-7
- aincr2=5.0d-8
+ aincr=1.0d-6
+ aincr2=5.0d-7
call cartder
write (iout,'(a)') '**************** dx/dalpha'
write (iout,'(a)')
nf=0
nfl=0
call zerograd
- aincr=1.0D-7
- print '(a)','CG processor',me,' calling CHECK_CART.'
+ aincr=1.0D-5
+ print '(a)','CG processor',me,' calling CHECK_CART.',aincr
nf=0
icall=0
call geom_to_var(nvar,x)
! call intcartderiv
! call checkintcartgrad
call zerograd
- aincr=1.0D-6
- write(iout,*) 'Calling CHECK_ECARTINT.'
+ aincr=2.0D-5
+ write(iout,*) 'Calling CHECK_ECARTINT.',aincr
nf=0
icall=0
call geom_to_var(nvar,x)
! Obtaining the gamma derivatives from sine derivative
if (phi(i).gt.-pi4.and.phi(i).le.pi4.or. &
phi(i).gt.pi34.and.phi(i).le.pi.or. &
- phi(i).gt.-pi.and.phi(i).le.-pi34) then
+ phi(i).ge.-pi.and.phi(i).le.-pi34) then
call vecpr(dc_norm(1,i-1),dc_norm(1,i-2),vp1)
call vecpr(dc_norm(1,i-3),dc_norm(1,i-1),vp2)
call vecpr(dc_norm(1,i-3),dc_norm(1,i-2),vp3)
projects / unres4.git / commitdiff