Page 386 - Handbook Of Multiphase Flow Assurance

Computer code (Makogon, 1994, 1997) 385

c prepare the rotation matrix
c call Initrotate(A1,A2,A3,v1,v2,v3,alpha)
call Initrotate
c see which end of the chain are we rotating
if(ncenter.gt.nbackbone/2) then
do i=ncenter,nbackbone
do j=1,monocount(i)
x=xms(monodata(i,j))
y=yms(monodata(i,j))
z=zms(monodata(i,j))
xn(monodata(i,j))=x*r11+y*r21+z*r31+r41
yn(monodata(i,j))=x*r12+y*r22+z*r32+r42
zn(monodata(i,j))=x*r13+y*r23+z*r33+r43
end do
end do
c assign new coordinates to the non-moving part of the chain
do i=1,ncenter2
do j=1,monocount(i)
xn(monodata(i,j))=xms(monodata(i,j))
yn(monodata(i,j))=yms(monodata(i,j))
zn(monodata(i,j))=zms(monodata(i,j))
end do
end do
c the Origin atom was moved, so check if it has crossed the PBC
CALL PBC

else
do i=1,ncenter
do j=1,monocount(i)
x=xms(monodata(i,j))
y=yms(monodata(i,j))
z=zms(monodata(i,j))
xn(monodata(i,j))=x*r11+y*r21+z*r31+r41
yn(monodata(i,j))=x*r12+y*r22+z*r32+r42
zn(monodata(i,j))=x*r13+y*r23+z*r33+r43
end do
end do
c assign new coordinates to the non-moving part of the chain
do i=ncenter2,nbackbone
do j=1,monocount(i)
xn(monodata(i,j))=xms(monodata(i,j))
yn(monodata(i,j))=yms(monodata(i,j))
zn(monodata(i,j))=zms(monodata(i,j))
end do
end do