Analytical Prediction of Low-Frequency Fluctuations Inside a One-Dimensional Shock
NNoname manuscript No. (will be inserted by the editor)
Analytical Prediction of Low-Frequency FluctuationsInside a One-dimensional Shock
Saurabh S. Sawant · Deborah A. Levin · Vassilios Theofilis
Received: date / Accepted: date
Abstract
Linear instability of high-speed boundary layers is routinely ex-amined assuming quiescent edge conditions, without reference to the internalstructure of shocks or to instabilities potentially generated in them. Our recentwork has shown that the kinetically modeled internal nonequilibrium zone ofstraight shocks away from solid boundaries exhibits low-frequency molecularfluctuations. The presence of the dominant low frequencies observed using theDirect Simulation Monte Carlo (DSMC) method has been explained as a con-sequence of the well-known bimodal probability density function (PDF) of theenergy of particles inside a shock. Here, PDFs of particle energies are derivedin the upstream and downstream equilibrium regions, as well as inside shocks,and it is shown for the first time that they have the form of the non-centralchi-squared (NCCS) distributions. A linear correlation is proposed to relatethe change in the shape of the analytical PDFs as a function of Mach number,within the range 3 ≤ M ≤
10, with the DSMC-derived average characteristiclow-frequency of shocks, as computed in our earlier work. At a given Machnumber M = 7 .
2, varying the input translational temperature in the range89 ≤ T tr, / ( K ) ≤ S. S. SawantDepartment of Aerospace, the University of Illinois at Urbana-Champaign, IL, 61801, USA.E-mail: [email protected]. A. LevinDepartment of Aerospace, the University of Illinois at Urbana-Champaign, IL, 61801, USA.V. TheofilisSchool of Engineering, University of Liverpool, The Quadrangle, Brownlow Hill, L69 3GH,UK Escola Politecnica, Universidade S˜ao Paulo, Av. Prof. Mello Moraes 2231, CEP 5508-900, S˜ao Paulo-SP, Brasil a r X i v : . [ phy s i c s . f l u - dyn ] J a n Saurabh S. Sawant et al. derived to predict low-frequency oscillations in shocks. Our model can be usedto provide realistic physics-based boundary conditions in receptivity and linearstability analysis studies of laminar-turbulent transition in high-speed flows.
Keywords shocks · noncentral chi-squared distributions · low-frequencyfluctuations · DSMC · transition Supersonic and hypersonic laminar to turbulent flow transition strongly de-pends on the freestream noise, which includes atmospheric or tunnel-inducedturbulence, acoustic, vorticity and entropy fluctuations, and particulates [1].These freestream disturbances interact with a shock and excite the downstreamboundary layer through a process known as receptivity [2]. The importanceof understanding this process is evident from several experimental works thathave shown a strong correlation between the Reynolds number of transitionand freestream noise [3,4,5,6]. Numerical research on the process of receptiv-ity and natural transition has focused primarily on the development of modelsto account for freestream disturbances [7,8] and their interaction with shocksthrough shock-capturing [9] or shock-fitting [10] schemes. However, these worksdo not account for the internal dynamic structure of a shock-wave which, asour recent work [11] has shown, exhibits low-frequency molecular fluctuations.The important role of molecular fluctuations in the laminar to turbulent tran-sition is supported by the works of Fedorov and Tumin [12] and Luchini [13] onthe receptivity of boundary layers to equilibrium molecular fluctuations usingLandau-Lifshitz’s theory of fluctuating hydrodynamics [14]. Since freestreamnoise inevitably passes through a shock before interacting with a boundary-layer, low-frequency fluctuations in the internal translational nonequilibriumzone of a shock may not be ignored in the study of receptivity. This is trueeven when the level of freestream noise is minimal, i.e. only containing equi-librium molecular fluctuations, because the leading-edge shock interacts witha boundary-layer formation.Some examples of the notable computational fluid dynamic (CFD) litera-ture on modeling the process of natural transition and receptivity include thework of Balakumar and Chou [7], who introduced the experimentally measuredfreestream fluctuations by Marineau et al. [15] in a Mach 10 flow over a coneat the outer computational boundary of their simulation domain upstreamof the detached shock modeled using a shock-capturing scheme. Hader andFasel [8] simulated a Mach 6 flow over a flared cone geometry following theexperiments in the BAM6QT quiet tunnel at Purdue [16]. They constructed asimple free-stream noise model and introduced random pressure pertubationsin the boundary-layer downstream of the leading-edge shock to account forthe effect of freestream disturbances. Ma and Zhong [17,18,19] studied thereceptivity of supersonic boundary layer to four types of freestream distur-bances, fast and slow acoustic waves, vorticity waves, and entropy waves. The itle Suppressed Due to Excessive Length 3 unsteady interaction of the freestream disturbance and the shock was capturedusing a high-order shock-fitting scheme of Zhong [10].In our work on the low-frequency fluctuations in a one-dimensional (1-D)shock [11], the shock layer was simulated using the particle-based, high-fidelityDSMC method. The frequencies of fluctuations in shocks were found to be anorder of magnitude lower than those in the freestream for the examined Machnumber range, 2 ≤ M ≤
10. We showed that this difference results from thewell-known bimodal nature of the PDF of gas particles in the shock, as op-posed to their Maxwellian distribution in the freestream. Based on the shapeof the DSMC-derived energy distributions at M=7.2, a reduced-order two-bindynamic model was constructed to account for a large number of collision in-teractions of particles and fluxes from neighboring zones. The model correctlypredicted the order-of-magnitude differences in frequencies in the shock versusthe freestream. A Strouhal number was also defined based on the bulk ve-locity upstream of the shock and the shock-thickness based on the maximumdensity-gradient inside the shock to nondimensionalize the low-frequencies ob-tained from DSMC. It remained practically constant when the Mach numberwas varied by keeping the upstream temperature constant, however, anotherset of test cases at a constant Mach number of 7.2 and varying input temper-ature revealed that the Strouhal number decreased with decrease in upstreamtemperature.CFD simulations aimed at understanding receptivity do not model thedetails of PDFs of particles or the internal structure of a shock wave andtherefore, cannot properly account for low-frequency fluctuations. Yet it ispossible to incorporate the effect of characteristic fluctuations of a leading-edge shock by constructing simple models similar to that of Hader and Fasel [8]to understand the process of natural transition in quiet tunnels. To constructsuch models, a simple method is needed to estimate the low-frequency of ashock-wave generated at arbitrary input conditions, which may not have beensimulated explicity using DSMC. Also, the Strouhal number may not be thebest way to estimate the low-frequency, especially if the input temperature isnot the same as that considered in our earlier DSMC-simulations. In this work,we present a different approach to estimate the characteristic low-frequencyof a shock by correlating the shape of the analytically derived PDFs withthe DSMC-computed frequencies for a parameter space of 3 ≤ M ≤
10 and89 ≤ T tr, / ( K ) ≤ Saurabh S. Sawant et al. distributions in the upstream and downstream regions with weighting factorsobtained from the theory of Mott-Smith [23] (see Sect. 4) and provide a pythoncode to easily generate these PDFs (see Appendix A). The variation of the ana-lytically derived PDFs is studied as a function of Mach number and their shapecharacteristic is linearly correlated with the DSMC-derived low-frequency ofthe fluctuations inside a shock for the Mach number range, 3 ≤ M ≤
10 at T tr, = 710 K (see Sect. 5). At M =7.2 and 89 ≤ T tr, / ( K ) ≤ Fig. 1
Normalized σ xx in a 1-D, M =7.2 shock of argon. The internal structure of a 1-D, steady, Mach 7.2 argon shock simulatedusing the kinetic DSMC method is shown in figure 1, where X is the directionnormal to the shock and λ is the upstream mean-free-path, and subscripts 1and 2 are used to denote freestream and downstream values. The upstreamnumber density, n , and translational temperature, T tr, , are 10 m − and710 K, respectively, and the downstream subsonic conditions are imposed bythe Rankine-Hugoniot jump conditions. Figure 1 shows the X -directional nor-mal stress, σ xx , normalized by the upstream pressure p , where σ xx is non-zeroinside the shock, − < X/λ <
6, and zero in the upstream and downstream itle Suppressed Due to Excessive Length 5 equilibrium regions, i.e.
X/λ ≤ −
16 and
X/λ ≥
6, respectively. Numericalprobe locations F , P , and D also marked on the profile of σ xx in the freestream( X/λ = −
29) region, inside the shock where the gradient of velocity is maxi-mum (
X/λ = − X/λ = 29), respectively.The complete details of the simulation and profiles of other macroparametersthrough the shock are given in our recent work [11].To compare molecular fluctuations in the shocks versus that in the freestream,we defined the normalized X -directional molecular energy of particles, ξ x =( v x − u x ) β , where v x and u x are instantaneous and bulk velocity componentsof molecules. σ xx is related to ξ x through the first moment of the PDF of f ξ x defined by mean µ ξ x as, µ ξ x = (cid:90) ξ x f ξ x dξ x = σ xx β ρ (1)where f ξ x is the PDF of molecules in the energy space ξ x , ρ = nm is the localmass density, m = 6 . × − kg is the molecular mass, β = (2 RT tr ) − / is the the inverse most-probable-speed, T tr is the local average translationaltemperature, and R is the specific gas constant. The analytical form of f ξ x inequilibrium and nonequilibrium regions is described in Sects. 3 and 4, respec-tively. Reference [11] shows that based on equation 1, the fluctuations in σ xx are related to fluctuations in the number of particles in the energy space, ξ x . ( a ) ( b ) Fig. 2
Fluctuations in the number of particles in PDF of f ξ x obtained from DSMC as afunction of ξ x at four time-windows ( a ) at probe P inside the shock and ( b ) probe F in thefreestream. The differences in fluctuations in the shock and the freestream are shown inFigs. 2( a ) and 2( b ), respectively. The fluctuations are defined by the quantity (cid:2) N ( ξ x , w ) − (cid:104) N ( ξ x ) (cid:105) w (cid:3) , where N ( ξ x , w ) is the total number of DSMC particleswithin the energy space ξ x and ξ x + ∆ξ x , where ∆ξ x = 0 . w to w + 1 and (cid:104) N ( ξ x ) (cid:105) w denotes the number of particles in the Saurabh S. Sawant et al. same energy space averaged over all time-windows, W = 646. Each time-window is equal to 100 timesteps of 3 ns and w =0 is defined to be the end oftransient time of 6 µ s required to establish the shock structure from the initialRankine-Hugoniot jump conditions. Figure 2( a ) shows at probe P inside theshock, the fluctuations within the energy space − ± P , these valueswere found to be 37.5 and 21.4 kHz, respectively. The entire spectrum offrequencies obtained using the power spectral density (PSD) of the mean-subtracted, DSMC-derived data of σ xx are described in Ref. [11]. At probe F in the freestream, however, figure 2( b ) shows no presence of low-frequenciesand the change in the number of particles in freestream is much lower thanthat in the shock indicating smaller amplitude fluctuations. In this section the PDF of f ξ x is derived for a gas in local equilibrium usingthe local Maxwell-Boltzmann (MB) PDF [24] of particle velocities and it isshown that the energy PDF is of the class of an NCCS PDF. For completeness,PDFs of f ξ y , f ξ z , and f ξ of the Y -directional ( ξ y ), Z -directional ( ξ z ), and totalenergies of particles ( ξ ) are also derived. The first moment of f ξ y , f ξ z and f ξ are related to Y and Z -directional normal stresses, σ yy , σ zz , and pressure, p ,as µ ξ y = (cid:90) ξ y f ξ y dξ y = σ yy β ρµ ξ z = (cid:90) ξ z f ξ z dξ z = σ zz β ρµ ξ = (cid:90) ξf ξ dξ = 3 pβ ρ = 32 (2)Note that ξ y = v y β , ξ z = v z β , and ξ = ( v − u ) β , where v y , v z are theinstantaneous components in the Y and Z -directions, respectively, u y , u z arethe bulk velocity components in the Y and Z -directions, respectively, and v , u are the instantaneous and bulk speed of molecules, respectively.The local MB PDF for particles having instantaneous molecular and bulkvelocity vectors v and u , respectively, and an equilibrium translational tem-perature, T tr , is written as, f v = β π / exp (cid:16) − β ( v − u ) (cid:17) (3) f v may change spatially depending on the local values of u , T tr , and applies toparticles in equilibrium, i.e. it is valid for X/λ ≤ −
16 in the freestream and itle Suppressed Due to Excessive Length 7
X/λ ≥ u i = (cid:90) v i f v dv x dv y dv z , for i ∈ { x, y, z } T tr = m κ b (cid:90) v f v dv x dv y dv z (4)The PDF for directional velocity components can be obtained by integratingthe differential probability f v dv x dv y dv z over the remaining two componentsand written as, f v i = β √ π exp (cid:16) − β ( v i − u i ) (cid:17) ; for i ∈ { x, y, z } (5)Equation 5 can be expressed as a PDF of the Gaussian distribution havingmean (bulk velocity) u i , and variance (2 β ) − as, N v i ( u i , (2 β ) − ) = 1 (cid:112) π (2 β ) − exp (cid:32) − ( v i − u i ) β ) − (cid:33) This can be further scaled as a PDF of the Gaussian distribution of a scaledvariable ˜ v i , having unit variance and non-zero scaled mean ˜ u i as, N ˜ v i (˜ u i ,
1) = 1 √ π exp (cid:32) − (˜ v i − ˜ u i ) (cid:33) (6)where ˜ v i = √ v i β and ˜ u i = √ u i β Note that the differential probability of particles lying within a velocity space v i and v i + dv i is the same as them being within the scaled velocity space ˜ v i and ˜ v i + d ˜ v i , i.e., N v i ( u i , (2 β ) − ) dv i = N ˜ v i (˜ u i , d ˜ v i From the scaled Gaussian variables ˜ v x , ˜ v y , and ˜ v z , the PDF for the sum ofsquares of independent combinations of these variables, χ , can be constructedwhich are known as the noncentral chi-squared distributions [25], f χ = 12 exp (cid:18) − χ + η (cid:19) (cid:18) χη (cid:19) D − I D − (cid:0) √ ηχ (cid:1) (7)They have two parameters, the number of degrees of freedom D =1, 2, or3 depending on the number of independent Gaussian variables used in thesum, χ , and a noncentrality parameter η , given by the sum of squares of the Saurabh S. Sawant et al. respective means of the Gaussian variables (˜ u x , ˜ u y , and ˜ u z ). Note that I γ ( x )is the modified Bessel function of the first kind, defined as, I γ ( x ) = (cid:88) m ≥ m ! Γ ( m + γ + 1) (cid:18) x (cid:19) m + γ The mean, µ χ , and standard deviation, Σ χ , of f χ is given as, µ χ = D + ηΣ χ = (cid:112) D + 2 η ) (8)To derive the PDF of f ξ x , we choose χ = ˜ v x . Since only one independentGaussian variable is required, we have D = 1. The noncentrality parameter forthe given χ is η = ˜ u x . By substituting χ , D , and η into equation 7, we obtain, f ˜ v x = 12 exp (cid:32) − ˜ v x + ˜ u x (cid:33) (cid:18) ˜ u x ˜ v x (cid:19) I − (˜ u x ˜ v x ) (9)and using equation 8, the mean and standard deviation, µ ˜ v x = 1 + ˜ u x Σ ˜ v x = (cid:112) u x ) (10)The distribution function in equation 9 can be scaled by a factor of 0.5 toobtain a distribution f ( ˜ v x / ) by noting that, f ˜ v x d ˜ v x = f ( ˜ v x / ) d (cid:32) ˜ v x (cid:33) which means that the probability of particles having normalized energies withinenergy space ˜ v x and ˜ v x + d ˜ v x remains unchanged. Therefore, we obtain, f ( ˜ v x / ) = 2 f ˜ v x = exp (cid:32) − ˜ v x + ˜ u x (cid:33) (cid:18) ˜ u x ˜ v x (cid:19) I − (˜ u x ˜ v x ) (11)This PDF can be shifted by − ˜ u x / f ( (˜ v x − ˜ u x ) / ) = exp (cid:32) − ˜ v x (cid:33) (cid:32) ˜ u x (cid:112) ˜ v x − ˜ u x (cid:33) I − (cid:16) ˜ u x (cid:112) ˜ v x − ˜ u x (cid:17) (12)This is the PDF of particles’ normalized X -directional energy, ξ x , as seen fromthe following relation. ˜ v x − ˜ u x v x β − u x β = ξ x itle Suppressed Due to Excessive Length 9 Therefore, we can write equation 12 in terms of ξ x as, f ξ x = exp (cid:16) − ( ξ x + u x β ) (cid:17) (cid:18) u x β √ ξ x (cid:19) I − (cid:16) u x β (cid:112) ξ x (cid:17) Finally, the noncentrality parameter can be substituted back where η =˜ u x = 2 u x β to obtain, f ξ x = exp (cid:18) − ξ x + η (cid:19) (cid:18) η ξ x (cid:19) I − (cid:16)(cid:112) ξ x η (cid:17) (13)The aforementioned process of scaling the PDF of f ˜ v x in equation 9 scalesthe mean and standard deviation in equation 10 by a factor of half, whereasshifting it further by − ˜ u x / − u x β , shifts the mean by the same factor butkeeps the standard deviation unchanged. Therefore, the mean and standarddeviation of f ξ x can be written as, µ ξ x = 0 . Σ ξ x = (cid:112) η )2 (14)Next, to obtain the PDF of f ξ y , we follow the same procedure, where χ = ˜ v y and D = 1, however, note that η = u y = 0. For η = 0, the noncentral chi-squared distribution in equation 7 reduces to a chi-squared distribution andfor D = 1, it has the form, f χ = (cid:114) χ π exp (cid:18) − χ (cid:19) (15)as can be verified by first expanding the series of modified Bessel function I − in equation 7, simplifying, and then substituting η = 0. Additionally, by usingthe same equation 8, the constant values of mean and standard deviation ofone and 1.414 can be obtained. By substituting χ = ˜ v y in equation 15, we get, f ˜ v y = ˜ v y √ π exp (cid:32) − ˜ v y (cid:33) (16)To obtain the distribution for ξ y = v y β = ˜ v y /
2, the distribution f ˜ v y can bescaled by a factor of half as, f ξ y = f ( ˜ v y / ) = ˜ v y √ √ π exp (cid:32) − ˜ v y (cid:33) or in terms of ξ y as, f ξ y = f ( ˜ v y / ) = 2 (cid:114) ξ y π exp (cid:0) − ξ y (cid:1) (17) The effect of scaling the PDF by a factor of half reduces the mean and standarddeviation to 0.5 and √
2, respectively. Note that the PDF of f ξ z has the sameform as the PDF of f ξ y because of axial symmetry in the 1-D case.Similarly, to obtain the PDF of ξ = ( v − u ) β , we choose χ = ˜ v x +˜ v y +˜ v z =˜ v , for which D = 3 and η = ˜ u x + ˜ u y + ˜ u z = ˜ u . Note that for the 1-D case,˜ u y = ˜ u z = 0, which results in η = ˜ u x = ˜ u . By substituting χ , D , and η inequation 7, we obtain, f ˜ v = 12 exp (cid:32) − ˜ v + ˜ u (cid:33) (cid:18) ˜ v ˜ u (cid:19) I (˜ u ˜ v ) (18)By following an approach similar to that for f ξ x , the PDF of f ˜ v can be firstscaled by a factor of half to obtain the PDF of (cid:0) ˜ v / (cid:1) and then shifted by − ˜ u / (cid:0) (˜ v − ˜ u ) / (cid:1) = ξ , which can be written in termsof ξ and η = 2 u β as, f ξ = exp (cid:18) − ξ + η (cid:19) (cid:18) ξη (cid:19) I (cid:16)(cid:112) ξη (cid:17) (19)with mean and standard deviation of f ξ given as, µ ξ = 1 . Σ ξ = (cid:112) η )2 (20)This completes our derivation of the PDFs of f ξ x , f ξ y , f ξ z , and f ξ , whichare valid in the equilibrium regions of freestream and downstream. The meanand standard deviation of these PDFs at probes F in the freestream and D in the downstream are listed in Table 1. The theoretical estimates of meanand standard deviation (not shown) agree with values obtained from DSMC-derived PDFs within 1%, as expected from the good agreement in PDFs shownin Fig. 3. Equation No. Parameters Probes F Probes D µ ξ , Σ ξ µ ξ x , Σ ξ x η = 0 µ ξ y , Σ ξ y Table 1
The mean and standard deviation of f ξ , f ξ x , and f ξ y using DSMC.itle Suppressed Due to Excessive Length 11( a ) ( b ) Fig. 3
The theoretically-derived PDFs of f ξ x , f ξ y , and f ξ at ( a ) probe F and ( b ) probe D ,respectively. Legends for ( b ) are the same as ( a ). The DSMC-derived PDFs are not shownas they agree well with theoretical PDFs. Figures 3( a ) and 3( b ) show the PDFs of f ξ x , f ξ y , and f ξ obtained fromthe analytical expressions given by equations 13, 17, and 19 at probes F inthe freestream and D in the downstream, respectively. Since both of theseregions are in equilibrium, the theoretically derived distributions agree verywell with the numerically obtained distributions from DSMC. The differencesin the respective PDFs at two probe locations are due to the differences intheir macroscopic flow parameters. In the freestream, the PDFs of f ξ x and f ξ are nearly symmetric about zero, as seen from Fig. 3( a ), whereas in thedownstream they are almost completely to the right of zero, as seen fromFig. 3( b ). This difference is due to the large noncentrality parameter of η =2 u x, β =87.25 in the freestream, for which the NCCS PDF tends towards aGaussian, whereas the value in the downstream is only η = 2 u x, β =0.3538.Parameter η can be interpreted as the ratio of the bulk flow energy of thegas to its thermal energy defined using the most-probable-speed, β − . In thedownstream, these two energies are expected to be comparable as the kineticenergy of the gas is converted to the thermal energy by the shock.Also, the differences between the PDFs of f ξ x and f ξ are because of thedifferences in degree, D . In the freestream, this differences is only noticeable asa small shift between the symmetric distributions, whereas in the downstreamthe difference in their shapes is quite conspicuous. Note, however, that at bothlocations the values of µ ξ x and µ ξ are such that the normalized viscous stress, τ xx , is zero based on the relation, τ xx β ρ = − ( σ xx − p ) β ρ = − (cid:18) µ ξ x − µ ξ (cid:19) (21)which is also consistent with the DSMC-derived profiles of τ xx (see Ref. [11]).Additionally, in the freestream, Fig. 3( a ) shows that there is a differencebetween the shape of the PDFs of f ξ y and f ξ x , although they both have the same degree, D = 1. This difference is because of the difference in η , whichis zero for the former and 87.25 for the latter. In the downstream, these twodistributions are very close to each other, as shown in Fig. 3( b ), because thevalue of η is closer, i.e. 0 and 0.3538 for PDFs of f ξ y and f ξ x , respectively. Forboth equilibrium regions, the mean values µ ξ y and µ ξ are 0.5 and 1.5, respec-tively, as shown in Table 1, such that τ yy = − ( σ yy − p ) = 0, based on a similarrelation as in equation 21. Furthermore, notice from Table 1 that the standarddeviations Σ ξ and Σ ξ x are much larger in the freestream than downstream,because η >> η . Although they are not directly related to any macroscopicflow parameter, they suggest a larger spread of X -directional and overall nor-malized energies of particles in the freestream, where the bulk flow energy ismuch larger than the thermal energy, as opposed to downstream where thetwo energies are comparable. This makes sense because in the downstream,particles collide more and distribute their energy better than the freestream,which causes a reduction in the spread of their energy distribution. To derive the energy PDFs of f ξ x , f ξ y , and f ξ inside the nonequilibriumregion of the 1-D shock, − < X/λ <
6, we take the approach of Mott-Smith [23], who approximated the PDF inside the shock of particle velocitiesas a linear combination of equilibrium PDFs in the upstream and downstreamequilibrium regions, i.e., nf v = φ f v , + φ f v , (22)where n is the local number density. PDFs of f v , and f v , are obtained bysubstituting the upstream and downstream macroscopic flow conditions of( n , u x, , β ) and ( n , u x, , β ), respectively, in equation 3. The use of normal-ization condition (cid:82) f v d v = 1 leads to the local number density, n = φ + φ .The weights φ and φ are functions of X/λ and are given as [22], φ = n (cid:18) α (cid:16) Xλ (cid:17)(cid:19) φ = n exp (cid:18) α (cid:16) Xλ (cid:17)(cid:19) (cid:18) α (cid:16) Xλ (cid:17)(cid:19) (23)where, α is a constant given as, α = 53 √ π (cid:112) RT tr, u x, u x, − u x, u x, + u x, Using the above equations, at probe P , the prediction of number density, n = φ + φ , is 0.26% higher than the DSMC-computed value, where φ and itle Suppressed Due to Excessive Length 13 φ are 7 . × and 1 . × , respectively. By integrating the differ-ential probability f v dv x dv y dv z over the remaining two components, the PDFof directional velocities, v i for i ∈ { x, y, z } can be obtained inside a shock.By taking moments of it, the macroscopic flow parameters of directional bulkvelocities and directional temperatures can also be obtained as, u i = (cid:90) v i f v i dv i , for i ∈ { x, y, z } T tr,i = mκ b (cid:90) v i f v i dv i (24)Figures 4( a ) and 4( b ) show the respective comparison of the PDFs of nor-malized velocities, v x β and v y β , obtained from DSMC with the theoretical pre-diction of Mott-Smith bimodal model and the equilibrium distribution givenby equation 5 at probe P . Figure 4( a ) shows that the PDF of v x β obtainedfrom DSMC skews towards the high velocity region in contrast to the symmet-ric unimodal equilibrium distribution. The Mott-Smith model, gives a qual-itatively better approximation in terms of bimodality of the PDF, however,it overestimates the contribution of the upstream high velocity particles andshows a second distinct peak. This difference leads to a 1.95% lower predic-tion of the bulk velocity, u x , in comparison to DSMC at probe P , althougha smaller difference is expected further away from the center of the shock,where the degree of nonequilibrium is lower. Figure 4( b ) shows a good agree-ment between Mott-Smith and DSMC except for the region of low Y -velocity,0 . < v y β < .
25, which is overestimated, but nonetheless, the model cor-rectly predicts u y = u z = 0. The predicted values of directional temperatures, T tr,x and T tr,y = T tr,z , are higher and lower than DSMC by 1.35 and 7.2%,respectively, and the predicted overall temperature, computed as the aver-age of directional temperatures, and the value of β are 2.7% lower and 1.4%higher than DSMC, respectively. Since the Mott-Smith model gives a reason-able estimation of macroscopic flow parameters, a similar idea can be used toanalytically derive bimodal NCCS PDFs inside the 1-D shock.We start with normalizing the weights φ and φ by the total weight φ + φ and use these for constructing a bimodal PDF as, f ξ i = ψ f ξ i , + ψ f ξ i , where, i ∈ { x, y, z } f ξ = ψ f ξ, + ψ f ξ, (25)where ψ = φ φ + φ and ψ = φ φ + φ The PDFs of f ξ x , , f ξ y , , and f ξ, are generated from PDFs of ˜ v x , ˜ v y , and ˜ v given in equation 9, 16, and 18, respectively, by first scaling them by a factorof s = 0 . β /β and then shifting by the amounts − u x β , 0, and − u β ,respectively. The local values of u x , u , and β are obtained from equation 24.Similarly, to obtain PDFs of f ξ x , , f ξ y , , and f ξ, the scaling factor of s = . β /β is used, while the amount of respective shifting is the same, and theprocedure for generating f ξ z , and f ξ z , is the same as that for f ξ y , and f ξ y , ,respectively, because of axial symmetry. The final PDFs of f ξ,j , f ξ x ,j , and f ξ y ,j , where j = 1 and 2 correspond to upstream and downstream conditions,respectively, are given as, f ξ,j = f ( v x,j − u ) β = 12 s j exp (cid:32) − v x,j − u ) β s j + η j (cid:33) × (cid:32) v x,j − u ) β s j η j (cid:33) I (cid:115) η j ( v x,j − u ) β s j f ξ x ,j = f ( v x,j − u x ) β = 12 s j exp (cid:32) − v x,j − u x ) β s j + η j (cid:33) × (cid:32) s j η j v x,j − u x ) β (cid:33) I − (cid:115) η j ( v x,j − u x ) β s j f ξ y ,j = f ( v y β ) = 1 s j (cid:114) v y,j β π exp (cid:16) − v y,j β (cid:17) (26)where, v x,j and v y,j are instantaneous velocities of particles that follow PDFsof f v x ,j and f v y ,j in equation 22, for j = 1 ,
2. Also, η j = 2 u x,j β j = 2 u j β j , for j = 1 ,
2. Note that the mean of these PDFs can be easily calculated as, µ ξ = ψ s (3 + η ) + ψ s (3 + η ) − u x β µ ξ x = ψ s (1 + η ) + ψ s (1 + η ) − u x β µ ξ x = ψ s + ψ s (27)however, for higher moments the distributions have to be numerically inte-grated. A simple Python code to generate the bimodal distributions of energyusing Mott-Smith fractions and its moments is given in Appendix A.Figures 4( d ), 4( e ), and 4( c ) show the respective comparison of the PDFsof f ξ x , f ξ y , and f ξ obtained from DSMC with the analytically derived bimodalNCCS PDFs using equation 25, denoted as ‘Mott-Smith’, and the equilibriumNCCS PDFs obtained from equations 13, 17, and 19 at local macroscopic flowconditions for probe P . The difference between PDFs obtained from DSMCand the equilibrium theory is expected due to nonequilibrium. It is small forPDF of f ξ y , indicating that the energy distribution in transverse directionsis not significantly affected. The DSMC-derived PDFs of f ξ and f ξ x exhibitinflection points at ξ = 0 .
49 and ξ x = 1 .
33, which indicate their bimodalnature.In our previous work [11], we have used the DSMC-derived PDF of f ξ x to show that the collisions between particles on two sides of the inflection itle Suppressed Due to Excessive Length 15( a ) ( b )( c ) ( d )( e ) Fig. 4
At probe P , the PDFs of normalized velocities ( a ) v x β and ( b ) v y β , normalizedenergies ( c ) ξ , ( d ) ξ x , ( e ) ξ y = ξ z . The DSMC-derived PDFs [11] are compared with theMott-Smith model and equilibrium theory. The legends shown in ( a ) also apply to ( b )–( e ). In( a ) and ( b ), the PDFs of normalized velocities are obtained by scaling velocity distributionsas, f v i β = β − f v i . point are responsible for the low-frequency fluctuations inside a shock. Fig-ure 4 shows that the Mott-Smith model approximates the bimodality of thePDF reasonably well except for the number of particles in the vicinity of theinflection points. Despite these differences, Table 2 shows that the mean andstandard deviation of the PDF of f ξ x obtained from the Mott-Smith modelare only 0.6% lower and 4.4% higher than the values obtained from the DSMC distribution, respectively, and the mean and standard deviation of the PDF of f ξ x are 1.3 and 5.8% higher than DSMC, respectively, whereas for the PDF of f ξ y , these quantities are 7.8 and 9.2% higher, respectively. We will ignore thesesmall differences and use the analytical PDF in the next section to correlateits shape with the low-frequency of fluctuations obtained from DSMC. Parameters DSMC Mott-Smith µ ξ , Σ ξ µ ξ x , Σ ξ x µ ξ y , Σ ξ y Table 2
At probe P inside a shock, the mean and standard deviation of f ξ , f ξ x , and f ξ y using DSMC and the Mott-Smith model. This section describes the correlation of the analytical bimodal NCCS PDFof f ξ x with the DSMC-derived characteristic average low-frequency, ν , in ar-gon shocks, obtained in our previous work [11]. Using these correlations, anestimate of ν is made for the Mach number range 3 ≤ M ≤
10 and inputtemperature range 89 ≤ T tr, ≤ σ xx inside shocks with the change of shape of the analyticalbimodal NCCS PDF of f ξ x for Mach numbers ranging from three to 10 at T tr, = 710 K. Figure 5( a ) shows the variation of the analytically-derived f ξ x at the location of maximum bulk velocity gradient in the shock ( X/λ =0) with Mach number. At Mach 2, the shape of the PDF is very similar toa regular chi-squared distribution, indicating only a small contribution fromhigh-energy upstream particles. However, as the Mach number increases, theircontribution increases leading to a prominant peak in the high-energy space, ξ x >
1. The DSMC-derived weighted average frequency is found to be linearlyproportional to the location of this peak. This is cosistent with the findingsin our previous work [11] that the value of low-frequency depends on thecollisions of low and high-energy particles from the subsonic downstream andsupersonic/hypersonic upstream, respectively. To determine the exact locationof the peak, the second gradient of the PDF is evaluated, f (cid:48)(cid:48) ξ x , as shown infigure 5( b ) for Mach 2, 4, and 10. At Mach 4 and 10, two inflection pointscan be seen, as defined by f (cid:48)(cid:48) ξ x = 0, whereas at Mach 2, such an inflectionpoint is absent. Between the two inflection points, f (cid:48)(cid:48) ξ x attains a local minima,denoted as ξ Mx , which corresponds to the respective peak location of f ξ x inthe high-energy space. It is seen from figure 5( c ) that the DSMC-derivedweighted-average low-frequency linearly changes with ξ Mx . The data-points itle Suppressed Due to Excessive Length 17 with increasing ξ Mx in figure 5( c ) correspond to Mach 3 to 10, which are alsolisted in Table 3.Note, however, that the shape of f ξ x does not change with changes in theupstream temperature, T tr, , at a given Mach number, because ξ x is normal-ized by β . Therefore, the change in low-frequency due to the variation infreestream temperature cannot be estimated using only the linear fit in fig-ure 5( c ). In this case, DSMC-derived ν is found to be linearly proportionalto the most-probable speed, β − , inside a shock calculated from the analyti-cal Mott-Smith velocity distribution function. Figure 5( d ) shows such a linearvariation at Mach 7.2 for data-points with increasing β − corresponding to T tr, / β − obtained from Mott-Smith with T tr, is also given in figure 5( d ). ( a ) ( b )( c ) ( d ) Fig. 5
At probe P , ( a ) analytical bimodal NCCS PDF, f ξ x , for Mach 2–10 and ( b ) secondderivative, f (cid:48)(cid:48) ξ x , for Mach 2, 4, and 10. Location of ξ Mx in the energy space is indicatedby verticle dashed-lines for Mach 4 and 10. f (cid:48)(cid:48) ξ x is smoothened using the Savitzky-Golayfilter provided in the SciPy [26] software. ( c ) For Mach 3–10, a linear variation of DSMC-derived ν with ξ Mx of the PDFs shown in ( a ). ( d ) At Mach 7.2, for cases with varied inputtemperature, linear variation of DSMC-derived ν in the shock with β − obtained fromMott-Smith bimodal velocity distribution. Also shown the variation of β − with T tr, .8 Saurabh S. Sawant et al. The aforementioned two approaches can be combined to determine ν be-tween Mach 3 to 10 and T tr, between 89 to 1420 K as, ν ( M ∗ , T ∗ tr, ) = ν ( M = 7 . , T ∗ tr, ) − dνdξ Mx [ ξ Mx ( M = 7 . − ξ Mx ( M ∗ )] (28)where M ∗ and T ∗ tr, are the desired Mach number and freestream temperature,respectively, and the slope dν/dξ Mx is obtained from figure 5( c ) to be equalto 36 kHz. For example, using equation 28, the predicted low-frequency at M =4 and T tr, =300 K is 6.4 kHz. The frequency at M =7.2 and T ∗ tr, =300 Kis 20.46 kHz, which is obtained by linear interpolation based on figure 5( d )using the Y -intercept and slope of the fitted line equal to 0.6867 kHz and dν/dβ − = 14 . m − , respectively. In comparison, the peak-frequency in thenormalized PSD of the DSMC-derived mean-subtracted instantaneous dataof σ xx ρ β − is 7 kHz, as shown in figure 6( a ), in close agreement with ourestimate. As discussed in Ref. [11], the PSD is obtained using the Welch’smethod [27,28] in the Scipy [26] software with two Hann-window weightedsegments of data sampled with a frequency of 100 MHz prior to the FastFourier Transform (FFT) such that the frequency resolution is 0.775 kHz.Considering the complexity of this full simulation, equation 28 is a usefulscaling relationship for extension to other Mach numbers and temperatures. M ξ Mx Table 3 ξ Mx for Mach numbers ranging from 3 to 10.( a ) ( b ) Fig. 6 ( a ) At M=4 and T tr, = 300 K, the PSD obtained from the mean-subtracted dataof σ xx ρ β − at the location of maximum bulk velocity gradient ( X/λ = 0). ( b ) Predictedcharacteristic average low-frequencies as a function of input temperature for Mach numbers3 to 10.itle Suppressed Due to Excessive Length 19 Finally, using this approach, the average low-frequencies are predicted forinput temperature ranging from 89 K to 1420 K for Mach numbers 3 to 10and shown in figure 6( b ). It is seen that the chacteristic average low-frequencyincreases with increase in temperature and Mach number, however, the maxi-mum range of frequencies is within 50 kHz for the parameter range examined.This figure can be used to read out the average low-frequency fluctuation for agiven input temperature and shock Mach number for researchers in hypersonictransition. This work demonstrates a simple approach to estimate the low-frequencyfluctuations in a 1-D shock structure of argon, using analytically derived PDFsof particle energies in local equilibrium. The latter compare well with DSMC-derived PDFs in the upstream and downstream equilibrium regions with re-spect to a one-dimensional shock. It is shown for the first time that they havethe form of NCCS PDFs. Then, using the Mott-Smith model, the bimodalNCCS PDFs inside the nonequilibrium regions of shocks are constructed as alinear combination of the equilibrium upstream and downstream PDFs. The bimodal NCCS PDFs at the location in the shock where the veloc-ity gradient is maximum (
X/λ = 0) are correlated with the low-frequenciesof shock fluctuations previously obtained from DSMC for 3 ≤ M ≤
10. Itis found that the weighted-averaged low-frequency is directly proportional tothe location of the peak in the PDF of f ξ x generated by the contribution ofparticles from the upstream. For cases where the Mach number is constantbut input temperature is varied, the low-frequencies are found to be propor-tional to the most-probable-speed at X/λ = 0 obtained from the Mott-Smithvelocity PDF. Using these linear functions, it is demonstrated that one canestimate the low-frequency of fluctuations for any arbitrary input conditionnot explicitly simulated in our previous work. In the future, a similar approachcan be used for other gases of practical importance, where the DSMC-derivedlow-frequencies can be obtained for a few selected input conditions and usingthose, a broader database of low-frequencies can be generated through linearcorrelations with analytical PDFs.The estimates provided for low-frequency unsteadiness generated at theshock can be used in receptivity and linear stabiliy analysis studies of laminar-turbulent transition of high-speed boundry layer flows, e.g. as boundary con-ditions in the analyses or to construct simplified models that account for theinteraction of leading-edge shock with the boundary layer. The same esti-mates can also be used to aid understanding of the changes in the spectrumof freestream noise after it interacts with the shock prior to excitation of theboundary-layer. The Appendix provides a python function to generate these PDFs.0 Saurabh S. Sawant et al.
The research conducted in this paper is supported by the Office of NavalResearch under the grant No. N000141202195 titled, “Multi-scale modeling ofunsteady shock-boundary layer hypersonic flow instabilities.” This work usedthe STAMPEDE2 supercomputing resources provided by the Extreme Scienceand Engineering Discovery Environment (XSEDE) at the Texas AdvancedComputing Center (TACC) through allocation TG-PHY160006.
A Python Code for Generating the Bimodal NCCS PDFs
This appendix gives the code snippet for generating bimodal NCCS distribution. The valuesof code variables ‘R’, ‘Ux1’, ‘Ttr1’, ‘beta1’, ‘Ux2’, ‘Ttr2’, ‘beta2’ are 208.243, 3572.24, 710,1 . × − , 944.74, 12120.6, 4 . × − , respectively. The values of code variables ‘Ux’,‘Ttr’, ‘beta’ are obtained from the Mott-Smith velocity distribution and equal to 2012.8,10149.8, 4 . × − , respectively, at probe P . The code variables ‘psi1’ and ‘psi2’ areMott-Smith fractions and equal to 0.40649 and 0.59351, respectively, at probe P . ∗ s i z e )s i z e 2=i n t ( p s i 2 ∗ s i z e )Var1 = R ∗ Ttr1s1 =0.5 ∗ ( Ttr1 / Ttr ) ∗ beta1Ux1Sq=Ux1 ∗ Ux1df1 = [ 1 , 1 , 3 ]nc1 = [ Ux1Sq/Var1 , 0 , Ux1Sq/Var1 ]Var2 = R ∗ Ttr2s2 =0.5 ∗ ( Ttr2 / Ttr )beta2Sq=beta2 ∗ beta2Ux2Sq=Ux2 ∗ Ux2df2 = [ 1 , 1 , 3 ]nc2 = [ Ux2Sq/Var2 , 0 , Ux2Sq/Var2 ]UxSq=Ux ∗ UxbetaSq=beta ∗ betaMeanSqBetaSq = [ UxSq ∗ betaSq , 0 , UxSq ∗ betaSq ]b i n s =3200bin edges NCCS = np . z e r o s ( b i n s +1)bin center NCCS = np . z e r o s ( b i n s )i =0 − MeanSqBetaSq [ i ] )itle Suppressed Due to Excessive Length 21mean , var , skew , kurt = ncx2 . s t a t s ( df2 [ i ] , nc2 [ i ] ,s c a l e=s2 , moments=’mvsk ’ )r2 = ncx2 . r v s ( df2 [ i ] , nc2 [ i ] , s i z e=s i z e 2 , s c a l e=s2 ,l o c= − MeanSqBetaSq [ i ] )r=np . c o n c a t e n a t e ( ( r1 , r2 ) )hist NCCS , bin edges NCCS = np . histogram ( r , b i n s=bins ,d e n s i t y=True )p r i n t ( ‘ mean , stddev ’ , r . mean ( ) , r . std ( ) )f i g , ax = p l t . s u b p l o t s ( 1 , 1)ax . h i s t ( r , d e n s i t y=True , b i n s=bins , h i s t t y p e =‘ s t e p f i l l e d ’ ,alpha =0.4)p l t . show ( )2 Saurabh S. Sawant et al.
References
1. S.P. Schneider, Progress in Aerospace Sciences (1-2), 1 (2004). DOI 10.1016/j.paerosci.2003.11.0012. M.V. Morkovin, Critical evaluation of transition from laminar to turbulent shear layerswith emphasis on hypersonically travelling bodies. Tech. Rep. AFFDL-TR-68-149, AirForce Flight Dynamics Laboratory (1969)3. J.L. Potter, AIAA Journal (10), 1907 (1968). DOI 10.2514/3.48994. S.R. Pate, AIAA Journal (6), 1082 (1971). DOI 10.2514/3.499195. S.P. Schneider, Journal of Spacecraft and Rockets (3), 323 (2001). DOI 10.2514/6.2000-22056. A. Wagner, E. Sch¨ulein, R. Petervari, K. Hannemann, S.R.C. Ali, A. Cerminara, N.D.Sandham, Journal of Fluid Mechanics , 495–531 (2018). DOI 10.1017/jfm.2018.1327. P. Balakumar, A. Chou, AIAA Journal pp. 193–208 (2018). DOI 10.2514/1.J0560408. C. Hader, H.F. Fasel, Journal of Fluid Mechanics , R3 (2018). DOI 10.1017/jfm.2018.3869. S. Lee, S.K. Lele, P. Moin, Journal of Fluid Mechanics , 225 (1997). DOI 10.1017/S002211209700510710. X. Zhong, Journal of Computational Physics (2), 662 (1998). DOI 10.1006/jcph.1998.601011. S.S. Sawant, D.A. Levin, V. Theofilis. A kinetic approach to studying low-frequencymolecular fluctuations in a one-dimensional shock (2020). arXiv: 2012.1459312. A. Fedorov, A. Tumin, AIAA Journal (7), 2335 (2017). DOI 10.2514/1.J05532613. P. Luchini, in Seventh IUTAM Symposium on Laminar-Turbulent Transition (Springer,2010), pp. 11–1814. L.D. Landau, E.M. Lifshitz,
Statistical Physics: Part 1 Volume 5 , 3rd edn. (PergamonPress, 1980)15. E.C. Marineau, G.C. Moraru, D.R. Lewis, J.D. Norris, J.F. Lafferty, H.B. Johnson, in (2015), p. 1737. DOI 10.2514/6.2015-173716. B. Chynoweth, C. Ward, R. Henderson, C. Moraru, R. Greenwood, A. Abney, S. Schnei-der, AIAA Paper (2014). DOI 10.2514/6.2014-007417. Y. Ma, X. Zhong, Journal of Fluid Mechanics , 31 (2003). DOI 10.1017/S002211200300478618. Y. Ma, X. Zhong, Journal of Fluid Mechanics , 79 (2003). DOI 10.1017/S002211200300479819. Y. Ma, X. Zhong, Journal of Fluid Mechanics , 63 (2005). DOI 10.1017/S002211200500383620. N.M. Laurendeau, Statistical Thermodynamics: Fundamentals and Applications (Cam-bridge University Press, 2005)21. W.G. Vincenti, C.H. Kruger,
Introduction to Physical Gas Dynamics (Wiley, New York,1965)22. G.A. Bird,
Molecular Gas Dynamics and the Direct Simulation of Gas Flows , 2nd edn.(Clarendon Press, 1994)23. H.M. Mott-Smith, Phys. Rev. , 885 (1951). DOI 10.1103/PhysRev.82.88524. F. Sharipov, Rarefied Gas Dynamics: Fundamentals for Research and Practice (JohnWiley & Sons, 2015)25. S. Andr´as, ´A. Baricz, Journal of Mathematical Analysis and Applications (2), 395(2008). DOI 10.1016/j.jmaa.2008.05.07426. SciPy-1.5.1, https: // docs. scipy. org/ doc/ scipy/ reference/ generated/ scipy. signal. welch.html
27. P. Welch, IEEE Transactions on Audio and Electroacoustics15