and discretize
UCES
Course: Methods and Analysis in UCES
Applied Area: irrotational and incompressible fluids, fluid flow in porous media, Darcy's law, ideal fluids
Model: partial differential equation
Method: finite difference, SOR
Assessment: existence of steam functions
Prerequisites:
Objectives:
General Information: We present two fluid flow models in 2D: flow in a saturated porous media and ideal fluids. Both these models are similar to steady state 2D heat diffusion.
Contact Person: R. E. White, NCSU,white@math.ncsu.edu
Revision Date: 8-3-94
Copyright ©1994, R. E. White. All rights reserved.
We present two fluid flow models in 2D: flow in a porous media and ideal fluids. Both these models are similar to steady state 2D heat diffusion. The porous media flow uses a potential formulation, and we apply this model to groundwater management. The ideal fluid flow model uses a stream function formulation, and we apply it to flow around an obstacle.
In both applications we will assume the velocity of the fluid has the form
The compressibility of the fluid can be quantified by
the divergence of the velocity. In 2D the divergence of (u,v) is
and discretize
| flow in and out of vol. (dx dx T) = | T dy (u(x+dx, y)
- u(x,y)) dt + T dx (v(x, y+dy) - v(x,y)) dt.
|
Divide by (dx dy T)dt and let dx and dy go to zero to get
(ux + vy ).
 |
| Figure: Incompressible and Irrotational 2D Fluid |
A fluid with no circulation or rotation can be described
by the curl of the velocity vector. In 2D the curl of (u, v) is
is
| momemtum = | A dy (v(x+dx, y)
- v(x,y)) - A dx (u(x, y+dy)
- u(x,y)). |
Divide by (A dy dx) and let dx and
dy go to zero to get
Application to Saturated 2D Fluid Flow in a Porous Medium. Consider a shallow saturated porous medium which is to have at least one well. We will assume the region is in the xy-plane and that the water moves towards the well in such a way that the velocity vector is in the xy-plane. At the top and bottom of the xy region we will assume there is no flow through these boundaries. However, assume there is ample supply from the left and right boundaries so that the pressure is fixed. The problem is to determine the flow rates of well(s), location of well(s) and number of wells so that there is still water to be pumped out.
If a cell does not contain a well and is in the interior, then
So, we have
 |
| Figure: Groundwater 2D Porous Flow |
Application to Ideal 2D Fluid Flow. The figure below depicts the flow about an obstacle. Because the fluid is not compressible it must significantly increase its speed in order to pass near the obstacle. This can cause severe erosion of the nearby soil. The problem is to determine these velocities.
Ideal 2D steady state fluid flow is defined to be incompressible and
irrotational. Thus, both 
, such that
x ,
y ) = (-v,u).xx -
yy = 0.
We call a stream function because
the curves 
),
y()) to a
constant are parallel to the velocity vectors (u,v). In order to see this,
let (x(), y())
= c, and use the chain rule:
Thus, (u,v) and the tangent to the curve given by
(x(), y())
= c
 |
| Figure: Ideal Fluid Flow Around an Obstacle |
Both models have a partial differential equation similar to that of the 2D heat diffusion model, but all three have different boundary conditions. For our present fluid flow problems, they are either a given function along part of the boundary, or they are a zero derivative for the remainder of the boundary.
Ideal Flow Around an Obstacle.
In both problems we will use the finite difference method coupled with the
SOR iterative method. For the (dx dy) cells in the interior this is exactly
as in the 2D heat diffusion problem. For the portions of the boundary where
the derivative is set equal zero on a half cell
)* u(i,j) +
*utemp.
In the following implementations observe where the extra lines of code are that reflect these derivative boundary conditions.
The groundwater model uses the following parameters:
| L = 5,000 | dx = h = 100 | xw = (iw-1)h | h = 100 |
| H = 1,000 | dy = h = 100 | yw = (jw-1)h | K = 10. |
A single well with a flow rate of -1000 was used in the first numerical experiment. The first output graphs are plots of the hydraulic head pressure as a function of x and y. Note the pressure near the well has dropped from 100 to about 30. The second experiment has two wells with the same flow rate. In this case the pressures are negative near both wells. This indicates that before any steady state solution was achieved, the wells went dry!
Maple Code for Fluid Flow in 2D Saturated Porous Medium.
Input Data:
> with(linalg):
> with(plots):
> K:=10;
> well:= -1000:
> iw:=6:
> jw:=16:
> eps:=.0001:
> nx:= 50:
> ny:=10:
> H:=1000:
> w:=1.7:
Procedure BVP2DPOR is Defined:
> BVP2DPOR:=proc(well,iw,jw,K,nx,ny,H,w,eps):
> u:=matrix(nx+1,ny+1,100):
> h:=H/ny:
> maxit:=400:
> tol:=eps*h*h:
> for m from 1 to maxit do
> numi:=0:
> j:= 1:
> for i from 2 to nx do
> utemp:= evalf(( 2*u[i,j+1]
+ u[i+1,j] + u[i-1,j])/4):
> utemp:= (1-w)*u[i,j] + w*utemp:
> error:= abs(utemp - u[i,j]):
> u[i,j]:=utemp:
> if error < tol then
numi:=numi +1
fi:
> od:
> for j from 2 to ny do
> for i from 2 to nx do
> utemp:= evalf(( u[i,j-1] + u[i-1,j]
+ u[i+1,j] + u[i,j+1])/4):
> if (i=iw and j=jw) then
utemp:= evalf(( u[i,j-1] + u[i-1,j]
+ u[i+1,j] + u[i,j+1] + well/K)/4):
> fi:
> utemp:= (1-w)*u[i,j] + w*utemp:
> error:= abs(utemp - u[i,j]):
> u[i,j]:=utemp:
> if error < tol then
numi:=numi +1
fi:
> od:
> od:
> j:= ny + 1:
> for i from 2 to nx do
> utemp:= evalf(( 2*u[i,j-1]
+ u[i+1,j] + u[i-1,j])/4):
> utemp:= (1-w)*u[i,j] + w*utemp:
> error:= abs(utemp - u[i,j]):
> u[i,j]:=utemp:
> if error < tol then
numi:=numi +1
fi:
> od:
> if numi = (nx-1)*(ny+1) then break fi:
> od:
> end:
Output Data:
> BVP2D(well,16,6,K,50,10,H,1.7,.0001):
> matrixplot(u,style=wireframe,axes=boxed);
> m;
27
> numi;
539
| |
| Figure: Well at (16,6) with Flow Rate 1000 |
| |
| Figure: Wells at (16,6) and (32,4) with Flow Rates of 1000 |
The code for the ideal flow around an obstacle has the parameters:
| L = 500 | dx = h = 20 | xp = (ip-1)h | u0 = 1 | ||||||
| H = 100 | dy = h = 20 | yp = (jp-1)h. |
 |
| Figure: Ideal Flow Around an Obstacle |
 |
| Figure: Ideal Flow Around a Bigger Obstacle |
Both these models have enough assumptions to rule out many real applications. For groundwater problems the soils are usually not fully saturated, and the hydraulic conductivity can be highly nonlinear and can vary with space according to the soil types. Often the soils are very heterogeneous, and the soil properties are unknown. Both models may require 3D calculations and irregular shaped domains. Fluid flow may often compressible and irrotational, for example, flow over an aircraft or weather prediction models. The good news is that the more complicated models have many subproblems which are similar to our present models from heat diffusion, fluid flow in saturated porous media and ideal fluid flow.
The existence of stream functions such that
x ,
y ) = (-v,u)
Let ux + vy = 0, let
Define the stream function to be the line integral of
x ,
y ) = (-v,u),x = -v = P.
The proof that y = u = Q
is similar, but one should use the other path
 |
| Figure: Two Paths to (x,y) |
Existence of Stream Function Theorem. If
u(x,y) and v(x,y) have continuous first order partial derivatives and
x ,
y ) = (-v,u)
and eps. Observe the number of iterations
required for convergence.
and eps. Observe the number of iterations
required for convergence.
