Page 433 - Handbook Of Multiphase Flow Assurance

432 10. Research methods in flow assurance

side=side*2.D0
xmaxh = unitside * Dsqrt(2.D0)
ymaxh = unitside * Dsqrt(1.5D0)
xmax = 2.D0 * xmaxh
ymax = 2.D0 * ymaxh
C now center a [111] surface
do i=1,npts
x1(i)=x1(i)-xmaxh/unitside*0.75D0
y1(i)=y1(i)-ymaxh/unitside*0.5D0
x2(i)=x2(i)-xmaxh/unitside*0.75D0
y2(i)=y2(i)-ymaxh/unitside*0.5D0
x3(i)=x3(i)-xmaxh/unitside*0.75D0
y3(i)=y3(i)-ymaxh/unitside*0.5D0
end do
endif

C check if geometries of water molecules are still intact
do i=1, npts
d1=(x1(i)-x2(i))**2+(y1(i)-y2(i))**2+(z1(i)-z2(i))**2
d1=dsqrt(d1) * unitside
d2=(x1(i)-x3(i))**2+(y1(i)-y3(i))**2+(z1(i)-z3(i))**2
d2=dsqrt(d2) * unitside
if(d1.gt.1.1.or.d2.gt.1.1.or.d1.lt.0.9.or.d2.lt.0.9) then
write(*,*)'Error in water geometry O-H bond length ',d1,d2
endif
end do

write(6,47)npts
47 format(1X,'A surface consisting of',I5,' water molecules was gener
%ated')
write(*,*)'Top surface was placed at z=0'
if (ap.eq.1) then
C write(6,49)side,side,C
write(6,49)xmax,ymax,C
write(6,48)xmaxh,ymaxh
else
C write(6,49)side*sqrt(2.),side*sqrt(1.5),C
write(6,49)xmax,ymax,C
write(*,*)'Center of the surface is at x=0, y=0'
endif
48 FORMAT('Center of the surface is at x=',F5.2,', y=',F5.2)
49 FORMAT(' Maximum surface size along axes x=',F6.2,', y=',F6.2,
_', z=',F6.2)
C code for writing the water box coordinates to a file surface.dat
C open(unit=2,file='surface.dat',status='unknown')