Shock tube simulation project

Aim: Shock Tube Simulation

Objectives:

• To set up a transient shock tube simulation
• To plot the pressure and temperature history in the entire domain
• To plot the cell count as a function of time
• To show the mesh generation using PARAVIEW.

Introduction:

The shock tube is an instrument used to replicate and direct blast waves at a sensor or a model to simulate actual explosions and their effects, usually on a smaller scale. Shock tubes can also be used to study aerodynamic flow under a wide range of temperatures and pressures that are difficult to obtain in other types of testing facilities. Shock tubes are also used to investigate compressible flow phenomena and gas-phase combustion reactions.

A shock wave inside a shock tube may be generated by a small explosion (blast-driven) or by the buildup of high pressures which cause the diaphragm to burst and a shock wave to propagate down the shock tube.

Theory:

A simple shock tube is a tube, rectangular or circular in cross-section, usually constructed of metal, in which gas at low pressure and gas at high pressure is separated using some form of a diaphragm. The diaphragm suddenly bursts open under predetermined conditions to produce a wave propagating through the low-pressure section. The shock that eventually forms increases the temperature and pressure of the test gas and induces a flow in the direction of the shock wave. Observations can be made in the flow behind the incident front or take advantage of the longer testing times and vastly enhanced pressures and temperatures behind the reflected wave.

The low-pressure gas, referred to as the driven gas, is subjected to the shock wave. The high-pressure gas is known as the driver gas. The corresponding sections of the tube are likewise called the driver and driven sections. The driver gas is usually chosen to have a low molecular weight, for safety reasons, with high speed of sound, but maybe slightly diluted to 'tailor' interface conditions across the shock. To obtain the strongest shocks the pressure of the driven gas is well below atmospheric pressure (a partial vacuum is induced in the driven section before detonation).

The test begins with the bursting of the diaphragm. Several methods are commonly used to burst the diaphragm.

• A mechanically-driven plunger is sometimes used to pierce it or an explosive charge may be used to burst it.
• Another method is to use diaphragms of plastic or metals to define specific bursting pressures. Plastics are used for the lowest burst pressures, aluminium and copper for somewhat higher levels and mild steel and stainless steel for the highest burst pressures. These diaphragms are frequently scored in a cross-shaped pattern to a calibrated depth to ensure that they rupture evenly, contouring the petals so that the full section of the tube remains open during the test time.
• Yet another method of rupturing the diaphragm utilizes a mixture of combustible gases, with an initiator designed to produce a detonation within it, producing a sudden and sharp increase in what may or may not be a pressurized driver. This blast wave increases the temperature and pressure of the driven gas and induces a flow in the direction of the shock wave but at a lower velocity than the lead wave.

The bursting diaphragm produces a series of pressure waves, each increasing the speed of sound behind them so that they compress into a shock propagating through the driven gas. This shock wave increases the temperature and pressure of the driven gas and induces a flow in the direction of the shock wave but at a lower velocity than the lead wave. Simultaneously, a rarefaction wave, often referred to as the Prandtl-Meyer wave, travels back into the driver gas.

Geometry Setup:

1. Import Geometry:

File> Import> Import STL> Browse to the location of geometry and import the .STL file

2. Diagnosis:

There are several intersections because of the diaphragm, which is created between the two parts in CAD modelling which is not required for the simulation.

The nonmanifold problems occurred due to the diaphragm as the vertices are shared by more than two triangles.

3. Delete Diaphragm:

The diaphragm can be viewed by hiding the triangles.

`triangle` as an entity> By `20^0` angle> Using cursor pick> Selecting triangles to be hidden> Right Click> Hide selection.

In the below images it can be seen that diaphragm vertices are shared by more than two triangles.

To avoid the nonmanifold problems, the diaphragm must be deleted.

Geometry> Repair> Delete> `triangle` as an entity> By `20^0` angle> Using cursor pick> Selecting triangles representing diaphragm> Apply.

4. Transform Normal:

Again Diagnosis, it shows normal orientation issues and it is because of the wrong direction of the normals. Normals must point inside the domain for simulation.

Normal orientation can be corrected by transforming normal.

Geometry> Transform> Normal> Propagate change from single triangle> `triangle` as an entity> Any> Selecting one triangle> Apply

Now diagnosis doesn't show any problem with the geometry.

5. Assigning Boundaries:

In converge every part of geometry is assumed to be a triangle. Selecting triangles and grouped to assign to specific boundaries.

Boundary> Flag> +> Create Boundaries> OK

a) High Pressure:

`triangle` as an entity selection-By Angle `20^0`> Using cursor pick option> Selecting the triangles as shown in below image> Assigning them to High-Pressure Boundary> Apply

a) Low Pressure:

`triangle` as an entity selection-By Angle `20^0`> Using cursor pick option> Selecting the triangles as shown in below image> Assigning them to Low-Pressure Boundary> Apply

Case Setup:

1. Application Type:

Selecting an application type as "Time Based".

2. Materials:

Gases are inside the shock tube.

Checking the solid simulation & Species option and unchecking the reaction mechanism option.

Confirming the permission to overwrite the gas.dat, mech.dat and therm.dat by loading predefined mixture Air.

a) Gas Simulation:

Keeping the default values for gas simulation > OK

b) Global Transport parameters:

Keeping all the default values for global transport properties > OK

c) Species:

Gas - Adding `N_2` and `O_2`

3. Simulation Parameters:

a) Run Parameters:

Solver: Transient

Simulation mode is Full hydrodynamic and Gas flow solver is Compressible.

Under Misc. tab disabling use shared memory option > OK.

b) Simulation Time Parameters:

Set the start time as 0 s, end time as 0.003 s, Initial and minimum time step as 1e-9 and Maximum time step as 1 second.

Keep all the other CFL values default > OK.

c) Solver parameters:

Keeping all settings as default > OK.

4. Regions and Initialization:

Grouping the boundaries to form volumetric zones called region. This defines a volumetric region with initial conditions like initial pressure, initial temperature, velocity and species concentration.

Enabling Events option.

a) Regions and Initialization:

Selecting Add to create a new region,

Add 2 Regions, namely the high-pressure region and low-pressure region.

i) High-Pressure Region:

For the High-Pressure Region setting the Stream ID as 0

Velocity - [0, 0, 0] m/s

Temperature as 300 K

Pressure as 600000 Pa.

Species `N_2` with Mass Fraction value 1.

ii) Low-Pressure Region:

For the Low-Pressure Region setting the Stream ID as 0

Velocity - [0, 0, 0] m/s

Temperature as 300 K

Pressure as 101325 Pa.

Species `O_2` with Mass Fraction value 1

b) Events:

When 2 regions are in contact with each other with no boundary in-between them, it is necessary to use events. The events flag creates triangles in the place between the 2 regions called as disconnect triangles.

Events can be cyclic, permanent or sequential. Based on how the event should occur.

In this case, the event can repeat itself in a sequence. Based on this setting, Converge creates disconnected triangles. These triangles are temporarily created to simulate the event. Triangles are created for when there should be no interaction between the 2 regions and these triangles are deactivated when there is supposed to be an interaction between the 2 regions.

At the start, 2 regions are not in contact and for that event is set as given below:

• Start second: 0 sec
• Region A: High pressure
• Region B: Low pressure
• Event: Close

After 0.001 seconds, to deactivate disconnected triangles and provide interaction between 2 regions the event is set as:

• Start second: 0.001 sec
• Region A: High pressure
• Region B: Low pressure
• Event: Open

5. Boundary:

i) High Pressure:

Boundary Type- WALL - Slip

Velocity Boundary Condition:

Wall motion Type - Stationary

Surface Movement: FIXED

Region: High-Pressure Region

ii) Low Pressure:

Boundary Type- WALL - Slip

Velocity Boundary Condition:

Wall motion Type - Stationary

Surface Movement: FIXED

Region: Low-Pressure Region

After assigning the boundary, enabling group colour by regions.

6. Physical Models:

• Turbulence Modeling:

Keep the Turbulence model as RNG k-epsilon

Keep all the settings default > OK

7. Grid Control:

Enabling Adaptive Mesh Refinement.

i) Base Grid:

Base Grid size = 0.001 m

ii) Adaptive Mesh Refinement:

Adaptive Mesh Refinement is a method of adapting the accuracy of a solution within certain sensitive or turbulent regions of simulation, dynamically and during the time the solution is being calculated. When solutions are calculated numerically, they are often limited to pre-determined quantified grids as in the Cartesian plane which constitutes the computational grid, or 'mesh'.

The available regions are High-Pressure Region and Low-Pressure Region which are changed to active regions.

• Minimum cells - 1
• Maximum cells - 200000
• Embedded type - Sub-Grid-Scale (SGS)
• Species scheme is used for AMR with Embed type -SGS (Sub Grid-Scale)
• Species Name - `N_2`
• Sub-grid criterion - 0.001
• Max embedding level as 3
• Timing-control type - Permanent

7. Output/Post-Processing

a) Post-Variable Selection

Using default settings for Typical Variables/Geometry.

Under Species/Passive checking the Mass Fraction box and adding `O_2` and `N_2`.

b) Output Files

Setting the time interval for writing 3D output data files - 1e-05 s

The time interval for writing restarting output to 1e-06 s

8. Export all the input files

Export all the input files into a separate folder

Files>Export>Export Input Files>setting up the desired location>OK

9. Copy and paste the mpiexec.exe file to the input files folder

Run Simulation

Cygwin - Cygwin is a collection of GNU and open-source tools that provide functionality similar to Linux distribution on windows.

To run the simulation, open Cygwin and navigate to the folder in which the input files are exported.

To run the simulation command with executable to be entered

$ mpiexec.exe -n 4 converge-intelmpi.exe restricted </dev/null> logfile.txt &

Results:

1. Meshing:

a) Base Mesh:

Base grid with grid size=1mm in dx,dy and dz.

b) Mesh with AMR:

Converge keep monitoring the curvature of the fluid property provided in the AMR and the mesh refinement takes place when the property curvature variation between the consecutive grid is more than the defined SGS value. This is done to capture the flow properties in the shock flow region accurately.

Due to AMR, the mesh is refined at the location where there more change in the slope of N2 mass fraction.

Animation:

From the meshing animation, we can see that the mesh refinement takes place at N2 - O2 interface and this mesh refinement keeps changing as the shock propagates back and forth.

2. Region Representation:

3. Mass Fraction `N_2` Contour:

As the diaphragm breaks the high-pressure region `N_2` pushes the `O_2` low-pressure region and then the low-pressure region reflects and pushes it towards `N_2` region. This shock waves keep moving back and forth as the waves have different pressures and due to which this occurs and after a while, if the simulation is run for more time the pressure eventually settles and converges to a value.

Animation:

4. Pressure Contour:

As the diaphragm is ruptured shock wave is generated & it enters driven region reaching the downstream end of the driven section. The shock wave produced is an oblique shock wave where the static properties increases at shock occurrence (after 0.001s) i.e. pressure, density and temperature in the `O_2` region (low-pressure region) increase. This also results in the expansion of `N_2` and hence the pressure decreases. The moment it hits the downstream end and gets reflected, their pressure values changes respectively and this flow pattern occurs till the shock intensity is decreased.

Animation:

Pressure Plot:

The High-pressure region is at 0.6 MPa and from the plot, it is clearly visible that simulation never reaches that pressure. As soon as the diaphragm breaks at 0.001 seconds the pressure increases and the high-pressure region pushes the low-pressure region again low-pressure region pushes high-pressure region and it is just like pushing and dissolving until reflected wave can occur. Generation of a shock wave in the low-pressure region leads to pressure rise and a slight drop in the plot is visible.

Both the regions generate a sine wave-like structure with their shock waves. This pressure eventually dissolves into steady value when a simulation is run for more time steps.

5. Temperature Contour:

Temperature contour also follows the same pattern as that of pressure i.e. temperature increases in `O_2` region when the shock wave is generated after the diaphragm is ruptured. When shock wave gets reflected then the temperature in `O_2` region decreases and vice-versa is the case for `N_2` region and same decreasing pattern occurs till shock intensity is decreased.

Animation:

Temperature Plot:

Initially, the temperature throughout the shock tube is 300 K. After the shock waves are generated, the temperature in both the regions start varying cyclically. The temperature in the low-pressure region varies between 300 K and 340 K. The temperature in the high-pressure region varies between 300 K and 200 K before attaining steady-state.

Velocity Contour:

The velocity increases sharp as diaphragm burst and increases forth and decreases back and vice versa, such waves create shock patterns like diamond shock pattern observed near the low-pressure region.

Animation:

Cell Count Plot:

The total cell count plot shows the variation of mesh sizes as a function of time during the entire simulation of 0.003 seconds. In the events option setting has done for the diaphragm to burst at 0.001 seconds as a result of that diaphragm bursting sudden spike is observed after 0.001 seconds i.e interaction of two regions. When the simulation begins the high-pressure region pushes the low-pressure region and the AMR gets activated, the variations of which can be observed in the plot. The sudden decrease in mesh size is when the reflected pressure pushes it back.

Conclusion:

The shock tube is mainly used to study the properties of the fluid at an elevated temperature and pressure.

The burst of the diaphragm creates compression waves which propagate in the driven section and expansion waves traversing in the driver section.

The shock wave when hits the wall gets refracted back and the mixing of the fluid occurs. The pressure gradients between driver and driven section are due to viscous drag and friction between shock wavefront and shock tube walls.

The gas behind the reflected wave is at rest relative to the shock tube. Usually, there is some mixing and interdiffusion of gases at the contact surface so that the temperature fall is less sudden than in the ideal case. In the case of IC engines, this shock tube data is important to find the ignition delay i.e, how long it takes to ignite fuel and air.