TITLE     'Flow through a Constricted Channel'		         { sim201.pde }
SELECT     errlim=1e-5     ngrid=1     spectral_colors
VARIABLES     phi                                	{ Velocity potential }
DEFINITIONS					{ SI units }
   Lx=1     Ly=1             		
   coef=0.5						{ Constriction coefficient }
   vx0=3.0                				{ Velocity at input end }
   p0=1e5              				{ Atmospheric pressure }
   dens=1e3            				{ Mass density }
   vx=dx(phi)     vy=dy(phi)      		{ Velocity components }
   v=vector( vx,vy)     vm=sqrt( vx^2+ vy^2)		{ Speed }
   p=p0+ 1/2*dens*(vx0^2-vm^2)         			{ Pressure }
   div_v=dx( vx)+ dy( vy)			{ Divergence, or div( v) }
   curl_z=dx( vy)- dy( vx)			{ Vorticity, or curl( v) }
EQUATIONS
   dxx( phi)+ dyy( phi)=0			{ Or div( grad( v)) }
BOUNDARIES
region 'domain' start 'outer' (0,Ly)
   natural( phi)=-vx0  line to (0,-Ly)				{ In }
   natural( phi)=0 line to (Lx,-Ly)  to (2*Lx,-Ly*coef)  to (3*Lx,-Ly*coef)
   value( phi)=0 line  to (3*Lx,Ly*coef)  				{ Out }
   natural( phi)=0 line to (2*Lx,Ly*coef) to (Lx,Ly) close
PLOTS
   contour( phi)     vector( v) norm     contour( vm) painted
   contour( p) painted     contour( p) zoom(1.5*Lx,0,  Lx,Ly)
   surface( p) zoom(1.5*Lx,0,  Lx,Ly)
   elevation( vm) on 'outer'    		{ Verify boundary conditions }
   contour( div_v)     contour( curl_z)
END