Week 5: Prandtl Meyer Shock problem

AIM :- To solve Prandtl Meyer Shock Problem.

OBJECTIVES :-

1. To understand the theory of formation of expansion fan due to infinite mach waves, originated at the sharp convex corner,
2. To understand the reason for the usage of Neumann BC at the outlet for supersonic flow,
3. To understand the effect of SGS (i.e. sub grid size) parameter on shock location,
4. To understand the effect of using subsonic Initial condition for supersonic flow.

THEORY AND INTRODUCTIONS :-

SHOCK FLOW BOUNDARY CONDITIONS :-

Boundary condition in fluid dynamics are the set of constraints to boundary value problems in computational fluid dynamics. These boundary conditions include inlet boundary conditions, outlet boundary conditions, wall boundary conditions, constant pressure boundary conditions, axisymmetric boundary conditions, symmetric bounday conditions and periodic or cyclic boundary conditions.

Transient problem require one more thing i.e. initial conditions where initial values of flow variable are specified at nodes in the flow domain. Various types of boundary conditions are used in CFD for different conditions and purposes and are discussed as follows :-

1. Inlet Boundary conditions :-

In inlet BC's the distribution of all flow variables needs tobe specified at inlet boundaries mainely flow velocity. This type of boundary conditions are commen and specified mostly where inlet flow velocity is known.

2. Outlet Boundary conditions :-

In outlet BC's, the distribution of all flow variable needs tobe specified, mainly flow velocity. This can be thought as a conjunction to inlet boundary conditions. This type of boundary conditions is commen and specified mostly where outlet velocity is known. The flow attains a fully developed state where no change occurs in the flow direction when the outlet is selected far away from the geometrical disturbances. In such region, an outlet could be outlined and the gradient of all variables could be equated to zero in the flow direction except pressure.

3. Wall Boundary conditions :-

The most commen boundary that comes upon in the confined fluid flow problems is the wall. This is also commenly known as no slip boundary conditions and this is the appropriate conditions for velocity components at the wall. The normal component could be set to zero straight away while the tangential component is set to the velocity of the wall.

Vnormal = 0,

Vtangential = Vwall

Heat transfer through the wall can be specified or if the walls are considered adiabatic, then heat transfer across the wall is set to zero.

Qadiabatic walls = 0.

4. Constant Pressure boundary condtions :-

This type of boundary conditions is used where boundary values of pressure are known and the exact details of the flow distribution are unknown. This includes pressure inlet and outlet conditions mainly. Typical example that utilize this boundary condition include buoyancy driven flows, internal flows with multiple outlets, free surface flows and external flows around objects. An example is flow outlet into atmosphere where pressure is atmospheric.

5. Axisymmetric Boundary conditions :-

In this BC's, the model is axisymmetric with respect to the main axis such that at the particular r=R, all theta's and each z= Z-slice, each flow variable has the same value. A good example is the flow in a circuilar pipe where the flow and pipe axes coincide.

Vr(R,theta,Z) = Constant

(r = R,theta,Z).

6. Symmetric Boundary conditions :-

In this boundary conditions, it is assumed that on the two sides of the boundary, same physical processes exist. All the variable have same value and gradient at the same distance from the boundary. It acts as a mirror that reflect all the flow distribution to the other side. The conditions at symmetric boundary are no flow across boundary and no scalar flux across boundary.

A good example is of a pipe flow with a symmetric obstacle in the flow. The obstacle divides the upper flow and lower flow as mirrored flow.

7. Periodic or Cyclic Boundary conditions :-

In periodic or cyclic boundary conditions arises from a different type of symmetry in a problem. If a component has a repeated pattern in flow distribution more than twice, thus violating the mirror image requirments required for symmetric boundary conditions. A good example would be swept vane pump. Where the mark area is repeated four times in r-theta co-ordinates. The cyclic symmetric areas should have the same flow variables and distribution and should satisfy that in every Z-slice.

In case of shock wave problems some typical boundary conditions are used as follows-

• Inlet boundary condition
• Outlet boundary condition
• Wall boundary condition
• Dirichlet Boundary conditon
• Neumann Boundary condition
• Mixed/Cauchy/Robin Boundary conditions

In Dirichlet BC, values of dependent variables are provided along the boundary of the domain. Also called as a fixed boundary conditions.

In Neumann BC, the value of the gradient of the dependent variable (flux) is provided along the boundary of the domain.

In Mixed/Cauchy/Robin BC, are all different types of combinations of the Neumann and Dirichlet Boundary conditions.

# SHOCK WAVE :-

Definitions :- 1) A large amplitude wave formed by the sudden compression of the medium through which the wave moves. Shock wave can be caused by explosions or by objects moving through a fluid at a speed greater than the speed of sound.

2) A region of abrupt change of pressure and density moving as a wave front at or above the velocity of sound, caused by an intense explosion or supersonic flow over a body.

# PHYSICS BEHIND SHOCK WAVES :-

When the speed of a source equals the speed of a sound (v = c) the wave fronts can not escape the source. The resulting pile of waves forms a large amplitude "sound barrier" that makes sustained flight at this speed difficult and risky.

The term "sound barrier" first came into use during world war-2. Fighter pilots engaged in high speed dives noticed several irregularities as flying speeds approach the speed of sound aerodynamic drag increased markedly, much more than normally associated with increased speed, while lift and maneuverability decreased in a similarly unusual manner. Pilots at the time mistakenly thought that these effects meant that supersonic flight was impossible, that somehow airoplanes would never travel faster than the speed of sound. They were wrong when the speed of a source exceeds the speed of sound (v>c) the wave fronts lag behind the source in a cone shaped region with the source at the vertex. The edge of the cone forms a supersonic wave front with an unusually large amplitude called a "shock wave". When a shock wave reaches an observer a "sonic boom" is heard.

Unlike ordinary sound waves, the speed of a shock waves varies with its amplitude. The speed of a shock wave is always greater than the speed of a sound in the fluid and decreases as the amplitude of the wave decreases. When the shock wave speed equals the normal speeds, the shock wave dies and is reduced to an ordinary sound wave.

The ratio of the speed of a moving object (v) to the speed of sound (c) in a fluid is known as the Mach number (Ma) in honor of Earnst Mach (1838-1916), the moravian physicist, psychologist and phylosopher who studied sound and ballistics

Ma = v/c

The Mach number is a dimensionless measure of speed commen in aerodynamics. Mach 0.5 is half the speed of sound, Mach 2 is twice the speed of sound, and so on. Speed less than the speed of sound have a Mach number between zero and one and are described as subsonic. Those greater than the speed of sound have Mach number greater than one are described as supersonic. Speeds approximately equal to the speed of sound have Mach numbers approximately equal to one and are described as transonic.

In my words,

Shocks are generated in supersonic flows, where the information about the future disturbances can reach the domain only through the packets of subsonic regions in the domain i.e. near boundary layer etc. By using Slip BC we further constrained the possibility of information flow about the future disturbance. So the flowing fluid in the domain can't able to predict the future disturbance, and so it cant make considerable change in its proparties to match with the proparties of the fluid at the region of disturbance. Therefore it encounters a sudden change in it's proparties while hitting the region of disturbance, which we observe as shocks.

In our case, at the convex corner infinite numer of points disturbances creates infinite number of Mach waves which fan out together and we represent it with the name "Expansion Fan".

REASON BEHIND THE BOUNDARY CONDITION USED FOR SHOCK FLOW PROBLEMS :-

# Physical Reason :-

As the future information can not be carried nside the domain (information travels as pressure waves, and the maxmum speed at which they can travel is the speed of sound) the proparties of the fluid in the domain can't be alterd to match the flow with the required outlet BC. So usage of Neumann BC indirectly says that the flow is completely dependent on the initial conditions of the domain.

# Mathmatical reason :-

As all the charecterastics lines travel out of the domain from the outlet boundary, the only way of defining the condition at the outlet BC is through Neumann.

Temperature based Adaptive Mesh refinement is used, with three different values of sub grid scale. Density based steady state solver is used.

# EFFECT OF SGS (i.e. SUB GRID SCALING) PARAMETER ON SHOCK LOCATION :-

As the sub grid scaling increases we get greater Mach number, SGS value increases means movement of flame increases.

Filtering enlarges the region where shock capturing effects are present. This is an intrinsic handicap of LES when shock capturing is used, shock capturing schemes necessarily introduce deccipation to smear out the shock and be able to capture it. Additional filter operations, such as those required in the test filtering of dynamic models will further enlarge the region affected by the shock capturing methods. This enlarge region will have a negative effect on the physical fidelity of resolved field, in particular immedietly downstream of the shock which in turn affects the computation of subgrid scale terms, where the focus of the model development in the present work is. For a quantification and modelling approach to estimate the error resulting from the damping of turbulance introduce by shock capturing scheme in DNS, plus some remarks on its influence in LES.

A careful treatment of those region perhapps turning off subgrid scale model at grid points where shock capturing is active and using skewed filters nearly might be required. Solution have been proposed for related flows, such as the use of deconvolution without shock capturing scheme or adaptive shock confining filters.

In our case,

The inlet condition was a flow at Mach number = 2 at 286.1k. The convex corner angle is 5.352 degrees. After the expansion the outlet mach number is 2.2 which matches with analytical calculations.

AMR is switch on with SGS temperature as the criteria for refinement.

GEOMETRY WITH SCALING :-

This is the CAD model, which is drag and drop here from skill-lync page. I scaling it for propar manner to take a gridding scale and sub grid size. When we work on this geometry one things always keep in mind that the normal of the 2D boundary always along Z-directions, which are inside or outside no matter. But, in this case normal inside the body because fluid present inside the body.

I created here two casses for different velocities as follows-

CASE-1 :-

For velocity = 680 m/s

Sub grid scaling for temperature = 0.05k

The set-up for the given challenge are as follows-

1. GAS SIMULATION :-

Here, I use Redlich-Kwong state of equation which is very popular and according to that there have critical tmperature and pressure as shown above.

2. GLOBAL TRANSPORT PROPARTIES :-

Here, have the value for Turbulent Prandtl number and Turbulent Schmidt number by default.

3. RUN PARAMETERS FOR SOLVER :-

Here I use steady state solver for full hydrodynamic simulation mode. Steady state solver because of there have no change with espect to time.

4. RUN PARAMETERS FOR MISC. :-

Here I use steady state monitor and uncheck the shared memory for getting contineous simulations.

5. STEADY STATE MONITOR :-

Here I consider average velocity for outlet boundary. For steady state I use per Cycle unit. For mass average flow I use following data-

Monitor start delay = 200 cyc

Sample size = 100 cyc

Minimum sample size = 0.0001 sec.

6. SIMULATION TIME PARAMETER :-

This is the most important steps for speeding the simulations. Here we can increase more time step for resulting the speed of the simulation increses. Some time no more compramise with the accuracy because the solution is most important to verify.

7. SOLVER PARAMETER FOR STEADY STATE SOLVER CONTROL :-

Here only one things keep in mind that maximum convection CFL limit is less than or equal to the simulation time parameter's, convection CFL limit.

8. SOLVER PARAMETER FOR NAVIER STROKE SOLVER :-

Here I use density based Navier strokes solver because of there have large velocity, and for large velocity pressure based solver not work proparly. As velocity increases, density decreases and vice-versa.

9. BOUNDARY NAMES :-

Here, I write boundary names and region name same as that of region and initialization.

10. INLET BOUNDARY :-

This is the inlet boundary, where pressure is zero normal gradient (NE) and velocity is 680 m/s also temperature is 300k and by adding air as a species.

11. OUTLET BOUNDARY :-

This is the outlet boundary and for this pressure and velocity boundary condition is Zero normal gradient (NE). Also by adding air as a fluid.

12. TOP AND BOTTOM WALLS :-

This is the Top and bottom wall boundary by which given following conditions-

Wall motion type -> Stationary

Surface movement -> Fixed,

UDF -> Slip (due to fluid flows along X-axis)

temperature, tke and eps boundary condition -> Zero normal gradient (NE).

13. REGION AND INITIALIZATION :-

Here have the velocity in X-directions, also given temperature and pressure values.

14. TURBULANCE MODELLING :-

Here, I use Reynold average Navier stroke (RANS), turbulent model and all other values are by default.

15. BASE GRID SIZE :-

Here, I take base grid size along X, Y, Z directions.

16. ADAPTIVE MESH REFINEMENT :-

This is the adaptive mesh refinement, where maximum cells are 2,00,000, available region change to Active region. Here I consider Temperature as a AMR group where I consider following data-

Embed type -> Sub grid scale,

Timing control type -> SEQUENTIAL,

Maximum embeding level -> 2,

Sub grid criteria -> 0.05,

Start cycles -> 5000,

End cycles -> 999999999.

Plotting results are as follows :-

1. MASS FLOW RATE FOR INLET :-

Here, mass flow rate decreases with respect to number of cycle increases because of large velocity.

2. MASS FLOW RATE FOR OUTLET :-

Here mass flow rate increases with respect to number of cycles.

3. PRESSURE PLOT :-

This is the pressure plot which is increases very slowly.This is because of close boundary.

4. PRESSURE CONTOUR FOR SURFACE WITH EDGE :-

At inlet have more pressure and at outlet have less pressure and at middle it crosses and we get glowing part. Here have velocity of light is more than the air velocity and thats why both part croses and glowing the middle part. Here for subsonic velocity because the value for Mach number is less here. Here, middle cell is smaller because of sub grid size and all other are due to base grid size.

5. PRESSURE CONTOUR FOR SURFACE :-

Here, showing clear picture of Mach number. Here mach number increases linearly from inlet and trailing edge.

Minimum pressure = 1900000 Pa

Maximum pressure = 2700000 Pa.

6. TEMPERATURE PLOT :-

First temperature increases and then decreases due to low velocity of fluid. Due to shock balance temperature automatically decreases as given above.

7. TEMPERATURE CONTOUR FOR SURFACE WITH EDGE :-

At inlet, temperature increases and at trailing edge it decreases.

At inlet maximum temperature = 270k

At outlet maximum temperature = 300k

8. TEMPERATURE CONTOUR FOR SURFACE :-

This is temperature contour which showing without meshing.

9. TOTAL CELLS :-

Here have total cells which are constant upto certain limit and then increases with respect to number of cycles. This is because of shock capturing from inlet to outlet. Total cell increases as more time steps available and more shocks are present.

10. VELOCIY MAGNITUDE PLOT :-

Here velocity magnitude first decreases and then increases as shown above due to the shock capturing capacity.

11. VELOCITY CONTOUR FOR SURFACE WITH EDGE :-

At inlet velocity is lower and at trailing edge velocity is higher so that the mach number increases linearly and meshing decreases due to sub grid size.

12. VELOCITY CONTOUR FOR SURFACE :-

It's velocity contour without meshing where,

At inlet velocity = 670 m/s,

At outlet velocity = 740 m/s.

From above plotting I created animation files to proove the set-up is work or not well.

Animation results are as follows :-

1. PRESSURE ANIMATION FOR SURFACE WITH EDGE :-

2. PRESSURE ANIMATION FOR SURFACE :-

3. TEMPERATURE ANIMATION FOR SURFACE WITH EDGE :-

4. TEMPERATURE ANIMATION FOR SURFACE :-

5. VELOCITY ANIMATION FOR SURFACE WITH EDGE :-

6. VELOCITY ANIMATION FOR SURFACE :-

p

CASE-2 :-

For velocity = 100 m/s,

Sub grid scaling for temperature = 0.05k

1. Inlet Boundary set-up :-

Here I changes inlet velocity from 680 m/s to 100 m/s for seeing sub grid size and shock effects and, all others settings are same as that of previous case.

2. Region and initialization set-up :-

Here, also I changes the velocity from 680 m/s to 100 m/s in X-directions and all others are same as that of previous case.

3. Mass flow rate for inlet :-

Here mass flow rate increases at certain number of cycle and the it goes to stable everywhere. This is due to the velocity decreases.

4. Mass flow rate for outlet :-

It also decreases with respect to number of cycles, due to the velocity decreases.

5. Total cell plot :-

Total cell constant for every cycle because there have no more shock occured due the velocity decreases.

6. Pressure plot :-

Here, pressure plot increases by very slow rate.

7. Temperature plot :-

Temperature plot also increases by very slow rate with respect to number of cycles.

8. Velocity magnitude plot :-

It's also increases with respect to the number of cycles.

9. Pressure contour for surface :-

For lower velocity pressure is same on overall region i.e. from 5.8 Pascal.

10. Temperature contour for surface :-

Temperature is higher and moderate on overall region as shown above.

11. Velocity contour for surface :-

Here also velocity is same on overall region i.e. 98 m/s.

From above results I created animation files as follows-

12. Pressure Animation :-

13. Temperature animation :-

14. Velocity animation :-

FINAL RESULTS :-

1. As we decrease the SGS temperature, the mesh becomes fine. The physical meaning of sub grid temperature is an error between the actual temperature and resolved temperature at each cell. When the error is larger than the specified value, converge refines the cell automatically. So when a lower value of SGS is specified, effectively the degree of freedom for error is less and hence more cells get refined.

2. From the temperature, pressure and velocity plots it can be observed that-

Temperature reduces to around 290k at the outlet (T2 > T1).

Pressure at the outlet is lower than the inlet (P2 > P1).

3. With fine mesh i.e. low value of SGS, the expansion of fan is observed more clearly and its width is also less.

4. When an initial condition 100m/s is provided, there is no shock formation in the domain as the flow at the inlet is subsonic. From the plots, it can be obtained that the velocity decreases and pressure increases after the corner as there is an increase in flow area. This is in accordance with Burnoulli's principle.

CONCLUSIONS :-

From above ploting results we can conclude that-

1. At lower velocity mass flow rate increases at inlet and decreases at outlet and vice versa.
2. At lower velocity pressure plot increases linearly from higher value and higher velocity pressure plot also increases but from comparitively lower values.
3. At lower velocity temperature plot increases only but for higher velocity temperature plot first increases and then decreases and vice versa.
4. At lower velocity velocity plot are increases only but for higher velocity temperature plot first decreases and then increases and vice versa.
5. Total cells are same for all velocity but at higher velocity it increases suddenly after beyond limit.\
6. For higher velocity pressure, temperature and velocity contour increases mach number linearly from start to trailing edge but for lower velocity pressure and velocity contours have lower value but for temperature contour there have higher and moderate temperature.
7. From velocity value and sub grid size value we can calculate and finded mach value.
8. At inlet, entering velocity is in the form of supersonic but upto trailing edge it changes to subsonic form due to velocity decreases at outlet.

As the SGS value is lower, the mesh refinement region becomes larger. It is due to reason that converge refines the mesh if the temperature difference is higher than the SGS value. So it helps to visualize the shock region better.

In this way, I set up the Prandtl Meyer shock Problem settings and finded the contour and plot by drawing plot over line also calculate total cell and mass flow rate at inlet and outlet end. For clarify this I create animation file which are showing that the clear results i.e. how mach number increases and decreases. Also comparing the plot by changing the velocity from higher to lower.

For shock problem velocity and sub griding size are main factor to change the results, if there have propar velocity and sub griding size then we should get good and accurate results. In this way I calculate and plotted everythings here to take the propar decissions.