Simulation of a 1D Super-sonic nozzle flow simulation using Macormack Method

Objective: Simulation of 1D supersonic nozzle flow simulation using the Macormack Method using Non-conservative and conservative forms of equations.

Non-conservative Flow:

In this non-conservative flow, a small control volume is taken inside the pipe flow but with a small change, the control volume here is moving with the instead of being fixed as in case of conservative flow.

McCormack method is used  for solving the governing equations. This concept is used in programming to give the exact solution:

(A)  Governing Equations:

(i) Continuity Equation:

{delrho'}/{delt'} = -rho'*{delv'}/{delx'} - rho'*v'*{del(ln A')}/{delx'} - v'*{delrho'}/{delx'}

(ii) Momentum Equation:

{delv'}/{delt'} = - v'*{delv'}/{delx'} -1/gamma*({delT'}/{delx'} + {T'}/{rho'}*{delrho'}/{delx'})

(iii) Energy Equation:

{delT'}/{delt'} = - v'*{delT'}/{delx'} -(gamma-1)*T'*({delv'}/{delx'} + v'*{del(lnA')}/{delx'})

(B)  Initial conditions of Profile at time, t =0:

Density, rho = 1-0.3146*x;

Temperature, T = 1- 0.2314*x

velocity, v = (0.1 + 1.09*x)*t^0.5

Area,A = 1 + 2.2*(x-1.5)^2

(C) Boundary conditions

Inlet:

v(1) = 2v(2) - v(3)

Outlet:

v(n) = 2v(n-1)-v(n-2)

rho(n) = 2rho(n-1)-rho(n-2)

T(n) = 2T(n-1)-T(n-2)

conservative Flow:

In this conservative flow, when a fluid is flowing through a pipe a small control volume is taken in some area inside that pipe and this control volume is fixed, so that the equation can be solved find out the properties of the flow through that point.

A) Governing Equations:

(i) Continuity Equation:

{delU_1}/{delt'} = -{delF_1}/{delx'}

(ii) Momentum Equation:

{delU_2}/{delt'} = -{delF_2}/{delx'} + J_2

(iii) Energy Equation:

{delU_3}/{delt'} = -{delF_3}/{delx'}

Solution vector:U

Flux vector:F

Source term:J_2

U_1= rho'A'

U_2=rho'A'v'

U_3=rho'({e'}/{gamma-1} + gamma/2v'^2)A'

F_1= rho'A'v'

F_2= rho'A'v' +(1/gamma)p'A'

F_3=rho'({e'}/{gamma-1} + gamma/2v'^2)v'A' +p'A'v'

J_2 = 1/gamma*p'*{dA'}/{dx'}

(B)  Initial conditions of Profile at time, t=0:

x(i)>=0 && x(i)<=0.5

rho(i) = 1

T(i) = 1

x(i)>=0.5 && x(i)<=1.5

rho(i) = 1-0.366*(x(i)-0.5);

T(i) = 1-0.167(x(i)-0.5);

x(i)>=1.5 && x(i)<=3.5

rho(i) = 0.634-0.3879(x(i)-1.5)

T(i) = 0.833-0.3507(x(i)-1.5)

Area,A = 1 + 2.2*(x-1.5)^2

Assuming steady state flow

v = 0.59/(rho*A)

(C) Boundary conditions

Inlet:

U_1(1)=rho(1)*A(1)

U_2(1)= 2U_2(2)-U_2(3)

v(1)= {U_2(1)}/{U_1(1)}

T(1) = 1

U_3(1) = U_1(1)*(({T(1)}/(gamma-1))+(gamma/2)*(v(1)^2))

Outlet:

U_1(n) = 2U_1(n-1) - U_1(n-2)

U_2(n) = 2U_2(n-1) - U_2(n-2)

U_3(n) = 2U_3(n-1) - U_3(n-2)

The flow field are Non-Dimensionalized, where the flow variable are referenced to their reservoir values.

Non-Dimensional Temperature:T'=T/T_o where T_o= reservoir Temperature

Non-Dimensional Density:rho' = rho/rho_owhere rho_o= reservoir Density

Non-Dimensional length:x'=x/L where L = length of nozzle

Non-Dimensional velocity:v'= v/a_o where a_o = speed of sound

Non-Dimensional Time:t'=t/(L/a_o)

Non-Dimensional AreaA'=A/(A^*) where A^*= Throat area

Calling a Non-conservative function:

function [v,rho,t,m,p,mach,v_t,rho_t,t_t,m_t,p_t,mach_t,k,time_execution_1] = non_conservative(n,x,dx,c,gamma,a,nt,throat)

tic
% Inputs
% Length of the domain
L = 3
% no of the grid points
n = 31
% x array along the length of the domain
x =linspace(0,L,n)
% Grid size
dx = x(2) - x(1)
% other parameters
gamma = 1.4
% calculate Initial profile, which are non dimensional
rho = 1-0.3146*x;
t = 1- 0.2314*x;   % t =temperature
v = (0.1 + 1.09*x).*t.^0.5;
a = 1 + 2.2*(x-1.5).^2; % Area, a

% Time steps
nt = 1600
% calculating the value of time step using courant criteria
for i = 1:n
del_t(i)= c*(dx/(t(i)^0.5+v(i)));
end
dt = min(del_t(i));

% outer time loop
for  k = 1:nt
rho_old = rho;
v_old = v;
t_old = t;
a_old = a;

% predictor method
for j = 2:n-1

dvdx = (v(j+1)-v(j))/dx;
drhodx = (rho(j+1)-rho(j))/dx;
dtdx = (t(j+1)-t(j))/dx;
dlogadx = (log(a(j+1)) - log(a(j)))/dx;

% continuity equation
drhodt_p(j) = -rho(j)*dvdx - rho(j)*v(j)*dlogadx - v(j)*drhodx;

% Momentum equation
dvdt_p(j) = -v(j)*dvdx - (1/gamma)*(dtdx + (t(j)/rho(j))*drhodx);

% Energy equation
dtdt_p(j) = -v(j)*dtdx - (gamma - 1)*t(j)*(dvdx + v(j)*dlogadx);

% solution update
rho(j) = rho(j) + drhodt_p(j)*dt;
v(j) = v(j) + dvdt_p(j)*dt;
t(j) = t(j) + dtdt_p(j)*dt;

end

% corrector method

for j = 2:n-1

dvdx = (v(j)-v(j-1))/dx;
drhodx = (rho(j)-rho(j-1))/dx;
dtdx = (t(j)-t(j-1))/dx;
dlogadx = (log(a(j)) - log(a(j-1)))/dx;

% continuity equation
drhodt_c(j) = -rho(j)*dvdx - rho(j)*v(j)*dlogadx - v(j)*drhodx;

% Momentum equation
dvdt_c(j) = -v(j)*dvdx - (1/gamma)*(dtdx + (t(j)/rho(j))*drhodx);

% Energy equation
dtdt_c(j) = -v(j)*dtdx - (gamma - 1)*t(j)*(dvdx + v(j)*dlogadx);
end

% compute the average time derivative
drhodt = 0.5*(drhodt_p + drhodt_c)
dvdt = 0.5*(dvdt_p + dvdt_c)
dtdt = 0.5*(dtdt_p + dtdt_c)

% Final solution update
for i = 2:n-1
rho(i) =rho_old(i) + drhodt(i)*dt;
v(i) = v_old(i) + dvdt(i)*dt;
t(i) = t_old(i) + dtdt(i)*dt;
end

% Apply the Boundary conditions

% Inlet
v(1) = 2*v(2) - v(3)

% outlet
v(n) = 2*v(n-1)-v(n-2);
rho(n) = 2*rho(n-1)-rho(n-2);
t(n) = 2*t(n-1)-t(n-2);

% Calculating solutions
% Mass flow rate
m = rho.*a.*v;
% pressure
p = rho.*t;
% Mach number
mach = v./((t).^0.5);

% Throat values
rho_t(k) = rho(throat);
v_t(k)=v(throat);
t_t(k)= t(throat);
p_t(k)=p(throat);
m_t(k)=m(throat);
mach_t(k)= mach(throat);

% plotting mass flow rate at different time steps
figure(1);
if k == 50
plot(x,m,'b','linewidth',2)
hold on
else
if k == 100
plot(x,m,'r','linewidth',2)
else
if k == 200
plot(x,m,'g','linewidth',2)
else
if k == 400
plot(x,m,'y','linewidth',2)
else
if k == 800
plot(x,m,'c','linewidth',2)
else
if k == 1600
plot(x,m,'m','linewidth',2)
end
end
end
end
end
end
grid on

xlabel('Non-Dimensional distance through nozzle(x)','fontweight','bold','fontsize',12)
ylabel('Non-Dimensional mass flow rate through nozzle(m)','fontweight','bold','fontsize',12)
title('Non-Dimensional mass flow rate(Non-conservative form) along x at different time steps','fontweight','bold','fontsize',12,'color','r')
legend({'50 time steps','100 time steps','200 time steps','400 time steps','800 time steps','1600 time steps'})

time_execution_1=toc
end

Calling a conservative function:

function [v,rho,t,m,p,mach,v_t,rho_t,t_t,m_t,p_t,mach_t,k,time_execution_2] = conservative(n,x,dx,c,gamma,a,nt,throat)

tic
% Inputs
% Length of the domain
L = 3
% no of the grid points
n = 31
% x array along the length of the domain
x =linspace(0,L,n)
% Grid size
dx = x(2) - x(1)
% other parameters
gamma = 1.4
% calculate Initial profile, which are non dimensional
for i = 1:n
if x(i)>=0 && x(i)<=0.5
rho(i) = 1;
t(i) = 1;
else
if x(i)>=0.5 && x(i)<=1.5
rho(i) = 1-0.366*(x(i)-0.5);
t(i) = 1-0.167*(x(i)-0.5);
else
if x(i)>=1.5 && x(i)<=3.5
rho(i) = 0.634-0.3879*(x(i)-1.5);
t(i) = 0.833-0.3507*(x(i)-1.5);
end
end
end
end

% Assuming steady state flow for the initial profile
v = 0.59./(rho.*a);

% Time steps
nt = 1600

% calculating the value of time step using courant criteria
for i = 1:n
del_t(i)= c*(dx/(t(i)^0.5+v(i)));
end
dt = min(del_t(i));

% solution vector parameters
p = rho.*t;
U1 = rho.*a;
U2 = rho.*a.*v;
U3 = rho.*((t./(gamma-1))+((gamma/2).*(v.^2))).*a

% outer Time loop
for k = 1:nt;
U1_old = U1;
U2_old = U2;
U3_old = U3;
% Defining flux variables
F1 = U2;
F2 = ((U2.^2)./U1) + ((gamma-1)/gamma)*(U3-(((gamma/2)*U2.^2)./U1));
F3 = ((gamma*U2.*U3)./U1) - (gamma*(gamma-1)/2)*((U2.^3)./(U1.^2));

% predictor step
for j=2:n-1

dF1dx = (F1(j+1)-F1(j))/dx;
dF2dx = (F2(j+1)-F2(j))/dx;
dF3dx = (F3(j+1)-F3(j))/dx;

% continuity equation
dU1dt_p(j) = -dF1dx
% Momentum equation
dU2dt_p(j) = -dF2dx + J2
% Energy equation
dU3dt_p(j) = -dF3dx

% solution update
U1(j) = U1(j) +  dU1dt_p(j)*dt;
U2(j) = U2(j) +  dU2dt_p(j)*dt;
U3(j) = U3(j) +  dU3dt_p(j)*dt;
end
% Decoding primitive terms to initialize the corrector step
rho = U1./a;
v = U2./U1;
t = (gamma-1)*((U3./U1)-((gamma/2)*((U2./U1).^2)));
p = rho.*t;
F1 = U2;
F2 = ((U2.^2)./U1) + ((gamma-1)/gamma)*(U3-(((gamma/2)*U2.^2)./U1));
F3 = ((gamma*U2.*U3)./U1)-(gamma*(gamma-1)/2)*((U2.^3)./(U1.^2))
% corrector step
for j=2:n-1
dF1dx = (F1(j)-F1(j-1))/dx;
dF2dx = (F2(j)-F2(j-1))/dx;
dF3dx = (F3(j)-F3(j-1))/dx;

% continuity equation
dU1dt_c(j) = -dF1dx
% Momentum equation
dU2dt_c(j) = -dF2dx + J2
% Energy equation
dU3dt_c(j) = -dF3dx
end

% Average derivatives
dU1dt = 0.5*(dU1dt_p + dU1dt_c);
dU2dt = 0.5*(dU2dt_p + dU2dt_c);
dU3dt = 0.5*(dU3dt_p + dU3dt_c);

% Final solution update
for i = 2:n-1
U1(i) = U1_old(i) +  dU1dt(i)*dt;
U2(i) = U2_old(i) +  dU2dt(i)*dt;
U3(i) = U3_old(i) +  dU3dt(i)*dt;
end

% Apply the Boundary conditions
% Inlet
U1(1)=rho(1)*a(1);
U2(1)= 2*U2(2)-U2(3);
v(1)= U2(1)/U1(1);
t(1) = 1;
U3(1) = U1(1)*((t(1)/(gamma-1))+((gamma/2)*(v(1)^2)));
% outlet
U1(n) = 2*U1(n-1) - U1(n-2);
U2(n) = 2*U2(n-1) - U2(n-2);
U3(n) = 2*U3(n-1) - U3(n-2);

rho = U1./a;
v = U2./U1;
t = (gamma-1)*((U3./U1) - ((gamma/2)*((U2./U1).^2)))

% caculating solutions
% Mass flow rate
m = rho.*v.*a;
% pressure
p = rho.*t;
% Mach number
mach = v./((t).^0.5);
% Throat values
rho_t(k) = rho(throat);
v_t(k)=v(throat);
t_t(k)= t(throat);
p_t(k)=p(throat);
m_t(k)=m(throat);
mach_t(k)= mach(throat);
% plotting mass flow rate at different time steps
figure(2);
if k == 50
plot(x,m,'b','linewidth',2)
hold on
else
if k == 100
plot(x,m,'r','linewidth',2)
else
if k == 200
plot(x,m,'g','linewidth',2)
else
if k == 400
plot(x,m,'y','linewidth',2)
else
if k == 800
plot(x,m,'c','linewidth',2)
else
if k == 1600
plot(x,m,'m','linewidth',2)
end
end
end
end
end
end
grid on

xlabel('Non-Dimensional distance through nozzle(x)','fontweight','bold','fontsize',12)
ylabel('Non-Dimensional mass flow rate through nozzle(m)','fontweight','bold','fontsize',12)
title('Non-Dimensional mass flow rate(conservative form) along x at different time steps','fontweight','bold','fontsize',12,'color','r')
legend({'50 time steps','100 time steps','200 time steps','400 time steps','800 time steps','1600 time steps'})

time_execution_2=toc
end

Main Matlab Program:

clear all
close all
clc

% Inputs
% Length of the domain
L = 3;
% Number of the grid points
n = 31;
% x array along the length of the domain
x = linspace(0,L,n);
% Grid spacing
dx = x(2) - x(1);
% Area,a
a = 1+2.2*(x-1.5).^2;
gamma = 1.4
% Exact node at the throat, numbers of nodes is an odd number
throat = ((n-1)/2)+1;

% courant number
c = 0.5

% No. of time steps
nt = 1600;

% non conservative form
[v1,rho1,t1,m1,p1,mach1,v_t1,rho_t1,t_t1,m_t1,p_t1,mach_t1,k1,time_execution_1]=non_conservative(n,x,dx,c,gamma,a,nt,throat);
% conservative form
[v2,rho2,t2,m2,p2,mach2,v_t2,rho_t2,t_t2,m_t2,p_t2,mach_t2,k2,time_execution_2]=conservative(n,x,dx,c,gamma,a,nt,throat);

% plots comparing nonconservative and conservative
figure(3)
s1=subplot(4,1,1);
plot(linspace(1,nt,nt), rho_t1);
hold on
plot(linspace(1,nt,nt), rho_t2);
hold on
plot(linspace(nt-600,nt+200,800),linspace(0.634,0.634,800),'color','g');
grid on
ylabel("rho'",'fontsize',12,'fontweight','bold')
h=legend(s1,{'Non-conservative Form', 'conservative Form','Analytical solution'})
set(h, 'Location','northeastoutside')
title('Timewise variation of flow parameters at throat','fontsize',12,'fontweight','bold')

s2=subplot(4,1,2)
plot(linspace(1,nt,nt), t_t1);
hold on
plot(linspace(1,nt,nt), t_t2);
hold on
plot(linspace(nt-600,nt+200,800),linspace(0.833,0.833,800),'color','g');
grid on
ylabel("T'",'fontsize',12,'fontweight','bold')
h=legend(s2,{'Non-conservative Form', 'conservative Form','Analytical solution'})
set(h, 'Location','northeastoutside')

s3=subplot(4,1,3)
plot(linspace(1,nt,nt), p_t1);
hold on
plot(linspace(1,nt,nt), p_t2);
hold on
plot(linspace(nt-600,nt+200,800),linspace(0.528,0.528,800),'color','g');
grid on
ylabel("pr'",'fontsize',12,'fontweight','bold')
h=legend(s3,{'Non-conservative Form', 'conservative Form','Analytical solution'})
set(h, 'Location','northeastoutside')
legend('Non-conservative Form', 'conservative Form','Analytical solution')

s4=subplot(4,1,4)
plot(linspace(1,nt,nt), mach_t1);
hold on
plot(linspace(1,nt,nt), mach_t2);
hold on
plot(linspace(nt-600,nt+200,800),linspace(1,1,800),'color','g');
grid on
xlabel("Number of time steps",'fontsize',12,'fontweight','bold')
ylabel("Ma'",'fontsize',12,'fontweight','bold')
h=legend(s4,{'Non-conservative Form', 'conservative Form','Analytical solution'})
set(h, 'Location','northeastoutside')

% variation in flow parameters along x axis
figure(4)
q1=subplot(4,1,1)
plot(x,rho1,'r')
hold on
plot(x,rho2,'b')
grid on
ylabel("rho'",'fontsize',12,'fontweight','bold')
h = legend(q1,{'Non conservative Form','conservative Form'});
set(h,'location','northeastoutside')
title('Variation of flow parameters along x axis','fontsize',12,'fontweight','bold')

q2=subplot(4,1,2)
plot(x,t1,'r')
hold on
plot(x,t2,'b')
grid on
ylabel("T'",'fontsize',12,'fontweight','bold')
h = legend(q2,{'Non conservative Form','conservative Form'});
set(h,'location','northeastoutside')

q3=subplot(4,1,3)
plot(x,p1,'r')
hold on
plot(x,p2,'b')
grid on
ylabel("pr'",'fontsize',12,'fontweight','bold')
h = legend(q3,{'Non conservative Form','conservative Form'});
set(h,'location','northeastoutside')

q4=subplot(4,1,4)
plot(x,mach1,'r')
hold on
plot(x,mach2,'b')
grid on
xlabel('Non-Dimensional Length','fontsize',12,'fontweight','bold')
ylabel("Ma'",'fontsize',12,'fontweight','bold')
h = legend(q4,{'Non conservative Form','conservative Form'});
set(h,'location','northeastoutside')

Results and conclusion:

Plot of Variation of flow parameters along x-axis: From the plots it is clear that the both the conservative and non-conservative form produce same plot after steady state, 1600 time steps are carried out by time marching approach, at which the steady state is almost obtained and the values of primitive variables are good in approximation with exact analytical solution.

Plot of Time variation of flow parameters at throat: •  From the plot it is clear that the non-conservative form lines are closer to the exact values of rho',T',pr', and Mach' when compared to conservative form of solution
• The Residue (difference between the numerical and exact solution) exists after 1600 time steps.
• Non-conservative form leads to smaller residue, and the quality index is determined by the rate of decay of the residue. So, non-conservative form obtains a better quality index when compared to the conservative form.
• From the plots, it is clear that the steady state is obtained with a good level of approximation after 400 time steps.

Plot of Non-Dimensional mass flow rate(Non-conservative form): • The mass flowrate obtained for the last three time steps (i.e., 400, 800, 1600) are almost overlying in the same region of space
• This plot substantiates our previous assumption from the observation of "Time variation of flow parameters at throat", i.e., after 400 time steps almost steady state is obtained .

Plot of Non-Dimensional mass flow rate(conservative form): • Even for this conservative form the mass flowrate obtained for the last three time steps (i.e., 400, 800, 1600) are almost overlying in the same region of space.
• The Initial conditions for the conservative form and non conservative form are different so the comparison is inappropriate.
• The steady state is attained after infinite steps, but from the plot after 400 the variation is too small with the exact solution.

Grid dependence:

• If the number of grid points is doubled, and there is only a marginal difference in the solution, then the calculations with the initial grid points are independent.
• In this case n =31, when modifying the value of 'n' the number of time steps should also be modified, so that the elapsed non dimensional time is same for both cases.
•  Grid independence test can also be carried out by comparing the values obtained by conservative and non-conservative forms at the throat, for n=31 and 61. It has to be noted that the non-conservative form results are more accurate when compared with the conservative form results.
• After doubling the number of grid points, the solutions of the conservative form are slightly closer to the exact value, but still attained the same accuracy of non-conservative form with half as many grid points i.e., n = 31. Therefore for practical purposes, n =31 can be treated as grid-independent.
 Flow parameters Non conservative form Conservative form Exact value n = 31 n = 61 n = 31 n = 61 rho' 0.639 0.638 0.652 0.639 0.634 t' 0.836 0.835 0.841 0.836 0.833 Pr' 0.534 0.532 0.548 0.533 0.528 mach' 0.999 1.0 0.978 0.998 1.00

WHICH METHOD IS FASTER?

For n = 31, c = 0.5

Non-Conservative simulation Time = 93.448

Conservative simulation Time = 171.8982

For n = 61, c = 0.5

Non-Conservative simulation Time = 170.474

Conservative simulation Time = 737.3677

Parametric Study of Gate Valve Shravankumar Nagapuri · 2019-10-22 03:14:10

Objectives: To perform a parametric study on the gate valve simulation. To obtain the mass flow rates at the outlet for 5 design points. To show the cut section view for different design points, and also to show the gate disc lift and fluid volume extraction. To s Read more

Symmetry vs Wedge vs HP equation Shravankumar Nagapuri · 2019-10-09 22:50:52

Objectives: To Write a Matlab program that can generate the computational blockMesh file  automatically for a wedge angle of 3 degrees and compare the results obtained for Symmetry BC and Wedge BC. To Write a Matlab program that takes an angle as input and gene Read more

Combustion simulation on the combustor model Shravankumar Nagapuri · 2019-10-05 16:23:22

1Q. Briefly explain about the possible types of combustion simulations in FLUENT. A:  Combustion or burning is a high temperature exothermic chemical reaction between a fuel and an oxidant accompanied by the production of heat, light and unburnt gases in the fo Read more

Cyclone separator challenge Shravankumar Nagapuri · 2019-10-01 17:50:02

Cyclone Separator: Principle and Working: A high speed rotating (air)flow is established within a cylindrical or conical container called a cyclone. Air flows in a spiral pattern, beginning at the top (wide end) of the cyclone and ending at the bottom (narrow) end bef Read more

Gearbox sloshing effect Shravankumar Nagapuri · 2019-09-28 15:37:58

Objective: To analyse the flow pattern of the fluid inside the gearbox, for two different clearances of the same geometry and each geometry flow is analyzed by two different fluids (Water and Oil) as lubricants. Gearbox sloshing effect:  Slosh refers Read more

Simulation of Flow through a pipe in OpenFoam Shravankumar Nagapuri · 2019-09-26 19:10:15

Objectives: To write a program in Matlab that can generate the computational mesh automatically for any wedge angle and Grading scheme. To calculate, length of the pipe using the entry length formula for laminar flow through a pipe To show that entry length is suf Read more

Conjugate Heat Transfer CHT Analysis on a graphics card Shravankumar Nagapuri · 2019-09-20 16:08:56

Objectives: To find out the maximum temperature attained by the processor. To Prove that the simulation has achieved convergence with appropriate images and plots. To Find out the heat transfer coefficient at appropriate areas of the model. To Identify potential h Read more

Finite Volume Method FVM Literature Review Shravankumar Nagapuri · 2019-09-17 13:14:28

Finite volume Method: The finite volume method (FVM) is a method for representing and evaluating partial differential equations in the form of algebraic equations. Similar to the finite difference method (FDM) or finite element method (FEM), values are calculated at Read more

Rayleigh Taylor instability Shravankumar Nagapuri · 2019-09-14 17:48:30

Q1. What are some practical CFD models that have been based on the mathematical analysis of Rayleigh Taylor waves? In your own words, explain how these mathematical models have been adapted for CFD calculations  Kelvin–Helmholtz instability: The Kelvin Read more

BlockMesh Drill down challenge Shravankumar Nagapuri · 2019-09-11 18:17:33

Objectives: To Generate the BlockMesh file for the given Geometry using OpenFOAM and use the icoFoam solver to simulate the flow through a backward-facing step. To create multiple meshes and will be comparing the results obtained from each mesh. To show the velocit Read more