<jupyter_start><jupyter_text>Shock tube and interface interactionThis test case is a classical application of a shock tube facility to study shock/interface interaction. Most interactions between locally uniform states are theoretically available to computation: waves speeds and intermediate states are computed. The $x/t$ diagram is then provided at the end. The limitation of these computation is given by expansion/discontinuities interaction where non-uniform flows cannot be simply described. Shock generationThe first Riemann problem is given the initial interaction between 2 zones with a pressure ratio. The driver gas `q20` is air. The driven gas `q10` (which will propagate the shock wave) is air too. Gas states will be denoted `qxn` where `x` is the initial zone number and 'n' the successive states in time.<jupyter_code>Pinit  = 1.e5   # reference pressure
Pratio = 2.38   # pressure ratio (left over right)
Tisoth = 293.15 # same temperature for both gas
r_air = 287.1
r_He  = 2078.
gam_air = 1.4
gam_He  = 5./3.
#
# driver gas (q20)
q20 = uq.unsteady_state(rho = Pinit*Pratio/(r_air*Tisoth),
                        u   = 0., 
                        p   = Pinit*Pratio, gamma=gam_air)
# driven gas (q10)
q10 = uq.unsteady_state(rho = Pinit/(r_air*Tisoth),
                        u   = 0., 
                        p   = Pinit, gamma=gam_air)
# Riemann problem
pb1 = riem.riemann_pb(q20, q10)
q21 = pb1.qstarL()
q11 = pb1.qstarR()
print("shock compression ratio: {:.3f}".format(q11.p/q10.p))
print("shock Mach number      : {:.3f}".format(pb1.right_fastest()/q10.asound()))<jupyter_output><empty_output><jupyter_text>Shock / Interface interactionThe previous computed shock wave will reach an interface between `q10` (driven gas) and `q10`. The Riemann problem when the shock interacts is defined by `q11` (shock downstream flow) and `q10`.<jupyter_code># Helium 
q00 = uq.unsteady_state(rho = Pinit/(r_He*Tisoth),
                        u   = 0., 
                        p   = Pinit, gamma=gam_He)
# Riemann problem
pb2 = riem.riemann_pb(q11, q00)
q12 = pb2.qstarL()
q01 = pb2.qstarR()
print("transmitted shock compression ratio: {:.3f}".format(q01.p/q00.p))
print("interface speed                    : {:.3f}".format(q12.u))<jupyter_output><empty_output><jupyter_text>Shock reflexionThe transmitted shock wave will bounce or reflect on the wall. Even if the reflected shock wave could have been computed directly, the trick here is to defined a symmetric state and a new Riemann problem. Helium is stopped and compressed to `q02` state.<jupyter_code># define symmetric (ghost state) 
q01gh = q01.copysymmetric()
# Wall interaction
pbw = riem.riemann_pb(q01, q01gh)
q02 = pbw.qstarL()
print("shock compression ratio: {:.3f}".format(q02.p/q01.p))
print("shock Mach number      : {:.3f}".format(pbw.left_fastest()/q01.asound()))<jupyter_output><empty_output><jupyter_text>Reshock of interfaceThe reflected shock will reach the interface with a new complex Riemann problem: the transmitted shock wave will change the air from `q12` to `q13`. A reflected right wave will be produced: it can be an expansion or a shock depending on the initial condition and correlated to the new interface velocity.<jupyter_code># Riemann problem
pb3 = riem.riemann_pb(q12, q02)
q13 = pb3.qstarL()
q03 = pb3.qstarR()
#print("interface speed        : {:.3f}".format(q13.u))<jupyter_output><empty_output><jupyter_text>Position of interactionsGiven all waves velocities, all interaction positions can be computed and drawn in $x/t$ diagram.<jupyter_code># positions
x_int   = 2.5 
x_Rwall = 3.
x_Lwall = -1.
#
def xt_intersect(x1, t1, u1, x2, t2, u2):
    ti = (-x2+x1 + t2*u2-t1*u1)/(u2-u1)
    return x1+u1*(ti-t1), ti
#
# shock / interface 
t_shock_int  = x_int / pb1.right_fastest()
#print("shock / interface   interaction at (x,t)=(%.3f,%.5f)"%(x_int, t_shock_int))
# shock / wall
t_shock_wall = (x_Rwall-x_int) / pb2.right_fastest() + t_shock_int
#print("shock / wall        interaction at (x,t)=(%.3f,%.5f)"%(x_Rwall, t_shock_wall))
# interface / reshock
x_reshock, t_reshock = xt_intersect(x_int,   t_shock_int,  q12.u,
                                    x_Rwall, t_shock_wall, pbw.left_fastest())
#print("interface / reshock interaction at (x,t)=(%.3f,%.5f)"%(x_reshock, t_reshock))
t_max = 1.2*t_reshock<jupyter_output><empty_output><jupyter_text>x/t diagram<jupyter_code>import matplotlib.pyplot as plt
%matplotlib inline
#
plt.figure(figsize=(15,8))
plt.axis([-1, 3, 0, t_max])
plt.xlabel(r'$x$', fontsize=18)
plt.ylabel(r"$t$", fontsize=18)
plt.title(u"shocktube in x/t diagram", fontsize=18)
#
# plot shocks
plt.plot([0, x_int,       x_Rwall,      x_reshock, x_reshock+pb3.left_fastest()*(t_max-t_reshock)], 
         [0, t_shock_int, t_shock_wall, t_reshock, t_max],
         color="red", linewidth=2)
# plot interface
plt.plot([x_int, x_int,       x_reshock, x_reshock+q13.u*(t_max-t_reshock)],
         [0,     t_shock_int, t_reshock, t_max],
         color="green", linewidth=2, linestyle="dashed")
plt.plot([0, q11.u*t_max],
         [0, t_max],
         color="green", linewidth=2, linestyle="dashed")
# plot expansion
t_exp_wall = x_Lwall/q20.left_acoustic()
plt.plot([0, x_Lwall,    x_Lwall+q21.right_acoustic()*(t_max-t_exp_wall)],
         [0, t_exp_wall, t_max],
         color="orange", linewidth=1)
plt.plot([0, x_Lwall],
         [0, x_Lwall/q21.left_acoustic()],
         color="orange", linewidth=1)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.12424409015942826, 'fraction_numerical': 0.035733919736118745, 'mean_word_length': 3.908273381294964, 'pattern_counts': {'":': 0, '<': 25, '<?xml version=': 0, '>': 25, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '6', 'dir': 'jupyter-structured-clean-dedup', 'id': '580304', 'license': 'MIT', 'n_tokens_mistral': 2104, 'n_tokens_neox': 1935, 'n_words': 522, 'path': 'docs/examples/instance.riemann.ipynb', 'repo_name': 'PierreMignerot/aerokit', 'size': '5457'}