TITLE     'Rectangular Duct, Pressure Driven'                     { sim251.pde }
SELECT     errlim=5e-3     ngrid=2     spectral_colors     stages=2   	
COORDINATES     cartesian3	
VARIABLES     vx     vy     vz     p
DEFINITIONS
   Lx=1.0     Ly=2.0     Lz=10.0     visc=1e4     delp=1e2
   vz_a=delp/Lz/(2*visc)*(Lx^2- x^2)	{ Solution for channel, p.235 }
   dens=1e3     Re=dens*globalmax( vz)*2*Ly/visc	
   v=vector( vx, vy, vz)     vm=magnitude( v)
   unit_x=vector(1,0,0)     unit_y=vector(0,1,0)     unit_z=vector(0,0,1)
   nx=normal( unit_x)     ny=normal( unit_y)     nz=normal( unit_z)
{ Natural boundary condition for p }
   natp=   if stage=2 then
      visc*[ nx*div( grad( vx))+ ny*div( grad( vy))+ nz*div( grad( vz))]
      else 0
EQUATIONS						{ For Re<<1 }
   vx:		dx( p)- visc*div( grad( vx))=0
   vy:		dy( p)- visc*div( grad( vy))=0
   vz:		dz( p)- visc*div( grad( vz))=0
   p:		div( grad( p))- 1e4*visc/Lx^2*div( v)=0
EXTRUSION						{ Parallel surfaces }
   surface 'bottom' z=0	
   layer 'liquid'						{ Layer containing liquid }
   surface 'top' z=Lz
BOUNDARIES
   surface 'bottom'  natural(vx)=0    natural(vy)=0   natural(vz)=0
      value(p)=delp
   surface 'top'    natural(vx)=0    natural(vy)=0   natural(vz)=0
   value(p)=0
region 'domain'						{ Curve to be extruded }
   start 'outer' (-Lx,-Ly)   value( vx)=0   value( vy)=0   value( vz)=0	
   natural( p)=natp   line to (Lx,-Ly)  to (Lx,Ly)  to (-Lx,Ly)  close
MONITORS
   contour( vz) painted on x=0     contour( vz) painted on z=5.0
PLOTS
   contour( p) painted on x=0 report(Re)
   contour( vz) painted on x=0   contour( vz) painted on y=0
   contour( vz) painted on z=0   contour( vz) painted on z=5.0
   contour( vz) painted on z=Lz
   elevation( p) from (0,0,0) to (0,0,Lz)
   elevation( vz, vz_a) from (-Lx,0,Lz) to (Lx,0,Lz)
   elevation( vz, vz_a) from (0,0,0) to (0,0,Lz)
END