NPARC Alliance Validation Archive
Validation Home   >   Archive   >   Shock Tube   >   Study #1

Shock Tube: Study #1

Figure 1 is described in the surrounding text
Figure 1. The density contours (normalized by state 4) for the flow in a shock tube.

Introduction

This study is an example study demonstrating the computational process for simulating the unsteady, inviscid flow in a shock tube. The computational results are compared with analytic results as obtained from the methods decribed in the text by Anderson.

Download tar File

All of the archive files of this validation case are available in the Unix compressed tar file stube01.tar.gz. The files can then be accessed by the commands:

gunzip stube01.tar.gz
tar -xvof stube01.tar

Grid

The grid is a single-block, planar grid generated using the Fortran program stubegr.f. The grid is evenly spaced with a grid density of 100 axial grid points and 7 radial grid points. A Plot3d file (formatted, single-block, whole, two-dimensional) is named stube.x.fmt.

The Plot3d grid file is converted to the common grid file (*.cgd) format using the CFCNVT utility. This is done using the command:

cfcnvt < cfcnvt.x.com

where the file cfcnvt.x.com is a command input file (standard Fortran input unit). A common grid file named stube.cgd is created.

Initial Conditions

The initial flow conditions are stationary flow with different pressure and density values on either side of the diaphragm. The flow conditions are specified in the WindUS data input file, stube01.dat, described later. The diaphram is located at x = 0.5 ft (grid index I = 50). A 10/1 pressure ratio is specified with the static pressure and temperature values given in Table 1, below.
Table 1. Initial conditions. 
Region Pressure (psia) Temperature (R)
1 1.0 416.0
4 10.0 520.0

Boundary Conditions

The INVISCID WALL boundary condition is applied at the J1 and JMAX boundaries of the grid. The FROZEN boundary condition is applied at the I1 and IMAX boundaries of the grid.

The boundary conditions are specified using GMAN. The procedure using the graphics interface is as follows:

Start GMAN.
gman

Read in the common grid file by typing at the GMAN prompt.

file stube.cgd

Specify that units are feet / slugs / seconds.

units fss
Switch to graphics mode.
swi
Display the grid.
Pick SHOW (upper right panel)
Pick SHOW SURFACES
Pick PICK K-PLANE
At this point the grid will be displayed since there is only one grid plane, k=1. Specify the boundary conditions.
Pick BOUNDARY COND. (left panel)
GMAN defaults to zone 1 since it is the only zone. The sequence is to select the boundary and specify the boundary condition type on that boundary. The initial boundary condition type is UNDEFINED. GMAN also defaults to picking a zone and picking zone 1. It is now waiting for the boundary to be specified.
Pick I1 (lower left panel)
At this point, the I1 grid plane is displayed. Since this case is planar, a line is displayed. To select the boundary condition type, execute the following menu choices
Pick MODIFY BNDY
Pick CHANGE ALL
Pick FROZEN
At this point, the bottom panel should contain the text "7 points were changed". The boundary condition specification is now saved for this boundary by the following menu choices
Pick BOUNDARY COND. (top left panel)
Pick YES-UPDATE FILE
The common grid file stube.cgd is modified to include this boundary condition specification. The boundary condition specification can be checked by selecting
Pick IDENTIFY PNTS. (left panel)
Pick FROZEN
The individual grid points will show up in color to indicate that those points are specified as FROZEN. The boundary conditions for the remaining three boundaries are set in a similar manner. For the IMAX boundary,
Pick BOUNDARY COND. (top left panel)
Pick PICK ZONE/BNDY
Pick PICK IMAX
Pick MODIFY BNDY
Pick CHANGE ALL
Pick FROZEN
Pick BOUNDARY COND.
Pick YES-UPDATE FILE
For the J1 boundary,
Pick BOUNDARY COND. (top left panel)
Pick PICK ZONE/BNDY
Pick PICK J1
Pick MODIFY BNDY
Pick CHANGE ALL
Pick INVISCID WALL
Pick BOUNDARY COND.
Pick YES-UPDATE FILE
For the JMAX boundary,
Pick BOUNDARY COND. (top left panel)
Pick PICK ZONE/BNDY
Pick PICK JMAX
Pick MODIFY BNDY
Pick CHANGE ALL
Pick INVISCID WALL
Pick BOUNDARY COND.
Pick YES-UPDATE FILE
Exit GMAN.
Pick TOP (top left panel)
Pick END
Pick YES-TERMINATE
The common grid file stube.cgd is now ready.

As GMAN operates, a journal file gman.jou is created which can later be used as a command file to re-run GMAN. Here this file has been renamed gman.com and GMAN can be re-run as:

gman < gman.com

Creating a command file from scratch represents an alternative to the graphical approach within GMAN.

Computation Strategy

The computation is performed using the time-marching capabilities of WIND to march in a time-accurate manner to simulate the unsteady flow. A constant time step is specified to be used for each time step. The number of time steps is specified such that the desired time interval is spanned.

Input Parameters and Files

The input data file for WIND is stube01.dat. The keyword inviscid indicates that the inviscid flow is assumed. The axisymmetric keyword indicates that the flow domain is assumed to be axisymmetric with a waterline of 0.0 feet and a degree of rotation of 5 degrees. The arbitrary inflow keyword block indicates that the static conditions described in Table 1, above, are to be used, with the higher pressure in the upstream half of the shock tube. The stages keyword indicates that the four-stage Jameson-type of Runge-Kutta time integration algorithm is to be used. The implicit none keyword indicates that the time integration is to be explicit. The cycles and iterations keywords indicate that five cycles with 100 time steps in each cycle is to be run. The timesteps keyword indicates that the actual time step in seconds is specified. The history keyword indicates that the static pressure at grid point (50,4,1) is to be printed out to the time history data file (stube.cth) every 5 time steps.

The input data file stube01.dat works only for Wind Version 5, and Wind-US versions 1 and 2. The input data file for version 4.0 of WIND is stube01.v4.dat. The only difference in the files is the form history keyword input, which changed between versions 4.0 and 5.0.

The input data file stube.spawn.dat demonstrates the use of the spawn keyword to write out the flow solution at the end of each cycle. This is useful in visualizing the unsteady flow. More information on using this option can be found in the documentation of the spawn keyword.

Computation

Version 2 of the Wind-US flow solver was run. For this computation, only one processor is used. The list file is stube01.lis and contains the output from the computation and lists the residual information for all of the iterations. The flow data file is stube01.cfl and contains the final solution.

Convergence

Information on the convergence history of the L2 residual can be obtained from the list file stube.lis using the utility RESPLT,

resplt < resplt.nsl2.com

The file resplt.nsl2.com is a command file containing the inputs for RESPLT to output the file nsl2.gen containing the residual history data.

The program CFPOST can now be used to read in nsl2.gen and plot the convergence history.

cfpost < cfpost.nsl2.com
Fig. 2 shows a plot of the residual history. The residual reaches an asymptotic value; however, notice that the drop in residual is about half of an order of magnitude. This is because the problem is unsteady and the solution is suppose to change over an iteration, which is now truly a time step.

Figure 2 is described in the surrounding text
Figure 2. Plot of the L2 residual history.

Post-Processing

The CFPOST utility was used to compute and plot the axial distributions of the pressure, density, and Mach number along the shock tube.

The following CFPOST files compute the pressure, density and Mach number, respectively: cfpost.p.com, cfpost.rho.com, and cfpost.mach.com. They also generate corresponding GENPLOT files of the results named p.gen, rho.gen, and mach.gen. So, for example, to create the pressure plot in Figure 3, the following procedure is used. First the pressure is calculated and the p.gen GENPLOT file is written by executing CFPOST:
cfpost < cfpost.p.com

Then another CFPOST command file, cfpost.p_comp.com is used to plot these Wind-US results with the theoretical results contained in the GENPLOT file p.theory.gen.

cfpost < cfpost.p_comp.com
The results are shown in Figure 3, below. Similarly, the density along the shock tube is calculated as follows.
cfpost < cfpost.rho.com

Then the CFPOST command file, cfpost.rho_comp.com is used to plot these Wind-US results for density with the theoretical results contained in the GENPLOT file rho.theory.gen.

cfpost < cfpost.rho_comp.com

These results are shown in Figure 4, below. The Mach number in the shock tube is calculated using the same procedure.

cfpost < cfpost.mach.com

Then the CFPOST command file, cfpost.mach_comp.com is used to plot these Wind-US results for mach number with the theoretical results contained in the GENPLOT file mach.theory.gen.

cfpost < cfpost.mach_comp.com

These results are shown in Figure 5, below.

Comparisons of the Results

The distributions of pressure, density, and Mach number along the shock tube at the final time are compared to the theoretical results and plotted in Figs. 3, 4, and 5. The comparisons show very good aggreement with minimal overshoot.

Figure 3 is described in the surrounding text
Figure 3. Comparison of the pressure at the final time.

Figure 4 is described in the surrounding text
Figure 4. Comparison of the density at the final time.

Figure 5 is described in the surrounding text
Figure 5. Comparison of the Mach number at the final time.

Log and Contact Information

This study was created on May 15, 1998 and updated on October 28, 2002 by John W. Slater. It was updated again on March 2, 2009 by Julianne Dudek. You may contact us at:
John Slater
NASA John H. Glenn Research Center, MS 5-12
21000 Brookpark Road
Brook Park, Ohio 44135
Phone: (216) 433-8513
e-mail: John.W.Slater@nasa.gov
Julianne Dudek
NASA John H. Glenn Research Center, MS 5-12
21000 Brookpark Road
Brook Park, Ohio 44135
Phone: (216) 433-2188
e-mail: Julianne.C.Dudek@nasa.gov

Last Updated: Tuesday, 2-Mar-2009 11:49:14 EDT