- include 'COMMON.VECTORS'
- dimension uyder(3,3,2),uzder(3,3,2)
-C Compute the local reference systems. For reference system (i), the
-C X-axis points from CA(i) to CA(i+1), the Y axis is in the
-C CA(i)-CA(i+1)-CA(i+2) plane, and the Z axis is perpendicular to this plane.
- do i=1,nres-1
- if (i.eq.nres-1) then
-C Case of the last full residue
-C Compute the Z-axis
- call vecpr(dc_norm(1,i),dc_norm(1,i-1),uz(1,i))
- costh=dcos(pi-theta(nres))
- fac=1.0d0/dsqrt(1.0d0-costh*costh)
-c write (iout,*) 'fac',fac,
-c & 1.0d0/dsqrt(scalar(uz(1,i),uz(1,i)))
- fac=1.0d0/dsqrt(scalar(uz(1,i),uz(1,i)))
- do k=1,3
- uz(k,i)=fac*uz(k,i)
- enddo
-C Compute the derivatives of uz
- uzder(1,1,1)= 0.0d0
- uzder(2,1,1)=-dc_norm(3,i-1)
- uzder(3,1,1)= dc_norm(2,i-1)
- uzder(1,2,1)= dc_norm(3,i-1)
- uzder(2,2,1)= 0.0d0
- uzder(3,2,1)=-dc_norm(1,i-1)
- uzder(1,3,1)=-dc_norm(2,i-1)
- uzder(2,3,1)= dc_norm(1,i-1)
- uzder(3,3,1)= 0.0d0
- uzder(1,1,2)= 0.0d0
- uzder(2,1,2)= dc_norm(3,i)
- uzder(3,1,2)=-dc_norm(2,i)
- uzder(1,2,2)=-dc_norm(3,i)
- uzder(2,2,2)= 0.0d0
- uzder(3,2,2)= dc_norm(1,i)
- uzder(1,3,2)= dc_norm(2,i)
- uzder(2,3,2)=-dc_norm(1,i)
- uzder(3,3,2)= 0.0d0
-C Compute the Y-axis
- do k=1,3
- uy(k,i)=fac*(dc_norm(k,i-1)-costh*dc_norm(k,i))
- enddo
- facy=fac
- facy=1.0d0/dsqrt(scalar(dc_norm(1,i),dc_norm(1,i))*
- & (scalar(dc_norm(1,i-1),dc_norm(1,i-1))**2-
- & scalar(dc_norm(1,i),dc_norm(1,i-1))**2))
- do k=1,3
-c uy(k,i)=facy*(dc_norm(k,i+1)-costh*dc_norm(k,i))
- uy(k,i)=
-c & facy*(
- & dc_norm(k,i-1)*scalar(dc_norm(1,i),dc_norm(1,i))
- & -scalar(dc_norm(1,i),dc_norm(1,i-1))*dc_norm(k,i)
-c & )
- enddo
-c write (iout,*) 'facy',facy,
-c & 1.0d0/dsqrt(scalar(uy(1,i),uy(1,i)))
- facy=1.0d0/dsqrt(scalar(uy(1,i),uy(1,i)))
- do k=1,3
- uy(k,i)=facy*uy(k,i)
- enddo
-C Compute the derivatives of uy
- do j=1,3
- do k=1,3
- uyder(k,j,1)=2*dc_norm(k,i-1)*dc_norm(j,i)
- & -dc_norm(k,i)*dc_norm(j,i-1)
- uyder(k,j,2)=-dc_norm(j,i)*dc_norm(k,i)
- enddo
-c uyder(j,j,1)=uyder(j,j,1)-costh
-c uyder(j,j,2)=1.0d0+uyder(j,j,2)
- uyder(j,j,1)=uyder(j,j,1)
- & -scalar(dc_norm(1,i),dc_norm(1,i-1))
- uyder(j,j,2)=scalar(dc_norm(1,i),dc_norm(1,i))
- & +uyder(j,j,2)
- enddo
- do j=1,2
- do k=1,3
- do l=1,3
- uygrad(l,k,j,i)=uyder(l,k,j)
- uzgrad(l,k,j,i)=uzder(l,k,j)
- enddo
- enddo
- enddo
- call unormderiv(uy(1,i),uyder(1,1,1),facy,uygrad(1,1,1,i))
- call unormderiv(uy(1,i),uyder(1,1,2),facy,uygrad(1,1,2,i))
- call unormderiv(uz(1,i),uzder(1,1,1),fac,uzgrad(1,1,1,i))
- call unormderiv(uz(1,i),uzder(1,1,2),fac,uzgrad(1,1,2,i))
- else
-C Other residues
-C Compute the Z-axis
- call vecpr(dc_norm(1,i),dc_norm(1,i+1),uz(1,i))
- costh=dcos(pi-theta(i+2))
- fac=1.0d0/dsqrt(1.0d0-costh*costh)
- fac=1.0d0/dsqrt(scalar(uz(1,i),uz(1,i)))
- do k=1,3
- uz(k,i)=fac*uz(k,i)
- enddo
-C Compute the derivatives of uz
- uzder(1,1,1)= 0.0d0
- uzder(2,1,1)=-dc_norm(3,i+1)
- uzder(3,1,1)= dc_norm(2,i+1)
- uzder(1,2,1)= dc_norm(3,i+1)
- uzder(2,2,1)= 0.0d0
- uzder(3,2,1)=-dc_norm(1,i+1)
- uzder(1,3,1)=-dc_norm(2,i+1)
- uzder(2,3,1)= dc_norm(1,i+1)
- uzder(3,3,1)= 0.0d0
- uzder(1,1,2)= 0.0d0
- uzder(2,1,2)= dc_norm(3,i)
- uzder(3,1,2)=-dc_norm(2,i)
- uzder(1,2,2)=-dc_norm(3,i)
- uzder(2,2,2)= 0.0d0
- uzder(3,2,2)= dc_norm(1,i)
- uzder(1,3,2)= dc_norm(2,i)
- uzder(2,3,2)=-dc_norm(1,i)
- uzder(3,3,2)= 0.0d0
-C Compute the Y-axis
- facy=fac
- facy=1.0d0/dsqrt(scalar(dc_norm(1,i),dc_norm(1,i))*
- & (scalar(dc_norm(1,i+1),dc_norm(1,i+1))**2-
- & scalar(dc_norm(1,i),dc_norm(1,i+1))**2))
- do k=1,3
-c uy(k,i)=facy*(dc_norm(k,i+1)-costh*dc_norm(k,i))
- uy(k,i)=
-c & facy*(
- & dc_norm(k,i+1)*scalar(dc_norm(1,i),dc_norm(1,i))
- & -scalar(dc_norm(1,i),dc_norm(1,i+1))*dc_norm(k,i)
-c & )
- enddo
-c write (iout,*) 'facy',facy,
-c & 1.0d0/dsqrt(scalar(uy(1,i),uy(1,i)))
- facy=1.0d0/dsqrt(scalar(uy(1,i),uy(1,i)))
- do k=1,3
- uy(k,i)=facy*uy(k,i)
- enddo
-C Compute the derivatives of uy
- do j=1,3
- do k=1,3
- uyder(k,j,1)=2*dc_norm(k,i+1)*dc_norm(j,i)
- & -dc_norm(k,i)*dc_norm(j,i+1)
- uyder(k,j,2)=-dc_norm(j,i)*dc_norm(k,i)
- enddo
-c uyder(j,j,1)=uyder(j,j,1)-costh
-c uyder(j,j,2)=1.0d0+uyder(j,j,2)
- uyder(j,j,1)=uyder(j,j,1)
- & -scalar(dc_norm(1,i),dc_norm(1,i+1))
- uyder(j,j,2)=scalar(dc_norm(1,i),dc_norm(1,i))
- & +uyder(j,j,2)
- enddo
- do j=1,2
- do k=1,3
- do l=1,3
- uygrad(l,k,j,i)=uyder(l,k,j)
- uzgrad(l,k,j,i)=uzder(l,k,j)
- enddo
- enddo
- enddo
- call unormderiv(uy(1,i),uyder(1,1,1),facy,uygrad(1,1,1,i))
- call unormderiv(uy(1,i),uyder(1,1,2),facy,uygrad(1,1,2,i))
- call unormderiv(uz(1,i),uzder(1,1,1),fac,uzgrad(1,1,1,i))
- call unormderiv(uz(1,i),uzder(1,1,2),fac,uzgrad(1,1,2,i))
- endif
- enddo
- do i=1,nres-1
- do j=1,2
- do k=1,3
- do l=1,3
- uygrad(l,k,j,i)=vblinv*uygrad(l,k,j,i)
- uzgrad(l,k,j,i)=vblinv*uzgrad(l,k,j,i)
- enddo
- enddo
- enddo
- enddo
- return
- end
-C-----------------------------------------------------------------------------
- subroutine check_vecgrad
- implicit real*8 (a-h,o-z)
- include 'DIMENSIONS'
- include 'DIMENSIONS.ZSCOPT'
- include 'COMMON.IOUNITS'
- include 'COMMON.GEO'
- include 'COMMON.VAR'
- include 'COMMON.LOCAL'
- include 'COMMON.CHAIN'
- include 'COMMON.VECTORS'
- dimension uygradt(3,3,2,maxres),uzgradt(3,3,2,maxres)
- dimension uyt(3,maxres),uzt(3,maxres)
- dimension uygradn(3,3,2),uzgradn(3,3,2),erij(3)
- double precision delta /1.0d-7/
- call vec_and_deriv
-cd do i=1,nres
-crc write(iout,'(2i5,2(3f10.5,5x))') i,1,dc_norm(:,i)
-crc write(iout,'(2i5,2(3f10.5,5x))') i,2,uy(:,i)
-crc write(iout,'(2i5,2(3f10.5,5x)/)')i,3,uz(:,i)
-cd write(iout,'(2i5,2(3f10.5,5x))') i,1,
-cd & (dc_norm(if90,i),if90=1,3)
-cd write(iout,'(2i5,2(3f10.5,5x))') i,2,(uy(if90,i),if90=1,3)
-cd write(iout,'(2i5,2(3f10.5,5x)/)')i,3,(uz(if90,i),if90=1,3)
-cd write(iout,'(a)')
-cd enddo
- do i=1,nres
- do j=1,2
- do k=1,3
- do l=1,3
- uygradt(l,k,j,i)=uygrad(l,k,j,i)
- uzgradt(l,k,j,i)=uzgrad(l,k,j,i)
- enddo
- enddo
- enddo
- enddo
- call vec_and_deriv
- do i=1,nres
- do j=1,3
- uyt(j,i)=uy(j,i)
- uzt(j,i)=uz(j,i)
- enddo
- enddo
- do i=1,nres
-cd write (iout,*) 'i=',i
- do k=1,3
- erij(k)=dc_norm(k,i)
- enddo
- do j=1,3
- do k=1,3
- dc_norm(k,i)=erij(k)
- enddo
- dc_norm(j,i)=dc_norm(j,i)+delta
-c fac=dsqrt(scalar(dc_norm(1,i),dc_norm(1,i)))
-c do k=1,3
-c dc_norm(k,i)=dc_norm(k,i)/fac
-c enddo
-c write (iout,*) (dc_norm(k,i),k=1,3)
-c write (iout,*) (erij(k),k=1,3)
- call vec_and_deriv
- do k=1,3
- uygradn(k,j,1)=(uy(k,i)-uyt(k,i))/delta
- uygradn(k,j,2)=(uy(k,i-1)-uyt(k,i-1))/delta
- uzgradn(k,j,1)=(uz(k,i)-uzt(k,i))/delta
- uzgradn(k,j,2)=(uz(k,i-1)-uzt(k,i-1))/delta
- enddo
-c write (iout,'(i5,3f8.5,3x,3f8.5,5x,3f8.5,3x,3f8.5)')
-c & j,(uzgradt(k,j,1,i),k=1,3),(uzgradn(k,j,1),k=1,3),
-c & (uzgradt(k,j,2,i-1),k=1,3),(uzgradn(k,j,2),k=1,3)
- enddo
- do k=1,3
- dc_norm(k,i)=erij(k)
- enddo
-cd do k=1,3
-cd write (iout,'(i5,3f8.5,3x,3f8.5,5x,3f8.5,3x,3f8.5)')
-cd & k,(uygradt(k,l,1,i),l=1,3),(uygradn(k,l,1),l=1,3),
-cd & (uygradt(k,l,2,i-1),l=1,3),(uygradn(k,l,2),l=1,3)
-cd write (iout,'(i5,3f8.5,3x,3f8.5,5x,3f8.5,3x,3f8.5)')
-cd & k,(uzgradt(k,l,1,i),l=1,3),(uzgradn(k,l,1),l=1,3),
-cd & (uzgradt(k,l,2,i-1),l=1,3),(uzgradn(k,l,2),l=1,3)
-cd write (iout,'(a)')
-cd enddo
- enddo
- return
- end
-C--------------------------------------------------------------------------
- subroutine set_matrices
- implicit real*8 (a-h,o-z)
- include 'DIMENSIONS'
- include 'DIMENSIONS.ZSCOPT'
- include 'COMMON.IOUNITS'
- include 'COMMON.GEO'
- include 'COMMON.VAR'
- include 'COMMON.LOCAL'
- include 'COMMON.CHAIN'
- include 'COMMON.DERIV'
- include 'COMMON.INTERACT'
- include 'COMMON.CONTACTS'
- include 'COMMON.TORSION'