Page 426 - Handbook Of Multiphase Flow Assurance

Computer code (Makogon, 1994, 1997) 425

C check the cell in direction of equation 2: x-2y+z=1
xyz=x1(i)-y1(i)-y1(i)+z1(i)
if(abs(xyz-dnint(xyz)).lt.1d-10) xyz=dnint(xyz)
if(xyz.GT.1.D0) remflag(i)=1
C check the cell in direction of equation 3: x-2y+z=-2
xyz=x1(i)-y1(i)-y1(i)+z1(i)
if(abs(xyz-dnint(xyz)).lt.1d-10) xyz=dnint(xyz)
if(xyz.LE.-2.D0) remflag(i)=1
C check the cell in direction of equation 4: -2x+y+z=-2
xyz=-x1(i)-x1(i)+y1(i)+z1(i)
if(abs(xyz-dnint(xyz)).lt.1d-10) xyz=dnint(xyz)
if(xyz.LE.-2.D0) remflag(i)=1
C check the cell in direction of equation 5: x+y+z=1
xyz=x1(i)+y1(i)+z1(i)
if(abs(xyz-dnint(xyz)).lt.1d-10) xyz=dnint(xyz)
if(xyz.GT.1.D0) remflag(i)=1
C check the cell in direction of equation 6: x+y+z=-2
xyz=x1(i)+y1(i)+z1(i)
if(abs(xyz-dnint(xyz)).lt.1d-10) xyz=dnint(xyz)
if(xyz.LE.-2.D0) remflag(i)=1
end do

C move the unremoved molecules to the front of coordinate arrays
nrem=0
do i=1,npts*12
if(remflag(i).eq.0) then
x1(i-nrem)=x1(i)
y1(i-nrem)=y1(i)
z1(i-nrem)=z1(i)
x2(i-nrem)=x2(i)
y2(i-nrem)=y2(i)
z2(i-nrem)=z2(i)
x3(i-nrem)=x3(i)
y3(i-nrem)=y3(i)
z3(i-nrem)=z3(i)
else
nrem=nrem+1
endif
enddo
C modify the number of molecules in [111] slice
npts=npts*12-nrem

C reorient the basic [111] block to align the top surface with
C the x-y plane at position z=0 and with center at x=0, y=0
do i=1, npts