X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=blobdiff_plain;f=source%2Funres%2Fsrc_MD-NEWSC%2Fchainbuild.F;fp=source%2Funres%2Fsrc_MD-NEWSC%2Fchainbuild.F;h=45a1a53837e28da1c8198f61d4d9973d13076818;hb=7308760ff07636ef6b1ee28d8c3a67a23c14b34b;hp=0000000000000000000000000000000000000000;hpb=9a54ab407f6d0d9d564d52763b3e2136450b9ffc;p=unres.git

diff --git a/source/unres/src_MD-NEWSC/chainbuild.F b/source/unres/src_MD-NEWSC/chainbuild.F
new file mode 100644
index 0000000..45a1a53
--- /dev/null
+++ b/source/unres/src_MD-NEWSC/chainbuild.F
@@ -0,0 +1,274 @@
+      subroutine chainbuild
+C 
+C Build the virtual polypeptide chain. Side-chain centroids are moveable.
+C As of 2/17/95.
+C
+      implicit real*8 (a-h,o-z)
+      include 'DIMENSIONS'
+      include 'COMMON.CHAIN'
+      include 'COMMON.LOCAL'
+      include 'COMMON.GEO'
+      include 'COMMON.VAR'
+      include 'COMMON.IOUNITS'
+      include 'COMMON.NAMES'
+      include 'COMMON.INTERACT'
+      logical lprn
+C Set lprn=.true. for debugging
+      lprn = .false.
+C
+C Define the origin and orientation of the coordinate system and locate the
+C first three CA's and SC(2).
+C
+      call orig_frame
+*
+* Build the alpha-carbon chain.
+*
+      do i=4,nres
+	call locate_next_res(i)
+      enddo     
+C
+C First and last SC must coincide with the corresponding CA.
+C
+      do j=1,3
+	dc(j,nres+1)=0.0D0
+        dc_norm(j,nres+1)=0.0D0
+	dc(j,nres+nres)=0.0D0
+        dc_norm(j,nres+nres)=0.0D0
+        c(j,nres+1)=c(j,1)
+        c(j,nres+nres)=c(j,nres)
+      enddo
+*
+* Temporary diagnosis
+*
+      if (lprn) then
+
+      call cartprint
+      write (iout,'(/a)') 'Recalculated internal coordinates'
+      do i=2,nres-1
+	do j=1,3
+	  c(j,maxres2)=0.5D0*(c(j,i-1)+c(j,i+1))
+        enddo
+        be=0.0D0
+        if (i.gt.3) be=rad2deg*beta(i-3,i-2,i-1,i)
+        be1=rad2deg*beta(nres+i,i,maxres2,i+1)
+        alfai=0.0D0
+        if (i.gt.2) alfai=rad2deg*alpha(i-2,i-1,i)
+        write (iout,1212) restyp(itype(i)),i,dist(i-1,i),
+     &  alfai,be,dist(nres+i,i),rad2deg*alpha(nres+i,i,maxres2),be1
+      enddo   
+ 1212 format (a3,'(',i3,')',2(f10.5,2f10.2))
+
+      endif
+
+      return
+      end
+c-------------------------------------------------------------------------
+      subroutine orig_frame
+C
+C Define the origin and orientation of the coordinate system and locate 
+C the first three atoms.
+C
+      implicit real*8 (a-h,o-z)
+      include 'DIMENSIONS'
+      include 'COMMON.CHAIN'
+      include 'COMMON.LOCAL'
+      include 'COMMON.GEO'
+      include 'COMMON.VAR'
+      cost=dcos(theta(3))
+      sint=dsin(theta(3))
+      t(1,1,1)=-cost
+      t(1,2,1)=-sint 
+      t(1,3,1)= 0.0D0
+      t(2,1,1)=-sint
+      t(2,2,1)= cost
+      t(2,3,1)= 0.0D0
+      t(3,1,1)= 0.0D0
+      t(3,2,1)= 0.0D0
+      t(3,3,1)= 1.0D0
+      r(1,1,1)= 1.0D0
+      r(1,2,1)= 0.0D0
+      r(1,3,1)= 0.0D0
+      r(2,1,1)= 0.0D0
+      r(2,2,1)= 1.0D0
+      r(2,3,1)= 0.0D0
+      r(3,1,1)= 0.0D0
+      r(3,2,1)= 0.0D0
+      r(3,3,1)= 1.0D0
+      do i=1,3
+        do j=1,3
+          rt(i,j,1)=t(i,j,1)
+        enddo
+      enddo
+      do i=1,3
+        do j=1,3
+          prod(i,j,1)=0.0D0
+          prod(i,j,2)=t(i,j,1)
+        enddo
+        prod(i,i,1)=1.0D0
+      enddo   
+      c(1,1)=0.0D0
+      c(2,1)=0.0D0
+      c(3,1)=0.0D0
+      c(1,2)=vbld(2)
+      c(2,2)=0.0D0
+      c(3,2)=0.0D0
+      dc(1,0)=0.0d0
+      dc(2,0)=0.0D0
+      dc(3,0)=0.0D0
+      dc(1,1)=vbld(2)
+      dc(2,1)=0.0D0
+      dc(3,1)=0.0D0
+      dc_norm(1,0)=0.0D0
+      dc_norm(2,0)=0.0D0
+      dc_norm(3,0)=0.0D0
+      dc_norm(1,1)=1.0D0
+      dc_norm(2,1)=0.0D0
+      dc_norm(3,1)=0.0D0
+      do j=1,3
+        dc_norm(j,2)=prod(j,1,2)
+	dc(j,2)=vbld(3)*prod(j,1,2)
+	c(j,3)=c(j,2)+dc(j,2)
+      enddo
+      call locate_side_chain(2)
+      return
+      end
+c-----------------------------------------------------------------------------
+      subroutine locate_next_res(i)
+C
+C Locate CA(i) and SC(i-1)
+C
+      implicit real*8 (a-h,o-z)
+      include 'DIMENSIONS'
+      include 'COMMON.CHAIN'
+      include 'COMMON.LOCAL'
+      include 'COMMON.GEO'
+      include 'COMMON.VAR'
+      include 'COMMON.IOUNITS'
+      include 'COMMON.NAMES'
+      include 'COMMON.INTERACT'
+C
+C Define the rotation matrices corresponding to CA(i)
+C
+#ifdef OSF
+      theti=theta(i)
+      if (theti.ne.theti) theti=100.0     
+      phii=phi(i)
+      if (phii.ne.phii) phii=180.0     
+#else
+      theti=theta(i)      
+      phii=phi(i)
+#endif
+      cost=dcos(theti)
+      sint=dsin(theti)
+      cosphi=dcos(phii)
+      sinphi=dsin(phii)
+* Define the matrices of the rotation about the virtual-bond valence angles
+* theta, T(i,j,k), virtual-bond dihedral angles gamma (miscalled PHI in this
+* program), R(i,j,k), and, the cumulative matrices of rotation RT
+      t(1,1,i-2)=-cost
+      t(1,2,i-2)=-sint 
+      t(1,3,i-2)= 0.0D0
+      t(2,1,i-2)=-sint
+      t(2,2,i-2)= cost
+      t(2,3,i-2)= 0.0D0
+      t(3,1,i-2)= 0.0D0
+      t(3,2,i-2)= 0.0D0
+      t(3,3,i-2)= 1.0D0
+      r(1,1,i-2)= 1.0D0
+      r(1,2,i-2)= 0.0D0
+      r(1,3,i-2)= 0.0D0
+      r(2,1,i-2)= 0.0D0
+      r(2,2,i-2)=-cosphi
+      r(2,3,i-2)= sinphi
+      r(3,1,i-2)= 0.0D0
+      r(3,2,i-2)= sinphi
+      r(3,3,i-2)= cosphi
+      rt(1,1,i-2)=-cost
+      rt(1,2,i-2)=-sint
+      rt(1,3,i-2)=0.0D0
+      rt(2,1,i-2)=sint*cosphi
+      rt(2,2,i-2)=-cost*cosphi
+      rt(2,3,i-2)=sinphi
+      rt(3,1,i-2)=-sint*sinphi
+      rt(3,2,i-2)=cost*sinphi
+      rt(3,3,i-2)=cosphi
+      call matmult(prod(1,1,i-2),rt(1,1,i-2),prod(1,1,i-1))
+      do j=1,3
+        dc_norm(j,i-1)=prod(j,1,i-1)
+        dc(j,i-1)=vbld(i)*prod(j,1,i-1)
+        c(j,i)=c(j,i-1)+dc(j,i-1)
+      enddo
+cd    print '(2i3,2(3f10.5,5x))', i-1,i,(dc(j,i-1),j=1,3),(c(j,i),j=1,3)
+C 
+C Now calculate the coordinates of SC(i-1)
+C
+      call locate_side_chain(i-1)
+      return
+      end
+c-----------------------------------------------------------------------------
+      subroutine locate_side_chain(i)
+C 
+C Locate the side-chain centroid i, 1 < i < NRES. Put in C(*,NRES+i).
+C
+      implicit real*8 (a-h,o-z)
+      include 'DIMENSIONS'
+      include 'COMMON.CHAIN'
+      include 'COMMON.LOCAL'
+      include 'COMMON.GEO'
+      include 'COMMON.VAR'
+      include 'COMMON.IOUNITS'
+      include 'COMMON.NAMES'
+      include 'COMMON.INTERACT'
+      dimension xx(3)
+
+c      dsci=dsc(itype(i))
+c      dsci_inv=dsc_inv(itype(i))
+      dsci=vbld(i+nres)
+      dsci_inv=vbld_inv(i+nres)
+#ifdef OSF
+      alphi=alph(i)
+      omegi=omeg(i)
+      if (alphi.ne.alphi) alphi=100.0
+      if (omegi.ne.omegi) omegi=-100.0
+#else
+      alphi=alph(i)
+      omegi=omeg(i)
+#endif
+      cosalphi=dcos(alphi)
+      sinalphi=dsin(alphi)
+      cosomegi=dcos(omegi)
+      sinomegi=dsin(omegi) 
+      xp= dsci*cosalphi
+      yp= dsci*sinalphi*cosomegi
+      zp=-dsci*sinalphi*sinomegi
+* Now we have to rotate the coordinate system by 180-theta(i)/2 so as to get its
+* X-axis aligned with the vector DC(*,i)
+      theta2=pi-0.5D0*theta(i+1)
+      cost2=dcos(theta2)
+      sint2=dsin(theta2)
+      xx(1)= xp*cost2+yp*sint2
+      xx(2)=-xp*sint2+yp*cost2
+      xx(3)= zp
+cd    print '(a3,i3,3f10.5,5x,3f10.5)',restyp(itype(i)),i,
+cd   &   xp,yp,zp,(xx(k),k=1,3)
+      do j=1,3
+        xloc(j,i)=xx(j)
+      enddo
+* Bring the SC vectors to the common coordinate system.
+      xx(1)=xloc(1,i)
+      xx(2)=xloc(2,i)*r(2,2,i-1)+xloc(3,i)*r(2,3,i-1)
+      xx(3)=xloc(2,i)*r(3,2,i-1)+xloc(3,i)*r(3,3,i-1)
+      do j=1,3
+	xrot(j,i)=xx(j)
+      enddo
+      do j=1,3
+        rj=0.0D0
+        do k=1,3
+          rj=rj+prod(j,k,i-1)*xx(k)
+        enddo
+        dc(j,nres+i)=rj
+        dc_norm(j,nres+i)=rj*dsci_inv
+        c(j,nres+i)=c(j,i)+rj
+      enddo
+      return
+      end