Page 385 - Handbook Of Multiphase Flow Assurance

384 10. Research methods in flow assurance

call Initrotate
c
do i=1,natom
xn(i)=xms(i)*r11+yms(i)*r21+zms(i)*r31+r41
yn(i)=xms(i)*r12+yms(i)*r22+zms(i)*r32+r42
zn(i)=xms(i)*r13+yms(i)*r23+zms(i)*r33+r43
end do
C
RETURN
END
C
C-------------------------------------------------------------------
C
SUBROUTINE NEWPIV
C this subroutine is used to rotate any number of 3D points about axis
c passing through 2 points with coordinates (a1,a2,a3) and (v1,v2,v3).
c angle of rotation ALPHA is in radians
implicit real*8 (a-h,o-z)
real*4 ran
integer*2 ncenter,ncenter2

common /matrix/ r11,r12,r13,r21,r22,r23,r31,r32,r33,r41,r42,r43
COMMON /types / iatomtype(1000,2),nbondfound,natom,ngroups,ncall
COMMON /chain / monodata(100,100),monocount(99),nbackbone
COMMON /NWPOS / xn(1000), yn(1000), zn(1000)
COMMON /STEPS / Deltrn, Delrot, Delpiv, Acctrn, Accrot, Accpiv
COMMON /SPACE / xms(1000), yms(1000), zms(1000)
COMMON /center/ ncenter,ncenter2
common /rotate/ A1,A2,A3,v1,v2,v3,alpha

C identify call type (pivot move)
ncall=3
c choose angle of rotation; and two backbone atoms as rotation axis
ran=roulet()
alpha=DELPIV*ran
ran=(roulet()+1.)*0.5
ncenter=int(nbackbone*ran)+1
ncenter2=ncenter+1
if(ncenter.gt.nbackbone/2) ncenter2=ncenter-1
a1=xms(monodata(100,ncenter))
a2=yms(monodata(100,ncenter))
a3=zms(monodata(100,ncenter))
v1=xms(monodata(100,ncenter2))
v2=yms(monodata(100,ncenter2))
v3=zms(monodata(100,ncenter2))