Page 434 - Handbook Of Multiphase Flow Assurance

Computer code (Makogon, 1994, 1997) 433

C do i=1, npts
C write(2,50)i,x1(i),y1(i),z1(i),x2(i),y2(i),z2(i),x3(i),y3(i),z3(i)
C end do
C50 format(I3,9F8.4)
C***************************************************************
C END OF THE HYDRATE SURFACE CODE
Nwatrs=npts
C move the water coordinates to the water array and scale the coordinates
C ... positions of oxygen atoms
DO 100 I = 1, Nwatrs
xw(I,3) = x1(i) * unitside
yw(I,3) = y1(i) * unitside
zw(I,3) = z1(i) * unitside
C ... positions of hydrogens
xw(I,1) = x2(I) * unitside
yw(I,1) = y2(I) * unitside
zw(I,1) = z2(I) * unitside
C
xw(I,2) = x3(I) * unitside
yw(I,2) = y3(I) * unitside
zw(I,2) = z3(I) * unitside
100 CONTINUE
C code for writing the water box coordinates to a file surface.dat
C open(unit=2,file='surface.dat',status='unknown')
C do i=1, npts
C write(2,50)i,xw(i,3),yw(i,3),zw(i,3),xw(i,1),yw(i,1),zw(i,1),
C _ xw(i,2),yw(i,2),zw(i,2)
C end do
C50 format(I3,9F8.4)
C
C some code to print the surface positions visually on screen
do 110 i=1,80
do 110 j=1,25
110 screen(i,j)=' '
if (ap.eq.1.or.ap.eq.3)then
C construct a [100] or a [110] surface
do 115 i=1,Nwatrs
if(zw(i,3).ge.-5.D0) then
C for oxygens
nx=+nint(xw(i,3)*2.D0)
ny=25-nint(yw(i,3))
C and for hydrogens
lx=nint(xw(i,2)*2.D0)
ly=25-nint(yw(i,2))
mx=nint(xw(i,1)*2.D0)