2 ***********************************************************************
3 * This subroutine calculates the derivatives of the consecutive virtual
4 * bond vectors and the SC vectors in the virtual-bond angles theta and
5 * virtual-torsional angles phi, as well as the derivatives of SC vectors
6 * in the angles alpha and omega, describing the location of a side chain
7 * in its local coordinate system.
9 * The derivatives are stored in the following arrays:
11 * DDCDV - the derivatives of virtual-bond vectors DC in theta and phi.
12 * The structure is as follows:
14 * dDC(x,2)/dT(3),...,dDC(z,2)/dT(3),0, 0, 0
15 * dDC(x,3)/dT(4),...,dDC(z,3)/dT(4),dDC(x,3)/dP(4),dDC(y,4)/dP(4),dDC(z,4)/dP(4)
16 * . . . . . . . . . . . . . . . . . .
17 * dDC(x,N-1)/dT(4),...,dDC(z,N-1)/dT(4),dDC(x,N-1)/dP(4),dDC(y,N-1)/dP(4),dDC(z,N-1)/dP(4)
21 * dDC(x,N-1)/dT(N),...,dDC(z,N-1)/dT(N),dDC(x,N-1)/dP(N),dDC(y,N-1)/dP(N),dDC(z,N-1)/dP(N)
23 * DXDV - the derivatives of the side-chain vectors in theta and phi.
24 * The structure is same as above.
26 * DCDS - the derivatives of the side chain vectors in the local spherical
27 * andgles alph and omega:
29 * dX(x,2)/dA(2),dX(y,2)/dA(2),dX(z,2)/dA(2),dX(x,2)/dO(2),dX(y,2)/dO(2),dX(z,2)/dO(2)
30 * dX(x,3)/dA(3),dX(y,3)/dA(3),dX(z,3)/dA(3),dX(x,3)/dO(3),dX(y,3)/dO(3),dX(z,3)/dO(3)
34 * dX(x,N-1)/dA(N-1),dX(y,N-1)/dA(N-1),dX(z,N-1)/dA(N-1),dX(x,N-1)/dO(N-1),dX(y,N-1)/dO(N-1),dX(z,N-1)/dO(N-1)
36 * Version of March '95, based on an early version of November '91.
38 ***********************************************************************
41 include 'COMMON.IOUNITS'
43 include 'COMMON.CHAIN'
44 include 'COMMON.DERIV'
46 include 'COMMON.LOCAL'
47 include 'COMMON.INTERACT'
48 double precision drt(3,3,maxres),rdt(3,3,maxres),dp(3,3),
49 &temp(3,3),prordt(3,3,maxres),prodrt(3,3,maxres)
50 double precision xx(3),xx1(3),alphi,omegi,xj,dpjk,yp,xp,xxp,yyp
51 double precision cosalphi,sinalphi,cosomegi,sinomegi,theta2,
52 & cost2,sint2,rj,dxoiij,tempkl,dxoijk,dsci,zzp,dj,dpkl
54 double precision fromto(3,3)
56 double precision fromto(3,3,maxdim)
57 c common /przechowalnia/ fromto
59 integer i,ii,j,jjj,k,l,m,indi,ind,ind1
60 * get the position of the jth ijth fragment of the chain coordinate system
61 * in the fromto array.
63 indmat(i,j)=((2*(nres-2)-i)*(i-1))/2+j-1
64 call chainbuild_extconf
68 * calculate the derivatives of transformation matrix elements in theta
97 * 3/10/2020 Adam: The fromto array to be created only for smaller
98 * systems; for large ones its elements to be calculated on the fly.
100 * generate the matrix products of type r(i)t(i)...r(j)t(j)
104 write(iout,*) i,i+1,ind
112 fromto(k,l,ind)=temp(k,l)
115 c write(iout,'(7hfromto 9f10.5)')((fromto(k,l,ind),l=1,3),k=1,3)
118 write(iout,*) i,j+1,ind
119 c write(iout,'(7htemp 9f10.5)')((temp(k,l),l=1,3),k=1,3)
120 c write(iout,'(7hrt 9f10.5)')((rt(k,l,j),l=1,3),k=1,3)
125 dpkl=dpkl+temp(k,m)*rt(m,l,j)
131 c write(iout,'(7hfromto 9f10.5)')((fromto(k,l,ind),l=1,3),k=1,3)
141 * Calculate derivatives.
147 * Derivatives of DC(i+1) in theta(i+2)
149 c write (iout,*) "theta i",i
150 c write(iout,'(7hprod 9f10.5)')((prod(k,l,i),l=1,3),k=1,3)
151 c write(iout,'(7hrdt 9f10.5)')((rdt(k,l,i),l=1,3),k=1,3)
152 c write(iout,*) "vbld",vbld(i+2)
157 dpjk=dpjk+prod(j,l,i)*rdt(l,k,i)
160 prordt(j,k,i)=dp(j,k)
163 c dcdv(j,ind1)=vbld(i+1)*dp(j,1)
164 dcdv(j,ind1)=vbld(i+2)*dp(j,1)
166 c write(iout,'(7hdcdv 3f10.5)')(dcdv(k,ind1),k=1,3)
168 * Derivatives of SC(i+1) in theta(i+2)
170 xx1(1)=-0.5D0*xloc(2,i+1)
171 xx1(2)= 0.5D0*xloc(1,i+1)
175 xj=xj+r(j,k,i)*xx1(k)
182 rj=rj+prod(j,k,i)*xx(k)
186 c write (iout,*) "dxdv",(dxdv(j,ind1),j=1,3)
188 * Derivatives of SC(i+1) in theta(i+3). The have to be handled differently
189 * than the other off-diagonal derivatives.
194 dxoiij=dxoiij+dp(j,k)*xrot(k,i+2)
196 dxdv(j,ind1+1)=dxoiij
198 c write(iout,*)ind1+1,(dxdv(j,ind1+1),j=1,3)
200 * Derivatives of DC(i+1) in phi(i+2)
206 dpjk=dpjk+prod(j,l,i)*drt(l,k,i)
209 prodrt(j,k,i)=dp(j,k)
211 c dcdv(j+3,ind1)=vbld(i+1)*dp(j,1)
212 dcdv(j+3,ind1)=vbld(i+2)*dp(j,1)
215 * Derivatives of SC(i+1) in phi(i+2)
218 xx(3)= xloc(2,i+1)*r(2,2,i)+xloc(3,i+1)*r(2,3,i)
219 xx(2)=-xloc(2,i+1)*r(3,2,i)-xloc(3,i+1)*r(3,3,i)
223 rj=rj+prod(j,k,i)*xx(k)
228 * Derivatives of SC(i+1) in phi(i+3).
233 dxoiij=dxoiij+dp(j,k)*xrot(k,i+2)
235 dxdv(j+3,ind1+1)=dxoiij
238 * Calculate the derivatives of DC(i+1) and SC(i+1) in theta(i+3) thru
239 * theta(nres) and phi(i+3) thru phi(nres).
245 c write(iout,*)'i=',i,' j=',j,' ind=',ind,' ind1=',ind1
246 call build_fromto(i+1,j+1,fromto)
247 c write(iout,'(7hfromto 9f10.5)')((fromto(k,l),l=1,3),k=1,3)
252 tempkl=tempkl+prordt(k,m,i)*fromto(m,l)
262 tempkl=tempkl+prordt(k,m,i)*fromto(m,l,ind)
268 c write(iout,'(7hfromto 9f10.5)')((fromto(k,l,ind),l=1,3),k=1,3)
269 c write(iout,'(7hprod 9f10.5)')((prod(k,l,i),l=1,3),k=1,3)
270 c write(iout,'(7htemp 9f10.5)')((temp(k,l),l=1,3),k=1,3)
271 * Derivatives of virtual-bond vectors in theta
273 c dcdv(k,ind1)=vbld(i+1)*temp(k,1)
274 dcdv(k,ind1)=vbld(j+2)*temp(k,1)
276 c write(iout,'(7hdcdv 3f10.5)')(dcdv(k,ind1),k=1,3)
277 * Derivatives of SC vectors in theta
281 dxoijk=dxoijk+temp(k,l)*xrot(l,j+2)
283 dxdv(k,ind1+1)=dxoijk
285 c write(iout,'(7htheta 3f10.5)')(dxdv(k,ind1),k=1,3)
287 *--- Calculate the derivatives in phi
294 tempkl=tempkl+prodrt(k,m,i)*fromto(m,l)
304 tempkl=tempkl+prodrt(k,m,i)*fromto(m,l,ind)
311 c dcdv(k+3,ind1)=vbld(i+1)*temp(k,1)
312 dcdv(k+3,ind1)=vbld(j+2)*temp(k,1)
317 dxoijk=dxoijk+temp(k,l)*xrot(l,j+2)
319 dxdv(k+3,ind1+1)=dxoijk
325 write (iout,'(a)') '****************** ddc/dtheta'
330 write (iout,'(2i4,3e14.6)') i,j,(dcdv(k,ii),k=1,3)
334 write (iout,'(a)') '******************* ddc/dphi'
339 write (iout,'(2i4,3e14.6)') i,j,(dcdv(k+3,ii),k=1,3)
344 write (iout,'(a)') '**************** dx/dtheta'
349 write (iout,'(2i4,3e14.6)') i,j,(dxdv(k,ii),k=1,3)
353 write (iout,'(a)') '***************** dx/dphi'
358 write (iout,'(2i4,3e14.6)') i,j,(dxdv(k+3,ii),k=1,3)
364 * Derivatives in alpha and omega:
372 if(alphi.ne.alphi) alphi=100.0
373 if(omegi.ne.omegi) omegi=-100.0
378 cd print *,'i=',i,' dsci=',dsci,' alphi=',alphi,' omegi=',omegi
383 temp(1,1)=-dsci*sinalphi
384 temp(2,1)= dsci*cosalphi*cosomegi
385 temp(3,1)=-dsci*cosalphi*sinomegi
387 temp(2,2)=-dsci*sinalphi*sinomegi
388 temp(3,2)=-dsci*sinalphi*cosomegi
389 theta2=pi-0.5D0*theta(i+1)
393 cd print *,((temp(l,k),l=1,3),k=1,2)
397 xxp= xp*cost2+yp*sint2
398 yyp=-xp*sint2+yp*cost2
401 xx(2)=yyp*r(2,2,i-1)+zzp*r(2,3,i-1)
402 xx(3)=yyp*r(3,2,i-1)+zzp*r(3,3,i-1)
406 dj=dj+prod(k,l,i-1)*xx(l)
416 subroutine build_fromto(i,j,fromto)
419 include 'COMMON.CHAIN'
420 include 'COMMON.IOUNITS'
422 double precision fromto(3,3),temp(3,3),dp(3,3)
423 double precision dpkl
426 * generate the matrix products of type r(i)t(i)...r(j)t(j) on the fly
428 c write (iout,*) "temp on entry"
429 c write (iout,'(3f10.5)') ((temp(k,l),l=1,3),k=1,3)
440 fromto(k,l)=temp(k,l)
450 dpkl=dpkl+temp(k,m)*rt(m,l,j-1)
462 c write (iout,*) "temp upon exit"
463 c write (iout,'(3f10.5)') ((temp(k,l),l=1,3),k=1,3)