source/unres/src_MD/bond_move.f

   1       subroutine bond_move(nbond,nstart,psi,lprint,error)
   2 C Move NBOND fragment starting from the CA(nstart) by angle PSI.
   3       implicit real*8 (a-h,o-z)
   4       include 'DIMENSIONS'
   5       integer nbond,nstart
   6       double precision psi
   7       logical fail,error,lprint
   8       include 'COMMON.GEO'
   9       include 'COMMON.CHAIN'
  10       include 'COMMON.VAR'
  11       include 'COMMON.REFSYS'
  12       include 'COMMON.IOUNITS'
  13       include 'COMMON.MCM'
  14       dimension x(3),e(3,3),rot(3,3),trans(3,3)
  15       error=.false.
  16       nend=nstart+nbond
  17       if (print_mc.gt.2) then
  18       write (iout,*) 'nstart=',nstart,' nend=',nend,' nbond=',nbond
  19       write (iout,*) 'psi=',psi
  20       write (iout,'(a)') 'Original coordinates of the fragment'
  21       do i=nstart,nend
  22         write (iout,'(i5,3f10.5)') i,(c(j,i),j=1,3)
  23       enddo
  24       endif
  25       if (nstart.lt.1 .or. nend .gt.nres .or. nbond.lt.2 .or.
  26      & nbond.ge.nres-1) then
  27         write (iout,'(a)') 'Bad data in BOND_MOVE.'
  28         error=.true.
  29         return
  30       endif
  31 C Generate the reference system.
  32       i2=nend
  33       i3=nstart
  34       i4=nstart+1
  35       call refsys(error)
  36 C Return, if couldn't define the reference system.
  37       if (error) return
  38 C Compute the transformation matrix.
  39       cospsi=dcos(psi)
  40       sinpsi=dsin(psi)
  41       rot(1,1)=1.0D0
  42       rot(1,2)=0.0D0
  43       rot(1,3)=0.0D0
  44       rot(2,1)=0.0D0
  45       rot(2,2)=cospsi
  46       rot(2,3)=-sinpsi
  47       rot(3,1)=0.0D0
  48       rot(3,2)=sinpsi
  49       rot(3,3)=cospsi
  50       do i=1,3
  51         e(1,i)=e1(i)
  52         e(2,i)=e2(i)
  53         e(3,i)=e3(i)
  54       enddo
  55
  56       if (print_mc.gt.2) then
  57       write (iout,'(a)') 'Reference system and matrix r:'
  58       do i=1,3
  59         write(iout,'(i5,2(3f10.5,5x))')i,(e(i,j),j=1,3),(rot(i,j),j=1,3)
  60       enddo
  61       endif
  62
  63       call matmult(rot,e,trans)
  64       do i=1,3
  65         do j=1,3
  66           e(i,1)=e1(i)
  67           e(i,2)=e2(i)
  68           e(i,3)=e3(i)
  69         enddo
  70       enddo
  71       call matmult(e,trans,trans)
  72
  73       if (lprint) then
  74       write (iout,'(a)') 'The trans matrix:'
  75       do i=1,3
  76         write (iout,'(i5,3f10.5)') i,(trans(i,j),j=1,3)
  77       enddo
  78       endif
  79
  80       do i=nstart,nend
  81         do j=1,3
  82           rij=c(j,nstart)
  83           do k=1,3
  84             rij=rij+trans(j,k)*(c(k,i)-c(k,nstart))
  85           enddo
  86           x(j)=rij
  87         enddo
  88         do j=1,3
  89           c(j,i)=x(j)
  90         enddo
  91       enddo
  92
  93       if (lprint) then
  94       write (iout,'(a)') 'Rotated coordinates of the fragment'
  95       do i=nstart,nend
  96         write (iout,'(i5,3f10.5)') i,(c(j,i),j=1,3)
  97       enddo
  98       endif
  99
 100 c     call int_from_cart(.false.,lprint)
 101       if (nstart.gt.1) then
 102         theta(nstart+1)=alpha(nstart-1,nstart,nstart+1)
 103         phi(nstart+2)=beta(nstart-1,nstart,nstart+1,nstart+2)
 104         if (nstart.gt.2) phi(nstart+1)=
 105      &                beta(nstart-2,nstart-1,nstart,nstart+1)
 106       endif
 107       if (nend.lt.nres) then
 108         theta(nend+1)=alpha(nend-1,nend,nend+1)
 109         phi(nend+1)=beta(nend-2,nend-1,nend,nend+1)
 110         if (nend.lt.nres-1) phi(nend+2)=
 111      &                beta(nend-1,nend,nend+1,nend+2)
 112       endif
 113       if (print_mc.gt.2) then
 114       write (iout,'(/a,i3,a,i3,a/)')
 115      & 'Moved internal coordinates of the ',nstart,'-',nend,
 116      & ' fragment:'
 117       do i=nstart+1,nstart+2
 118         write (iout,'(i5,2f10.5)') i,rad2deg*theta(i),rad2deg*phi(i)
 119       enddo
 120       do i=nend+1,nend+2
 121         write (iout,'(i5,2f10.5)') i,rad2deg*theta(i),rad2deg*phi(i)
 122       enddo
 123       endif
 124       return
 125       end