Solutions Manual for Computational Fluid Dynamics for Mechanical Engineering, 1st Edition By George Qin All Chapters (1-14)

Solutions Manual for Computational Fluid
Dynamics for Mechanical Engineering, 1e by
George Qin (All Chapters)


Chapter 1
1. Show that Equation (1.14) can also be written as
𝜕𝑢 𝜕𝑢 𝜕𝑢 𝜕2𝑢 𝜕2𝑢 1 𝜕𝑝
+𝑢 +𝑣 = 𝜈 ( 2 + 2) −
𝜕𝑡 𝜕𝑥 𝜕𝑦 𝜕𝑥 𝜕𝑦 𝜌 𝜕𝑥
Solution
Equation (1.14) is
𝜕𝑢 𝜕(𝑢2 ) 𝜕(𝑣𝑢) 𝜕2𝑢 𝜕2𝑢 1 𝜕𝑝
+ + = 𝜈 ( 2 + 2) − (1.13)
𝜕𝑡 𝜕𝑥 𝜕𝑦 𝜕𝑥 𝜕𝑦 𝜌 𝜕𝑥
The left side is
𝜕𝑢 𝜕(𝑢2 ) 𝜕(𝑣𝑢) 𝜕𝑢 𝜕𝑢 𝜕𝑢 𝜕𝑣
+ + = + 2𝑢 +𝑣 +𝑢
𝜕𝑡 𝜕𝑥 𝜕𝑦 𝜕𝑡 𝜕𝑥 𝜕𝑦 𝜕𝑦
𝜕𝑢 𝜕𝑢 𝜕𝑢 𝜕𝑢 𝜕𝑣 𝜕𝑢 𝜕𝑢 𝜕𝑢
= +𝑢 +𝑣 +𝑢( + ) = +𝑢 +𝑣
𝜕𝑡 𝜕𝑥 𝜕𝑦 𝜕𝑥 𝜕𝑦 𝜕𝑡 𝜕𝑥 𝜕𝑦
since
𝜕𝑢 𝜕𝑣
+ =0
𝜕𝑥 𝜕𝑦
due to the continuity equation.
2. Derive Equation (1.17).
Solution:
From Equation (1.14)
𝜕𝑢 𝜕(𝑢2 ) 𝜕(𝑣𝑢) 𝜕2𝑢 𝜕2𝑢 1 𝜕𝑝
+ + = 𝜈 ( 2 + 2) −
𝜕𝑡 𝜕𝑥 𝜕𝑦 𝜕𝑥 𝜕𝑦 𝜌 𝜕𝑥
Define
𝑢 𝑣 𝑥𝑖 𝑡𝑈 𝑝
𝑢̃ = , 𝑣̃ = , 𝑥̃𝑖 = , 𝑡̃ = , 𝑝̃ =
𝑈 𝑈 𝐿 𝐿 𝜌𝑈 2
Equation (1.14) becomes
𝑈𝜕𝑢̃ 𝑈 2 𝜕(𝑢̃2 ) 𝑈 2 𝜕(𝑣̃𝑢 ̃) 𝜈𝑈 𝜕 2 𝑢̃ 𝜕 2 𝑢̃ 𝜌𝑈 2 𝜕𝑝̃
+ + = 2 ( 2 + 2) −
𝐿 ̃ 𝐿𝜕𝑥̃ 𝐿𝜕𝑦̃ 𝐿 𝜕𝑥̃ 𝜕𝑦̃ 𝜌𝐿 𝜕𝑥̃
𝑈 𝜕𝑡
Dividing both sides by 𝑈 2 /𝐿, Equation (1.17) follows.

3. Derive a pressure Poisson equation from Equations (1.13) through (1.15):




Downloaded by: tutorsection | Want to earn $1.236
Distribution of this document is illegal extra per year?

, 𝜕2𝑝 𝜕2𝑝 𝜕𝑢 𝜕𝑣 𝜕𝑣 𝜕𝑢
+ = 2𝜌 ( − )
𝜕𝑥 2 𝜕𝑦 2 𝜕𝑥 𝜕𝑦 𝜕𝑥 𝜕𝑦
Solution:
𝜕𝑢 𝜕𝑣
+ =0 (1.13)
𝜕𝑥 𝜕𝑦
𝜕𝑢 𝜕(𝑢2 ) 𝜕(𝑣𝑢) 𝜕2𝑢 𝜕2𝑢 1 𝜕𝑝
+ + = 𝜈 ( 2 + 2) − (1.14)
𝜕𝑡 𝜕𝑥 𝜕𝑦 𝜕𝑥 𝜕𝑦 𝜌 𝜕𝑥
𝜕𝑣 𝜕(𝑢𝑣) 𝜕(𝑣 2 ) 𝜕2𝑣 𝜕2𝑣 1 𝜕𝑝
+ + = 𝜈 ( 2 + 2) − (1.15)
𝜕𝑡 𝜕𝑥 𝜕𝑦 𝜕𝑥 𝜕𝑦 𝜌 𝜕𝑦
Taking 𝑥-derivative of each term of Equation (1.14) and 𝑦-derivative of each term of Equation (1.15),
then adding them up, we have
𝜕 𝜕𝑢 𝜕𝑣 𝜕 2 (𝑢2 ) 𝜕 2 (𝑣𝑢) 𝜕 2 (𝑣 2 )
( + )+ + 2 +
𝜕𝑡 𝜕𝑥 𝜕𝑦 𝜕𝑥 2 𝜕𝑥𝜕𝑦 𝜕𝑦 2
2 2
𝜕 𝜕 𝜕𝑢 𝜕𝑣 1 𝜕2𝑝 𝜕2𝑝
= 𝜈 ( 2 + 2) ( + ) − ( 2 + 2)
𝜕𝑥 𝜕𝑦 𝜕𝑥 𝜕𝑦 𝜌 𝜕𝑥 𝜕𝑦
Due to continuity, we have
𝜕2𝑝 𝜕2𝑝 𝜕 2 (𝑢2 ) 𝜕 2 (𝑣𝑢) 𝜕 2 (𝑣 2 )
+ = −𝜌 [ + 2 + ]
𝜕𝑥 2 𝜕𝑦 2 𝜕𝑥 2 𝜕𝑥𝜕𝑦 𝜕𝑦 2
= −2𝜌(𝑢𝑥 𝑢𝑥 + 𝑢𝑢𝑥𝑥 + 𝑢𝑥 𝑣𝑦 + 𝑢𝑣𝑥𝑦 + 𝑢𝑥𝑦 𝑣 + 𝑢𝑦 𝑣𝑥 + 𝑣𝑦 𝑣𝑦 + 𝑣𝑣𝑦𝑦 )
𝜕 𝜕 𝜕𝑢 𝜕𝑣
= −2𝜌 [(𝑢𝑥 + 𝑢 + 𝑣 ) ( + ) + 𝑢𝑦 𝑣𝑥 + 𝑣𝑦 𝑣𝑦 ]
𝜕𝑥 𝜕𝑦 𝜕𝑥 𝜕𝑦
𝜕𝑢 𝜕𝑣 𝜕𝑣 𝜕𝑢
= −2𝜌(𝑢𝑦 𝑣𝑥 + 𝑣𝑦 𝑣𝑦 ) = −2𝜌(𝑢𝑦 𝑣𝑥 − 𝑢𝑥 𝑣𝑦 ) = 2𝜌 ( − )
𝜕𝑥 𝜕𝑦 𝜕𝑥 𝜕𝑦
4. For a 2-D incompressible flow we can define the stream function 𝜙 by requiring
𝜕𝜙 𝜕𝜙
𝑢= ; 𝑣=−
𝜕𝑦 𝜕𝑥
We also can define a flow variable called vorticity
𝜕𝑣 𝜕𝑢
𝜔= −
𝜕𝑥 𝜕𝑦
Show that
𝜕2𝜙 𝜕2𝜙
𝜔 = − ( 2 + 2)
𝜕𝑥 𝜕𝑦
Solution:
𝜕𝑣 𝜕𝑢 𝜕 𝜕𝜙 𝜕 𝜕𝜙 𝜕2𝜙 𝜕2𝜙
𝜔= − = (− ) − ( ) = − ( 2 + 2)
𝜕𝑥 𝜕𝑦 𝜕𝑥 𝜕𝑥 𝜕𝑦 𝜕𝑦 𝜕𝑥 𝜕𝑦




Downloaded by: tutorsection | Want to earn $1.236
Distribution of this document is illegal extra per year?

, Chapter 2
𝑑𝜙
1. Develop a second-order accurate finite difference approximation for ( 𝑑𝑥 ) on a non-uniform
𝑖
mesh using information (𝜙 and 𝑥 values) from mesh points 𝑥𝑖−1 , 𝑥𝑖 and 𝑥𝑖+1. Suppose 𝛿𝑥𝑖 =
𝑥𝑖+1 − 𝑥𝑖 = 𝛼𝛿𝑥𝑖−1 = 𝛼(𝑥𝑖 − 𝑥𝑖−1 ).

Solution:

Assume close to the 𝑖 𝑡ℎ point, 𝜙(𝑥) = 𝜙𝑖 + 𝑏(𝑥 − 𝑥𝑖 ) + 𝑐(𝑥 − 𝑥𝑖 )2 + 𝑑(𝑥 − 𝑥𝑖 )3 …
𝑑𝜙 𝑑𝜙
Then 𝑑𝑥 = 𝑏 + 2𝑐(𝑥 − 𝑥𝑖 ) + ⋯ and ( 𝑑𝑥 ) = 𝑏.
𝑖

Now 𝜙𝑖+1 = 𝜙(𝑥𝑖+1 ) = 𝜙𝑖 + 𝑏(𝑥𝑖+1 − 𝑥𝑖 ) + 𝑐(𝑥𝑖+1 − 𝑥𝑖 )2 + ⋯ = 𝜙𝑖 + 𝑏Δ𝑥𝑖 + 𝑐Δ𝑥𝑖2 + 𝑑Δ𝑥𝑖3 …
2 3
And 𝜙𝑖−1 = 𝜙(𝑥𝑖−1 ) = 𝜙𝑖 + 𝑏(𝑥𝑖−1 − 𝑥𝑖 ) + 𝑐(𝑥𝑖−1 − 𝑥𝑖 )2 + ⋯ = 𝜙𝑖 − 𝑏Δ𝑥𝑖−1 + 𝑐Δ𝑥𝑖−1 − 𝑑Δ𝑥𝑖−1 …
2
So Δ𝑥𝑖−1 𝜙𝑖+1 − Δ𝑥𝑖2 𝜙𝑖−1 = (Δ𝑥𝑖−1
2
− Δ𝑥𝑖2 )𝜙𝑖 + 𝑏Δ𝑥𝑖 Δ𝑥𝑖−1 (Δ𝑥𝑖 + Δ𝑥𝑖−1 ) + 𝑑Δ𝑥𝑖2 Δ𝑥𝑖−1
2
(Δ𝑥𝑖 +
Δ𝑥𝑖−1 ) + ⋯
2
Δ𝑥𝑖−1 𝜙𝑖+1 −Δ𝑥𝑖2 𝜙𝑖−1 −(Δ𝑥𝑖−1
2
−Δ𝑥𝑖2 )𝜙𝑖
And 𝑏 = Δ𝑥𝑖 Δ𝑥𝑖−1 (Δ𝑥𝑖 +Δ𝑥𝑖−1 )
− 𝑑Δ𝑥𝑖 Δ𝑥𝑖−1 + ⋯

𝑑𝜙
A 2nd order finite difference for ( 𝑑𝑥 ) is therefore
𝑖
2
𝑑𝜙 Δ𝑥𝑖−1 𝜙𝑖+1 − Δ𝑥𝑖2 𝜙𝑖−1 − (Δ𝑥𝑖−1
2
− Δ𝑥𝑖2 )𝜙𝑖 𝜙𝑖+1 + (α2 − 1)𝜙𝑖 − α2 𝜙𝑖−1
( ) =𝑏≈ =
𝑑𝑥 𝑖 Δ𝑥𝑖 Δ𝑥𝑖−1 (Δ𝑥𝑖 + Δ𝑥𝑖−1 ) α(α + 1)Δ𝑥𝑖−1

2. Use the scheme you developed for problem 1 to evaluate the derivative of 𝜙(𝑥) =
sin(𝑥 − 𝑥𝑖 + 1) at point 𝑖. Suppose Δ𝑥𝑖−1 = 0.02 and Δ𝑥𝑖 = 0.01. Compare your
numerical result with the exact solution, which is cos(1). Then halve both Δ𝑥𝑖−1 and Δ𝑥𝑖 ,
and redo the calculation. Is the scheme truly second-order accurate?
Solution:
clear; clc;
dxi = 0.01; dxim1 = 0.02; alpha = dxi/dxim1;
for iter = 1 : 2
x = [-dxim1,0,dxi];
phi = sin(x+1);




Downloaded by: tutorsection | Want to earn $1.236
Distribution of this document is illegal extra per year?

, dphidx = (phi(3)+(alpha^2-1)*phi(2)-alpha^2*phi(1))/(alpha*(alpha+1)*dxim1);
err(iter) = dphidx-cos(1);
dxi = dxi/2; dxim1 = dxim1/2;
end
err(1)/err(2)


It is truly 2nd-order accurate.

3. Reproduce the calculation presented in Section 2.1.2 Example: Laminar Channel Flow.
Solution:
%--------------------------------------------------------------%
% FULLY-DEVELOPED CHANNEL FLOW %
% By George Qin for "A Course of Computational Fluid Dynamics" %
%--------------------------------------------------------------%
tic
clear; clc;
% PARAMETERS
H = 1; N = 5;
% FACE LOCATIONS
yf = linspace(0,H,N+1);
% NODE LOCATIONS
y = 0.5*(yf(1:end-1)+yf(2:end));
% DELTA Y
dy = yf(2)-yf(1);
% THE THREE DIAGONAL VECTORS IN THE COEFFICIENT MATRIX
as = -(1/dy^2)*ones(1,N);
ap = (2/dy^2)*ones(1,N);
an = -(1/dy^2)*ones(1,N);
b = ones(1,N);
% SPECIAL VALUES AT THE BOUNDARIES (BOUNDARY CONDITIONS)
%ap(1) = 3/dy^2; as(1) = 0;
%ap(1) = 4/dy^2; an(1) = - (4/3)/dy^2; as(1) = 0;
ap(1) = 3/dy^2; as(1) = 0; b(1) = 3/4;
ap(end) = 1/dy^2; an(end) = 0;
% SOLVE THE SYSTEM OF LINEAR EQUATIONS WITH TDMA ALGORITHM
u = TDMA(as,ap,an,b); % this is in the appendix of the textbook
toc
% COMPARE WITH EXACT SOLUTION
u_exact = y.*(1-0.5*y);
err = (u_exact-u)./u_exact * 100;
% RECALCULATE EXACT SOLUTION VECTOR FOR PLOTTING
y_exact = 0:0.01:1;
u_exact = y_exact.*(1-0.5*y_exact);
plot(y_exact,u_exact,y,u,'ro')
xlabel('$y$','FontSize',20,'Interpreter','Latex')
ylabel('$u$','FontSize',20,'Interpreter','Latex')
h_legend=legend('Exact Solution','Numerical Solution');
set(h_legend,'FontSize',14,'Interpreter','Latex','Location','NorthWest')


4. Show that the method used in Section 2.1.2 Example: Laminar Channel Flow is first-
order accurate by using the global error estimate technique.
Solution:

These finite difference equations and their truncation errors are reproduced here:




Downloaded by: tutorsection | Want to earn $1.236
Distribution of this document is illegal extra per year?