Chapter 32, Problem 9P

The velocity of water flow through the porous media can be related to head by D'Arcy's law

$q_n=-K\frac{dh}{dn}$

where K is the hydraulic conductivity and $q_n$ is discharge velocity in the n direction. If $K=5\times {10}^{-4}\text{ cm/s}$, compute the water velocities for Prob. 32.8.

To determine

To calculate: The water flowsvelocity through the porous media for the Prob. 32.8, if the hydraulic conductivity is $K=5\times {10}^{-4}\text{ cm/s}$.

Answer to Problem 9P

Solution:

The water flow velocity at every node is,

-5.205E-04	-5.542E-04	-6.593E-04	-7.249E-04
-5.079E-04	-5.315E-04	-6.989E-04	-7.942E-04
-4.668E-04	-3.967E-04	-4.429E-04

Explanation of Solution

Given Information:

Write the expression for D’Arcy’s law.

$q_n=-K\frac{dh}{dn}$

Here, $K$ is the hydraulic conductivity and $q_n$ is the discharge velocity in $n$ direction.

The hydraulic conductivity is $K=5\times {10}^{-4}\text{ cm/s}$.

Formula used:

Consider the Laplace Equation,

$\frac{{\partial }^{2}h}{\partial x^2}+\frac{{\partial }^{2}h}{\partial y^2}=0$

Write the central difference approximation for the second derivative.

$\frac{h_{i+1,j}-2h_{i,j}+h_{i-1,j}}{\Delta x^2}+\frac{h_{i,j+1}-2h_{i,j}+h_{i,j-1}}{\Delta y^2}=0$

Calculation:

Refer to Figure P32.8, draw the nodal diagram.

Recall the Laplace Equation,

$\frac{{\partial }^{2}h}{\partial x^2}+\frac{{\partial }^{2}h}{\partial y^2}=0$

The central difference approximation applies for the second derivative in above Laplace equation,

$\frac{h_{i+1,j}-2h_{i,j}+h_{i-1,j}}{\Delta x^2}+\frac{h_{i,j+1}-2h_{i,j}+h_{i,j-1}}{\Delta y^2}=0$ …… (1)

At the node, $\left(i,j\right)=\left(0,0\right)$ the equation (1) reduces to

$\frac{h_{1,0}-2h_{0,0}+h_{-1,0}}{\Delta x^2}+\frac{h_{0,1}-2h_{0,0}+h_{0,-1}}{\Delta y^2}=0$ …... (2)

Approximate the all external nodes with a central finite difference,

$\begin{array}{c}\frac{dh}{dy}=\frac{h_{0,1}-h_{0,-1}}{2\Delta y}\\ h_{0,1}-h_{0,-1}=2\Delta y\frac{dh}{dy}\\ =2\Delta y\cdot 0\\ =0\end{array}$

Thus, $h_{0,1}=h_{0,-1}$ …… (3)

Now with a central finite difference, approximate the external nodes.

$\begin{array}{c}\frac{dh}{dx}=\frac{h_{1,0}-h_{-1,0}}{2\Delta x}\\ 2\Delta x\cdot \frac{dh}{dx}=h_{1,0}-h_{-1,0}\\ h_{-1,0}=h_{1,0}-2\Delta x\cdot \frac{dh}{dx}\\ =h_{1,0}-2\times 1\times 1\end{array}$

Solve further,

$h_{-1,0}=h_{1,0}-2$ …… (4)

Substitute (3) and (4) in (2).

$\begin{array}{c}\frac{h_{1,0}-2h_{0,0}+\left(h_{1,0}-2\right)}{\Delta x^2}+\frac{h_{0,1}-2h_{0,0}+h_{0,1}}{\Delta y^2}=0\\ \frac{2h_{1,0}-2h_{0,0}-2}{\Delta x^2}+\frac{2h_{0,1}-2h_{0,0}}{\Delta y^2}=0\\ 4h_{0,0}-2h_{1,0}-2h_{0,1}=-2\end{array}$

Similarly, at the node, $\left(i,j\right)=\left(0,1\right)$ the equation (1) can write,

$-h_{0,0}+4h_{0,1}-h_{0,2}-2h_{1,1}=-2$

Similarly, at the node, $\left(i,j\right)=\left(0,2\right)$ the equation (1) can write,

$-2h_{0,1}+4h_{0,2}-2h_{1,2}=-2$

Similarly, at the node, $\left(i,j\right)=\left(1,0\right)$ the equation (1) can write,

$-h_{0,0}+4h_{1,0}-h_{2,0}-2h_{1,1}=0$

Similarly, at the node, $\left(i,j\right)=\left(1,1\right)$ the equation (1) can write,

$-h_{0,1}-h_{1,0}-h_{1,2}+4h_{1,1}-h_{2,1}=0$

Similarly, at the node, $\left(i,j\right)=\left(1,2\right)$ the equation (1) can write,

$-h_{0,2}-2h_{1,1}+4h_{1,2}-h_{2,2}=0$

Similarly, at the node, $\left(i,j\right)=\left(2,0\right)$ the equation (1) can write,

$-2h_{1,0}+4h_{2,0}-2h_{2,1}=0$

Similarly, at the node, $\left(i,j\right)=\left(2,1\right)$ the equation (1) can write,

$-h_{2,0}-h_{1,1}+4h_{2,1}-h_{2,2}=20$

Similarly, at the node, $\left(i,j\right)=\left(2,2\right)$ the equation (1) can write,

$-h_{1,2}-2h_{2,1}+4h_{2,2}=20$

Thus, the system of all linear equations is,

$\begin{array}{c}4h_{0,0}-2h_{1,0}-2h_{0,1}=-2\\ -h_{0,0}+4h_{0,1}-h_{0,2}-2h_{1,1}=-2\\ -2h_{0,1}+4h_{0,2}-2h_{1,2}=-2\end{array}$

And,

$\begin{array}{c}-h_{0,0}+4h_{1,0}-h_{2,0}-2h_{1,1}=0\\ -h_{0,1}-h_{1,0}-h_{1,2}+4h_{1,1}-h_{2,1}=0\\ -h_{0,2}-2h_{1,1}+4h_{1,2}-h_{2,2}=0\end{array}$

And,

$\begin{array}{c}-2h_{1,0}+4h_{2,0}-2h_{2,1}=0\\ -h_{2,0}-h_{1,1}+4h_{2,1}-h_{2,2}=20\\ -h_{1,2}-2h_{2,1}+4h_{2,2}=20\end{array}$

Write all equation in matrix form.

$\left[\begin{array}{ccccccccc}4& -2& -2& 0& 0& 0& 0& 0& 0\\ -1& 0& 4& -1& 0& -2& 0& 0& 0\\ 0& 0& -2& 4& 0& 0& 0& 0& -2\\ -1& 4& 0& 0& -1& -2& 0& 0& 0\\ 0& -1& -1& 0& 0& 4& -1& 0& -1\\ 0& 0& 0& -1& 0& -2& 0& -1& 4\\ 0& -2& 0& 0& 4& 0& -2& 0& 0\\ 0& 0& 0& 0& -1& -1& 4& -1& 0\\ 0& 0& 0& 0& 0& 0& -2& 4& -1\end{array}\right]\left[\begin{array}{c}{h}_{00}\\ h_{10}\\ h_{01}\\ h_{02}\\ h_{20}\\ h_{11}\\ h_{21}\\ h_{22}\\ h_{12}\end{array}\right]=\left[\begin{array}{c}-2\\ -2\\ -2\\ 0\\ 0\\ 0\\ 0\\ 20\\ 20\end{array}\right]$

Use the MATLAB to solve the above equations, write the following code in MATLAB.

A =zeros(9,9);

A(1,1)=4;

A(1,2)=-2;

A(1,3)=-2;

A(2,1)=-1;

A(2,3)=4;

A(2,4)=-1;

A(2,6)=-2;

A(3,3)=-2;

A(3,4)=4;

A(3,9)=-2;

A(4,1)=-1;

A(4,2)=4;

A(4,5)=-1;

A(4,6)=-2;

A(5,2)=-1;

A(5,3)=-1;

A(5,6)=4;

A(5,7)=-1;

A(5,9)=-1;

A(6,4)=-1;

A(6,6)=-2;

A(6,8)=-1;

A(6,9)=4;

A(7,2)=-2;

A(7,5)=4;

A(7,7)=-2;

A(8,5)=-1;

A(8,6)=-1;

A(8,7)=4;

A(8,8)=-1;

A(9,7)=-2;

A(9,8)=4;

A(9,9)=-1;

%

b =[-2;-2;-2;0;0;0;0;20;20];

Sol = A\b

The output is,

Sol =

16.2279

17.1589

16.2969

16.3372

17.7853

17.3113

18.4117

18.5502

17.3775

$\left[\begin{array}{c}{h}_{00}\\ h_{10}\\ h_{01}\\ h_{02}\\ h_{20}\\ h_{11}\\ h_{21}\\ h_{22}\\ h_{12}\end{array}\right]=\left[\begin{array}{c} 16.2279 \\ 17.1589\\ 16.2969\\ 16.3372\\ 17.7853\\ 17.3113\\ 18.4117\\ 18.5502\\ 17.3775\end{array}\right]$

Thus, the distribution of head of the system is shown below.

16.3372	17.37748	18.55022	20
16.29691	17.31126	18.4117	20
16.22792	17.15894	17.78532

Now from the above table find the value of $\frac{dh}{dx}$ and $\frac{dh}{dy}$.

$\frac{dh}{dx}=-k^{\prime }\frac{h_{i+1,j}-h_{i-1,j}}{2\Delta x}$

And

$\frac{dh}{dy}=-k^{\prime }\frac{h_{i,j+1}-h_{i,j-1}}{2\Delta y}$

Calculate the value of $\frac{dh}{dx}$ for every node.

$\begin{array}{c}{\left(\frac{dh}{dx}\right)}_{0,2}=h_{1,2}-h_{0,2}\\ =17.33748-16.3372\\ =1.04029\end{array}$

And,

$\begin{array}{c}{\left(\frac{dh}{dx}\right)}_{0,1}=h_{2,1}-h_{0,1}\\ =17.31126-16.29691\\ =1.01435\end{array}$

Calculate all the value of $\frac{dh}{dx}$ the same way for every node. After calculating all the values of $\frac{dh}{dx}$ the following table will be obtained:

1.04029	1.10651	1.31126	1.44978
1.01435	1.05740	1.34437	1.58830
0.93102	0.77870	0.62638

Calculate the value of $\frac{dh}{dy}$ for every node.

$\begin{array}{c}{\left(\frac{dh}{dy}\right)}_{0,2}=h_{0,2}-h_{0,1}\\ =16.3372-16.29691\\ =0.04029\end{array}$

Calculate all the value of $\frac{dh}{dy}$ for every node same way. After calculating all the values of $\frac{dh}{dy}$ obtained the following table:

0.04029	0.06623	0.13852	0.00000
0.05464	0.10927	0.38245	0.00000
0.06898	0.15232	0.62638

Now calculate the value of $\frac{dh}{dn}$.

${\frac{dh}{dn}}_n=\sqrt{{\frac{dh}{dx}}_n^2+{\frac{dh}{dy}}_n^2}$

Substitute the values of $\frac{dh}{dx}$ and $\frac{dh}{dy}$ for every node then the value of $\frac{dh}{dn}$ to the corresponding nodes will be obtained.

For ${\left(\frac{dh}{dn}\right)}_{0,2}$,

$\begin{array}{c}{\left(\frac{dh}{dn}\right)}_{0,2}=\sqrt{{\left(\frac{dh}{dx}\right)}_{0,2}^2+{\left(\frac{dh}{dy}\right)}_{0,2}^2}\\ =\sqrt{{\left(1.04029\right)}^2+{\left(0.04029\right)}^2}\\ =1.04107\end{array}$

Calculate for every node the value of $\frac{dh}{dn}$ same way. After calculating all the values obtained the following table.

1.04107	1.10849	1.31855	1.44978
1.01582	1.06303	1.39771	1.58830
0.93357	0.79345	0.88583

Apply the D’Arcy law to find the discharge velocity in the n direction:

$q_n=-K\frac{dh}{dn}$

Here, $K=5\times {10}^{-4}\text{ cm/s}$ is the hydraulic conductivity.

Now Calculate for the velocity ($q_{0,2}$).

$\begin{array}{c}{q}_{0,2}=-K{\left(\frac{dh}{dn}\right)}_{0,2}\\ =-\left(5\times {10}^{-\text{4}}\text{ cm/s}\right)\times 1.04107\\ =-5.205\times {10}^{-5}\text{ cm/s}\end{array}$

Calculate the velocity for every node same way, and got the following table:

-5.205E-04	-5.542E-04	-6.593E-04	-7.249E-04
-5.079E-04	-5.315E-04	-6.989E-04	-7.942E-04
-4.668E-04	-3.967E-04	-4.429E-04