Newly identified principle for aerodynamic heating in hypersonic flows

Instability evolution in a transitional hypersonic boundary layer and its effects on aerodynamic heating are investigated over a 260 mm long flared cone. Experiments are conducted in a Mach 6 wind tunnel using Rayleigh-scattering flow visualization, fast-response pressure sensors, fluorescent temperature-sensitive paint (TSP) and particle image velocimetry (PIV). Calculations are also performed based on both the parabolized stability equations (PSE) and direct numerical simulations (DNS). Four unit Reynolds numbers are studied, 5.4, 7.6, 9.7 and $11.7\times 10^{6}~\text{m}^{-1}$ . It is found that there exist two peaks of surface-temperature rise along the streamwise direction of the model. The first one (denoted as HS) is at the region where the second-mode instability reaches its maximum value. The second one (denoted as HT) is at the region where the transition is completed. Increasing the unit Reynolds number promotes the second-mode dissipation but increases the strength of local aerodynamic heating at HS. Furthermore, the heat generation rates induced by the dilatation and shear processes (respectively denoted as $w_{\unicode[STIX]{x1D703}}$ and $w_{\unicode[STIX]{x1D714}}$ ) were investigated. The former item includes both the pressure work $w_{\unicode[STIX]{x1D703}1}$ and dilatational viscous dissipation $w_{\unicode[STIX]{x1D703}2}$ . The aerodynamic heating in HS mainly arose from the high-frequency compression and expansion of fluid accompanying the second mode. The dilatation heating, especially $w_{\unicode[STIX]{x1D703}1}$ , was more than five times its shear counterpart. In a limited region, the underestimated $w_{\unicode[STIX]{x1D703}2}$ was also larger than $w_{\unicode[STIX]{x1D714}}$ . As the second-mode waves decay downstream, the low-frequency waves continue to grow, with the consequent shear-induced heating increasing. The latter brings about a second, weaker growth of surface-temperature HT. A theoretical analysis is provided to interpret the temperature distribution resulting from the aerodynamic heating.

Newly identified principle for aerodynamic heating in hypersonic flows

Introduction
Aerodynamic heating is a key issue in hypersonic flow research due to its strong relevance to the safety of ultra-high-speed flight. Accurate knowledge about where and why, and with what strength, strong heating peaks occur at the wall is of crucial importance to the performance and safety of hypersonic vehicles (Sivolella 2014). Under real flight conditions, laminar-to-turbulence transition is one of the most important sources bringing about uncertain aerodynamic heating that might adversely impact a vehicle's thermal protection system (TPS) Greene & Hamilton 2006;McGinley et al. 2006). It is generally acknowledged that enhanced heat transfer is related to the evolution of flow disturbances and onset of turbulence, which in turn modifies the skin friction. Extensive studies have been performed using physical and numerical experiments as well as theoretical analyses (Hopkins & Inouye 1972;Brostmeyer & Nagamatsu 1984;Horvath et al. 2002;Wadhams et al. 2008;Berridge et al. 2010;Franko & Lele 2013). In fact, the region with a dramatic increase of skin friction and heat transfer (denoted herein as HT) is commonly identified as the onset of the turbulent state in experiments and numerical simulations.
An important mechanism in the transition process, absent in low Mach number flows, is the instability of dilatational, longitudinal waves (Emanuel 1992;Gad-el-Hak 1995), such as second-and higher instability modes, which was first identified based on linear stability analysis of hypersonic boundary layers (Mack 1969). Those modes behave as acoustic waves reflecting between the solid wall and the sonic line in hypersonic boundary layers (Fedorov 2011), and play an increasingly important role as the Mach number increases -as demonstrated in previous studies (Stetson & Kimmel 1992;Bountin, Shiplyuk & Sidorenko 2000;Fujii 2006;Hofferth et al. 2013;Zhang, Tang & Lee 2013;Casper, Beresh & Schneider 2014;Sivasubramanian & Fasel 2015;Zhang et al. 2015;Zhu et al. 2016). Recent investigations of hypersonic boundary layers have indicated the appearance of an additional peak in heat transfer (denoted herein as HS) within the transitional region, as well as a second rapid growth of heat transfer near the end of transition (HT). Schneider and his group (Berridge et al. 2010) were the first to observe this phenomenon on a Mach 6 flared cone under quiet flow conditions. Since 1998, this group has made significant contributions to re-establish a hypersonic quiet wind tunnel at Purdue University. They also pioneered the application of temperature-sensitive paint (TSP) in hypersonic experiments. In their TSP results published in 2010 (Berridge et al. 2010), hot streaks were found during the early stage of transition over a flared cone, corresponding to a peak value HS in the spanwise-averaged temperature. PCB piezoelectric pressure sensors were also applied within that region, and indicated an appearance of second-mode instability. Sivasubramanian & Fasel (2015) used direct numerical simulations (DNS) to simulate the experiment of Berridge et al. (2010). Streamwise hot streaks were also observed before the boundary layer became turbulent. Either spanwise-averaged heat transfer or that along a single streak exhibited both HS and HT along the streamwise direction. Note that the Stanton number at HS was twice that at HT, but the skin friction of the former was only half that of the latter, so that a new aerodynamic heating mechanism might be active. However, because the focus of Sivasubramanian & Fasel's work was the fundamental resonance of second-mode instability, no explanation of HS has been presented in their paper.
The behaviour of hypersonic aerodynamic heating seems to rely on the types of instability mode. Based on DNS, Franko & Lele (2013) investigated the transition mechanisms of a Ma = 6 planar hypersonic boundary layer. Three types of such 154 Y. Zhu, C. Lee, X. Chen, J. Wu, S. Chen and M. Gad-el-Hak mechanisms were identified, including first-mode oblique interaction, second-mode fundamental resonance and second-mode oblique interaction. For the first-mode oblique interaction mechanism, a distinct overshoot of heat transfer was observed near the end of the transitional region, HT, which was attributed to the generation of streamwise vortices and correspondingly increased wall shear stress. However, for the fundamental resonance dominated by two-dimensional second mode, an additional peak of heat transfer, HS, appeared where the second mode reached its maximum. The temperature values at HS are slightly higher than those at HT, while the associated skin friction of the former is no more than half of the latter. Franko & Lele (2013) stated that the appearance of HS was attributed to the second-mode evolution, but no physical explanation was given in terms of fundamental thermo-aerodynamics.
Understanding the physics of second-mode instability is clearly important. By applying inviscid theory to compressible boundary layer stability, Mack (1969) found a family of neutral curves, denoted as second and higher modes, whenever the phase speed of the disturbances is subsonic relative to the free stream. As explained by Morkovin (1987) and later by Mack (1990), the inviscid solutions of these higher modes were consistent with the acoustic behaviour depicted by Lighthill's equation. Mack (1990) suggested the terminology 'acoustic modes' to describe those modes because they were dominated by acoustic waves. When the waves reflected in an embedded region between the sonic line and solid wall, the flow was alternately compressed and expanded owing to the boundary condition at the sonic line. Not surprisingly, the consequent variations of density or velocity divergence could be captured experimentally using techniques such as schlieren or particle image velocimetry (PIV). Hofferth et al. (2013) have applied a focused schlieren technique on a 0.5 m long, 5 • half-angle flared cone in their Mach 6 Quiet Wind Tunnel at Texas A&M University. The density variation caused by second-mode instability resulted in light intensity variations on the schlieren image, which were acquired by a high-bandwidth photo-detector. Using this technique, Hofferth et al. (2013) determined the frequency and amplitude evolution of the second mode. Laurence, Wagner & Hannemann (2012) investigated a 1.1 m long, 7 • half-angle cone in a free-piston-driven, reflected-shock wind tunnel under high-entropy conditions. Employing a 1000 W short-arc Xe lamp together with a camera at frame rates of up to 1 MHz, Laurence et al. (2012Laurence et al. ( , 2016 obtained schlieren images of typical 'rope-like' structures of the second mode. The intensity variation in the image manifested the density variation caused by the compression and expansion process inherent in the second-mode instability. Very recently, Zhang et al. (2015) and Zhu et al. (2016) have conducted twodimensional PIV measurements, combined with PCB sensors, over a 0.26 m long, 5 • half-angle flared cone in the Mach 6 wind tunnel at Peking University. They found that the second-mode instability initially grew fast, reached its maximum and finally decayed to zero along the streamwise direction. However, low-frequency waves continued to grow until the boundary layer became turbulent. The PIV results clearly showed a strong compression and expansion process accompanying the evolution of second mode. The corresponding shear-and dilatation-induced viscous dissipation of energy were respectively evaluated; however, their connection to aerodynamic heating was not investigated. They then applied temperature-sensitive paint, as well as PCB sensors, to simultaneously measure both the surface-temperature distribution and the instability amplitude evolution. PIV measurements were also performed to evaluate the dissipation function of the flow. The recent work not only confirmed the previously observed heat-transfer peak HS that correlated with the second-mode-dominated Newly identified principle for aerodynamic heating in hypersonic flows 155 instability evolution, but also revealed a new physical mechanism to promptly generate this heating impact and its strong connection to the instability evolution. For the first time, a direct relation was identified between the local heating peak HS before transition was completed and the second-mode-induced dilatation process. The new discovery has been published as a rapid communication (Zhu et al. 2018), and soon after highlighted by Wolverton (2018) and Sun & Oran (2018) for its significant achievement from the physics perspective as well as in the engineering design of hypersonic flights. This paper is an expansion of the preliminary work reported in Zhu et al. (2018), methodically providing more experimental and numerical evidence. Compared to the rapid communication, the present expanded manuscript describes an in-depth investigation of the newly discovered heating mechanism. The importance of pressure work in aerodynamic heating is particularly elucidated. A theoretical analysis of heat transfer is also described to explain the relationship between aerodynamic heating and surface-temperature distribution.

Facility and model
The experiments are carried out in a Mach 6 wind tunnel (M6QT) at Peking University Tang 2014;Zhang 2014;Zhang et al. 2015;Zhu et al. 2016), which at present is one of three operational hypersonic quiet wind tunnels in the world (Schneider 2013). The tunnel is of the open-jet configuration with nozzle exit diameter of 160 mm. Longitudinal suction can be applied upstream of the nozzle throat and the wind tunnel runs at quiet or noisy conditions with the suction valve open or closed, respectively. For the quiet condition, the turbulence level is below 0.2 %. Owing to low free-stream disturbances and the upper limit of unit Reynolds number for the quiet condition, natural transition does not occur up to the end of the present model. Therefore, to obtain a whole picture of the transition process in the present work, the wind tunnel runs at noisy condition with the suction valve closed. To avoid liquefying the air, the flow is pre-heated to a nominal stagnation temperature of 433 K. During the typical run time of 30 s, the stagnation pressure holds nearly constant with variation of less than 3 %. Four flow conditions with Re unit = 5.4, 7.6, 9.7 and 11.7 × 10 6 m −1 are investigated, where Re unit is the free-stream unit Reynolds number. As indicated in Zhu et al. (2016), the stagnation temperature and normalized pitot pressure fluctuations of the free stream are, respectively, 433 K and 2.2 %.
The model used in the present study is a flared, instability-enhanced cone, schematically shown in figure 1. The full length is L = 260 mm. Its geometry consists of a 5 • half-angle circular conical profile for the first 100 mm of axial distance, followed by a tangent flare of radius 931 mm until the base of the cone at the 260 mm axial position. The first 50 mm of length, which can be removed and replaced with an arbitrary nose-tip shape, is made of stainless steel. In this work, the tip is sharp with a nominal 50 µm radius. The rest of the model is made of Bakelite. The origin of the coordinate is located at the cone tip, with x being the streamwise coordinate along the cone's surface, y the coordinate normal to that surface, z the transverse coordinate normal to the x-y plane and l the axial coordinate along the axis of the model.
The cone is installed along the centreline of the nozzle with zero angle of attack. The tip is 50 mm deep into the nozzle exit. The length of the flow field that is not affected by the reflected Mach waves is more than 400 mm, which is long enough to embed the whole length of the model. 2.2. Fast-response surface pressure measurements An array of surface-mounted PCB fast-response sensors is used to evaluate the evolution of instability waves in hypersonic wall-bounded flows (Fujii 2006;Zhang et al. 2013Zhang et al. , 2015Zhu et al. 2016). In this work, PCB 132A31 sensors are flush-mounted along one ray of the model, which is defined as the zero circumferential angle. The sensor is a piezoelectric quartz sensor with a flat frequency response between 11 KHz and 1 MHz. The diameter of the sensor head is 3.18 mm, but the effective sensing area is only approximately 0.581 mm 2 . The signals of the PCB sensors are first processed by an ICP signal conditioner and then recorded by Donghua 5939 data-acquisition system with a sample rate of 1 MHz and a duration of 100 ms.

Temperature-sensitive paint
To investigate the surface-temperature distribution, fluorescent temperature-sensitive paint (TSP) (Liu & Sullivan 2005) is painted all over the model surface (figure 1), with a thickness of less than 0.1 mm. The fluorescence is excited continuously by an LED light and recorded by a Photron FASTCAM CMOS high-speed camera with an image resolution of 1024 × 1024 pixels and sample rate of 100 Hz. A Nikkor Micro 60 mm lens is mounted ahead of the camera using long-wavelength-pass filter to remove the LED light. The camera's recording rate is 100 Hz, corresponding to an integral time of 10 ms on each complementary metal oxide semiconductor (CMOS) pixel. A black background image is sampled before recording and subtracted from each image.
The model is first placed in a thermal bath to calibrate the TSP's performance. Fluorescent images are then recorded through a quartz-glass window once the model Newly identified principle for aerodynamic heating in hypersonic flows 157 reaches thermal equilibrium. Intensity ratios are then calculated under different temperatures to calibrate the TSP images and the results are depicted in figure 1.
The TSP's response time is affected by a number of factors, including the fluorescence lifetime, the image exposure time and the paint's heat lag time. The first two factors are 1 µs and 0.01 s, respectively. The last factor can be evaluated from where δ, ρ, C and k are, respectively, the thickness, density, specific heat and thermal conductivity of the paint layer. For our TSP, δ = 0.1 mm, ρ = 1200 kg m −3 , C = 1464 J kg −1 K −1 and k = 0.19 W m −1 K −1 . Thus, the heat lag time is estimated to be 0.1 s, which is considered to be the TSP's response time.
2.4. Particle image velocimetry Due to rapid developments in CCD cameras since the first application of PIV in hypersonic flows (Moraitis & Riethmuller 1988), this technique is now more commonly used in high-speed flows (Scarano 2005 The key issues include the selection of seeding particles, their dispersion in the flow and the analysis of particle image recordings (Scarano 2005).
To seed the present flow, TiO 2 particles are applied for their high-temperature stability and high index of refraction. The nominal diameter of the selected particles is 90 nm. The particle tracking capability is assessed by Scarano (2005) as the ratio where τ p is the particle's relaxation time as it passes through a velocity step (e.g. an oblique shock wave) from u 1 to u 2 , where u 1 and u 2 are the initial and final velocities during a single particle's displacement. The particle's velocity is given by where t is the time. The flow time scale τ f is the slipping time of the flow, defined as where du is the velocity spatial variation, and dl is the distance for that variation.
In the present work, τ p is assessed experimentally to be 2 µs (Tang 2014). For a hypersonic boundary layer with thickness δ, dl is characterized by the wavelength of the second-mode wave λ, and du is no more than 5 % of the free-stream velocity U ∞ (Stetson et al. 1983;Stetson & Kimmel 1992). Thus, τ f λ/(0.05U ∞ ). In view of our experimental conditions, U ∞ = 850 m s −1 , and λ 2.5 mm. Therefore, τ f 58.5 µs, and τ 0.034. As discussed in Zhu et al. (2016), the injection of PIV particles somewhat dampens the amplitude of the high-frequency second mode. Therefore, the velocity fluctuations below the sonic line are slightly decreased.

158
Y. Zhu, C. Lee, X. Chen, J. Wu, S. Chen and M. Gad-el-Hak The particles are illuminated by a double-cavity Nd:YAG laser from Continuum, generating 6 ns duration light pulses at a wavelength of 532 nm, with a maximum energy of 2 J per pulse. The laser sheet is projected vertically from the top window of the test section and aligned with the cone's top centreline. The light sheet's thickness is 1 mm and width 100 mm. The time delay of the laser pulses is set at 1 µs and the sample rate is 2.5 Hz. The light scattered by the particles is viewed by a PCO SensiCam QE CCD camera from a side window, which is equipped with a Nikkor Micro 200 mm lens. The field of view is 34 × 26 mm 2 .
The algorithm to evaluate the PIV images has been discussed by Zhu et al. (2013). To conduct window correlation of the particle images, at least 10 particles are needed in each interrogation window (Adrian 1991). However, this is not always satisfied in laminar boundary layers. Furthermore, the low signal-to-noise ratio of particle images can decrease the Q-value of the correlation peak, making it submerged in the noise, or provides wrong peaks. Therefore, only the images that satisfy both standards are selected.
Finally, to solve the problem of very large in-plane displacement and its gradient on the PIV images -commonly encountered in hypersonic boundary layers -we have developed an improved image-preprocessing method (Zhu et al. 2013). PIV images are first rotated to make the wall boundary nearly horizontal. The evaluations start in 256 × 256 pixel coarse samples and then 128 × 128 pixel medium samples to capture the mean displacement fields of the time series. With the mean results as reference, instantaneous recordings are then evaluated with 64 × 64 pixel samples to resolve the fine structures in the boundary layer. Smaller samples are unavailable because of low particle density near the wall. Vectors are then post-processed by a global histogram filter and a Q-ratio of 2. Vectors could then be replaced if a second correlation peak proved to be suitable.

Rayleigh-scattering technique
The whole transition process is also investigated qualitatively using Rayleigh-scattering flow visualization (RSFV), which can offer clear information regarding the evolution of the structures in the boundary layer. This method was used by Smits & Lim (2000) and Auvity, Etz & Smits (2001) to investigate various hypersonic flows. CO 2 gas is injected upstream of the test section. The mass injection rate of carbon dioxide is no more than 5 % of the free-stream flow. Owing to the very low static temperature of the hypersonic flow, the CO 2 gas changes phase to become minute, solid particles in the free stream. The particles scatter when illuminated by the laser, so that this region is white on the grey-scale CCD image. The scattering is of the Rayleigh type because the particle diameters are much smaller than the laser wavelength. On the other hand, these CO 2 particles sublimate to gas near the wall because of the relatively high temperature there. This region is black on the CCD image. The line distinguishing the white and black regions is the condensation line. The set-up of the light source and camera is the same as in the PIV, except that the camera operates in the single-frame mode.
3. Theoretical and numerical approaches 3.1. Analysis based on parabolized stability equations We seek solutions in the form Newly identified principle for aerodynamic heating in hypersonic flows 159 where x, y and φ are body-oriented coordinates representing, respectively, streamwise, wall-normal and circumferential directions;Ψ ≡ (ρ,ū,v,w,T) denotes the basic state; Ψ mn ≡ (ρ,û,v,ŵ,T) mn is the shape function for the Fourier mode (m, n); α mn and m are, respectively, the associated streamwise (complex) wavenumber and radian frequency. All variables are non-dimensionalized by the edge value of the mean flow. Finally, n is an integer called the azimuthal wavenumber, which represents the number of waves around the circumference of the cone. The mode shape,Ψ mn , is assumed to be a slowly varying function of the streamwise coordinate x, so that its second derivative with respect to x is negligible. Substituting (3.1) into the Navier-Stokes equations and subtracting the mean-flow equation yields the governing equations for the shape function of each mode (m, n): where F mn is the nonlinear forcing term, and expressions of the operators A and B can be seen in appendix A.
The above system of equations is parabolic and thus requires boundary and initial conditions. The initial conditions are provided by the linear stability theory (LST). The boundary conditions arê  (Chang & Malik 1994;Li & Malik 1996). Both the curvature effects and external pressure-gradient effects have been taken into account. The wavenumber is updated at each station as follows:

160
Y. Zhu, C. Lee, X. Chen, J. Wu, S. Chen and M. Gad-el-Hak Iteration continues until the largest change is less than 10 −5 . A norm based on both kinetic and internal energy is also utilized, with which the parabolized stability equations (PSE) results have no distinct difference compared to the present ones.
In the present work, the mean flow for PSE calculations is obtained via an Euler boundary layer approach (Pruett & Chang 1998). The Euler equations are solved with FLUENT to obtain the inviscid flow field around the flared cone, and the non-similar boundary layer is then calculated using a spectral collocation method (Malik & Spall 1991). An adiabatic wall condition is used.
In the PSE simulations, a planar second mode ( f 2nd , n), where f 2nd = 350 kHz and n = 0, with an initial amplitude (maximum of temperature root-mean-square, r.m.s.) of 5 %, is seeded at x = 0.384L, which is the beginning of the flared region. The initial second-mode amplitude is chosen such that the second mode becomes saturated at approximately x = 0.65L, as observed in the experiments. In order to trigger early transition, an oblique low-frequency mode ( f low , n), where f low = 20 kHz and n = 50, with an amplitude of 0.1 % is also seeded. The initial mode shapes are provided by LST. The number of grid points along the wall-normal direction is 151. Four modes, the mean-flow distortion, the difference mode ( f 2nd − f low , −n), the summation mode ( f 2nd + f low , n) and the harmonic (2f 2nd , 0), are directly triggered via nonlinear interactions. These four modes along with the two seeded modes further generate higher-order modes of which nine modes are included here:

Direct numerical simulations
DNS has been performed on the supercomputer Tianhe-2 in the Guangzhou Supercomputing Centre. The basic laminar (steady) flow is first established including the bow shock and without any perturbation, and is then perturbed to obtain the time-dependent simulations of the transitional flow. The code Hoam-OpenCFD, developed by Li, Fu & Ma (2008), is used in these simulations. In this code, the Navier-Stokes equations (see appendix B) are solved numerically using a high-order finite-difference method.
(i) Convection terms are split by using Steger-Warming splitting, and are discretized with a seventh-order weighted essentially non-oscillatory (WENO) scheme. (ii) An eighth-order central finite-difference scheme is used for the viscous terms. (iii) A third-order total-variation-diminishing type of Runge-Kutta method is used for time advance.
The computational domain of the unsteady (transition) simulations in the circumferential direction spans from the 0 • meridian plane to the 45 • plane. The mesh number (streamwise × wall normal × circumference) is 4000 × 200 × 500, as shown in figure 2. Both the curvature effects and external pressure-gradient effects have been taken into account. The boundary conditions for simulations of the unsteady flow are as follows.  The amplitude , measured in terms of peak streamwise velocity fluctuations across the boundary layer, is equal to 1 % for each mode; the mode shapes, Ψ k ≡ (ρ,û,v,ŵ,T) k , are provided by LST; β k and k are, respectively, the corresponding spanwise wavenumber and radian frequency; and z is the spanwise coordinate in numerical space. The detailed parameters of the seeded modes are illustrated in table 1. (ii) A non-reflecting boundary condition is used on the outflow boundary in the streamwise direction. (iii) A periodic boundary is used in the circumference direction. (iv) For the wall boundary, an adiabatic wall temperature together with the assumption ∂p/∂y = 0 is used on the wall, and a second-order, one-side, finite-difference scheme is used to compute the wall pressure.  For the lowest unit Reynolds number, the second mode continuously grows until the end of the model whereas the low-frequency waves grow slowly. The amplitude of the latter is approximately 25 % that of the former. As Re unit increases to 7.6 × 10 6 m −1 , a saturation of second-mode amplitude occurs at x/L = 0.73. It then decays until the end of the model. Meanwhile, the low-frequency waves enjoy a faster growth, reaching a maximum value nearly the same as that of the second mode. As Re unit increases to 11.7 × 10 6 m −1 , the decay of the second mode is brought forward to x/L = 0.65, whereas the growth of the low-frequency waves is further promoted. For the highest unit Reynolds number, the annihilation of the second mode is further brought forward. By contrast, the maximum value of the low-frequency waves is nearly three times that of the second mode.
Comparison of instability growth is further made between PCB, PSE and DNS, as shown in figure 4. Herein, −α i represents the growth ratio, which can be calculated from the equation where A is the amplitude and α i is the imaginary part of the streamwise wavenumber α. As shown, DNS and PCB agree very well with each other whereas the PSE data, especially for low-frequency waves, are higher than the other two. As discussed by Chen et al. (2017), one possible reason is that PSE purely considers the two fastest-growing modes whereas both the experiment and DNS involve a considerable level of random disturbance in the oncoming stream. The latter two growth ratios are integrative results of a wide spectrum. For example, although the present PSE has considered three-dimensional waves, only the mode (20, 50) with frequency of 20 kHz and circumferential wavenumber of 50 is seeded, which might not completely  satisfy the complicated real flow condition. Despite this, the physical process depicted by PSE is substantially the same as the other two. Surface-temperature distributions are then calculated from the TSP results. As shown in figure 5, increasing the unit Reynolds number greatly promotes the aerodynamic heating as well as the appearance of local hot spots (indicated as HS) in the streamwise direction. For Re unit = 5.4 × 10 6 m −1 (figure 5a), the temperature rises slowly to the end of the model. However, no prominent local peak is observed in the temperature rise distribution. As Re unit increases to 7.6 × 10 6 m −1 (figure 5b), the aerodynamic heating becomes more pronounced, and a hot spot HS appears at x/L = 0.79. As Re unit increases to 9.7 × 10 6 m −1 (figure 5c), the first hot spot HS is advanced to x/L = 0.68 with a new hot spot HT appearing at x/L = 0.89. The latter is related to the onset of turbulence at the end of the transition process (Berridge et al. 2010;Franko & Lele 2013;Sivasubramanian & Fasel 2015). For the highest Re unit (figure 5d), the hot spot HS at x/L = 0.68 becomes stronger, with its peak value higher than that in HT. Streamwise streaks are also seen in the temperature rise distribution, consistent with the findings of Berridge et al.  Newly identified principle for aerodynamic heating in hypersonic flows 165 value varies along the spanwise position, indicating the three-dimensional development of the boundary layer. HS might not be observed clearly if the temperature probes had been arranged at the valleys. Evolution of second-mode and low-frequency waves arising from the PCB array at φ = 0 are compared with the TSP-measured surface-temperature rise along ten circumferential angles from φ = −18 • to 18 • , as shown in figure 6. Temperatures within PCB sensors are evaluated by cubic-spline interpolation. It can be clearly seen that the appearance of HS is strongly correlated with that of the second-mode peak. For the lowest Re unit (figure 6a), the second mode grows until the end and no HS is found. For Re unit = 7.6 × 10 6 m −1 (figure 6b), the second-mode peak is located at x/L = 0.71, which is approximately 0.7 L ahead of HS. For Re unit = 9.7 × 10 6 m −1 (figure 6c), the second-mode peak is located at x/L = 0.65, which is approximately 0.6-0.7 L ahead of HS. For the highest Reynolds number (figure 6d), the distance between the second-mode peak and HS decreases to 0.3 L. Comparison of instability evolutions and spanwise-averaged surface temperature is also made for PSE (figure 6e) and DNS (figure 6f ). Correlations between the peak of the second mode and the hot spot HS are also observed.
In general, as the unit Reynolds number increases, the second mode significantly decreases but HS quickly increases. This relationship seems logical because the aerodynamic heating should come from the dissipation of the flow mechanical energy. The heat production mechanism and its relationship with the flow structure near HS need to be further investigated.

Heat generation rate
An RSFV picture of flow structures and the DNS-calculated temperature field near HS are both presented in figure 7. Such a comparison is reasonable because the RSFV actually reflects the temperature distribution in the field of view (Chen et al. 2017). As shown, the boundary layer is initially laminar in the upstream region, but then regular rope-like structures appear at x/L = 0.55, which is typical of second-mode instability waves because their wavelength is twice the boundary layer thickness. The second-mode waves grow, saturate and finally decay to zero around x/L = 0.76. For experiments, the boundary layer becomes completely turbulent soon afterwards at approximately x/L = 0.85, whereas for DNS, onset of turbulence is delayed until nearly the end of the model (x/L = 1.0). Figure 8 presents the instantaneous flow field near HS. The velocity magnitude normalized by the free-stream velocity, the dilatation ω 2 and θ 2 are presented in, respectively, (a), (b) and (c). Good agreement has been achieved between the PIV and DNS results. As shown, a strong compression-expansion process accompanies the evolution of the second-mode instability (figure 8b), reflecting the acoustic nature of the latter (Mack 1990). The vorticity distribution both near the wall and the critical layer are also modified due to viscous effects, creating periodic vortical waves in the boundary layer (figure 8c). Such a mechanism has been discussed in detail in Zhu et al. (2016).
Time-averaged skin friction can then be calculated from PIV and DNS results. It would seem reasonable that the increase of wall skin friction due to the disturbance evolution is the major source of additional aerodynamic heating in the boundary layer. However, both for DNS and PIV (figure 9), the increase of skin friction at HS from laminar state is less than 15 % of that at HT, which cannot account for the enhanced heat transfer observed in figure 5. The inherent aerodynamic heating mechanism for HS needs to be further investigated.  6. Comparison of second-mode (red) and low-frequency waves (blue) with streamwise evolution of temperature at different spanwise angles (grey) and their spanwise average (green). Obtained from time-averaged PCB spectra and TSP results for (a) Re unit = 5.4 × 10 6 m −1 , (b) Re unit = 7.6 × 10 6 m −1 , (c) Re unit = 9.7 × 10 6 m −1 and (d) Re unit = 11.7 × 10 6 m −1 . Based on (e) DNS, and ( f ) PSE for Re unit = 9.7 × 10 6 m −1 . The subscripts 2nd and lf respectively represent the second-mode and low-frequency waves.
The heat generation rate per unit volume can be explicitly split into contributions of vorticity ω and dilatation θ (Lele 1994; Wu, Ma & Zhou 2006): where BI represents a boundary integral; w θ and w ω denote, respectively, the work done by dilatation and shear stress; µ = µ b + 4µ/3 is the longitudinal viscosity, with µ and µ b being, respectively, the shear and bulk viscosities. As stated above, the work done by the wall skin friction, reflected in the term w ω , is the major source of aerodynamic heating in low Mach number wall-bounded flows. For high Mach number flows, however, the term w θ , including both the pressure work and dilatational viscous dissipation, might be considerable. Its important role was recognized by Emanuel (1992) and Gad-el-Hak (1995). Associated dilatational viscous dissipation was highlighted in our previous work (Zhang et al. 2015;Zhu et al. 2016). Because of their different natures, we further split w θ into two parts: (4.6) w θ2 = µ θ 2 .
(4.7) Figure 10 presents the instantaneous streamwise evolution of w θ1 , w θ2 and w ω in the boundary layer near HS. As shown, the term w θ1 is alternatively positive and negative, representing an alternative creation and absorption of heat. Its time-averaged value not only depends on the amplitude of p and θ, but also the phase angle between them. The terms w θ2 and w ω are always positive because of their viscous dissipation nature. This means heat is created wherever the shear or dilatation process exists.
We next evaluate the evolution of time-averaged w θ1 , w θ1 and w ω in the flow field. Two problems are encountered, however. One is that w θ1 is difficult to obtain experimentally because a simultaneous measurement of the instantaneous pand θ-field is far beyond the measurement capabilities available today for hypersonic flows. To solve this problem, we assumed that p is nearly constant in the boundary layer and w θ1 is mainly created near the wall because p and θ reach their maximum values there (Morkovin 1987). Due to the no-slip condition, a balance between the normal and shear stress is achieved by The second problem is the exact value of µ b . The dilatational dissipation for a gas arises from an internal molecular relaxation process, which is, of course, a dispersion of frequency. As interpreted by Landau & Lifshitz (1959), the bulk viscosity that leads to attenuation of dilatational waves increases monotonically with frequency. A more accurate evaluation of the bulk viscosity coefficient in hypersonic boundary layers is necessary, which is an issue of crucial importance to the community. In this work, we use the most conservative evaluation as µ b = 2µ/3 for N 2 gas at 430 K. Figures 11 and 12 present the time-averaged w θ1 , w θ1 and w ω in the x-y plane near HS and HT, respectively. The values are normalized by their respective maximum value. PIV and DNS results agree reasonably well, which validates our previous assumptions. For HS, the contribution to aerodynamic heating mainly comes from the dilatation process, especially w θ1 , rather than its shear counterpart w ω . The value of w θ2 is also considerable in that region. The main energy of w θ1 and w θ2 lies at the bottom of the boundary layer, because the dilatation process mainly exists beneath the sonic line. Some tiny negative regions of w θ1 are also seen in DNS, but not in PIV because the evaluated pressure there might be far from reality. For HT, the aerodynamic heating mainly comes from the fast increase of w ω very near the wall, which agrees well with the explanation given by Franko & Lele (2013) for HT in that ultimate transition process. Time-averaged w θ 1 , w θ2 and w ω are then integrated along the y-direction in the boundary layer and plotted in figure 10, as well as their correlation with the instability amplitudes. The values are normalized by their respective maxima. It is clearly seen both in PIV and DNS that the dilatation-induced heat at HS, mainly from w θ1 , is more than five times that from w ω . The value of the latter increases to nearly 75 % that of w θ1 at HS for PIV but only 52 % for DNS because transition is not completed in DNS. Note that low-frequency waves continuously grow until the end of the model because it belongs to a vortical instability, as discussed by Chen et al. (2017). It benefits from Newly identified principle for aerodynamic heating in hypersonic flows Values are time-averaged along the normal direction. w θ 1 = −pθ (orange); w θ 2 = µ θ 2 (green); w ω = µω 2 (magenta). Re unit = 9.7 × 10 6 m −1 . SM, second mode; LFW, low-frequency waves. V1 and V2 stand for, respectively, the PIV view range x ∈ [157, 187] mm and [190,220] mm.
the increasing shear stress after HS. One interesting phenomenon, shown in both DNS and PIV, is that unlike w θ2 and w ω the peak width of w θ1 is much smaller than that of the second-mode amplitudes at HS. As stated earlier, this asynchrony reflects the fact the w θ 1-induced heating not only depends on the fluctuation amplitudes, but also the phase angle between them. The details of this mechanism will be discussed in future work.
As introduced above, one important phenomenon observed in prior studies is streamwise hot streaks, which were not only identified in the experiments (Berridge et al. 2010;Chou et al. 2011;Chynoweth et al. 2014) but also in the DNS (Sivasubramanian & Fasel 2015). However, such streaks are less clear in the present TSP results. One reason is that the noisy condition of the incoming stream makes the instantaneous spanwise flow structures rapidly vary in time. The effect of hot streaks is spanwise smoothed in the experimental data, as compared to DNS and PSE. To investigate this phenomenon accurately, we seed the most unstable spanwise wavenumber in the DNS so that those hot streaks are successfully repeated, as shown in figure 14(a). Figures 14(b), 14(c), and 14(d) show, respectively, the instantaneous y-direction-averaged w θ1 , w θ2 and w ω over the model surface. It can be seen that every streak is composed of second-mode waves, which is similar to the mechanism suggested by Lee & Wu (2008) (their figures 11 and 22) for incompressible flows. To verify the aerodynamic heating mechanism within these hot streaks, w θ1 , w θ2 and w ω are calculated at different circumferential positions, covering from one streak (φ = 19 • ) to its adjacent valley (φ = 10 • ), as shown in figure 15. No matter what the circumferential position is, the dilatation-induced w θ dominates the creation of aerodynamic heating. Especially along the streak, as shown in figure 15(d), w θ2 is larger than w ω even if it is underestimated. . DNS results of surface distribution of (a) mean temperature and instantaneous heat generation rates. Including (b) w ω = µω 2 , (c) w θ 1 = −pθ and (d) w θ 2 = µ θ 2 . Each quantity is normalized by its maximum value.

Heat transfer analysis
Theoretically, we may study the time-averaged surface-temperature distribution based on the steady-state energy equation (4.10) where ρ, C v , T, p and k are, respectively, the density, specific heat, temperature, pressure and thermal conductivity of the gas, and φ denotes the time-average of Φ.
For two-dimensional flow near the wall, the equation can be simplified to where the three terms represent, respectively, heat convection in x, conduction in y and heat sources caused by dissipation. In order to integrate (4.11), we linearize this Values are averaged along the time and the normal direction in the boundary layer. w θ 1 = −pθ (orange); w θ 2 = µ θ 2 (green); w ω = µω 2 (magenta). Re unit = 9.7 × 10 6 m −1 . equation by replacing u(y) with its mean value U, and assuming C v and k to be constants. The boundary conditions at an adiabatic wall and far from the wall are ∂T ∂y = 0, at y = 0; T = 0, at y = ∞. (4.12a,b) The initial temperature profile at x = 0 is assumed to be T 0 (y). The linearized problem can be solved using Laplace transformation, which yields the wall temperature T w = T(x, 0): (4.13) 174 Y. Zhu, C. Lee, X. Chen, J. Wu, S. Chen and M. Gad-el-Hak Neglecting the variation of the heat source in the y-direction, T 2 can be further simplified to (4.16) While T 1 represents the convection-conduction process without heat sources, T 2 stands for the integration of heat sources along the streamwise direction. The x-wise gradient of T w is given by that of T 1 and T 2 : (4.18) The sign of ∂T 1 /∂x depends on the initial distribution of T 0 (y). For a hypersonic boundary layer, dT 0 /dy < 0 and d 2 T 0 /dy 2 < 0, thus ∂T 1 /∂x < 0, which has the effect of cooling the wall, but its magnitude decreases with x. When φ reaches its maximum value, it begins to decay although T 2 still increases. Thus, if φ(x) decreases so rapidly that (4.19) at some downstream position, a peak value of T occurs there and the position slightly shifts downstream relative to the Φ-peak.
To verify (4.18), we compare the streamwise surface-temperature gradient and q = w θ1 + w θ2 + w ω for Re unit = 9.7 × 10 6 m −1 both from experiments and DNS, and second-mode instability amplitude for Re unit = 7.6 and 11.7 × 10 6 m −1 because no PIV data are available for the latter two. As shown in figure 16, the first peak positions of surface-temperature gradient agree well with that of q or second-mode amplitude, which justifies a posteriori the linearization of (4.11).

Conclusions
In this paper, a combined experimental, numerical and theoretical study of instability modes in a hypersonic boundary layer and their relevance to aerodynamic heating was conducted. The evolution of instabilities and wall temperature distributions were observed via physical and numerical experiments. The present study confirmed a previously observed heat-transfer peak that correlated with the second-mode-dominated transition process as well as revealing a new physical mechanism to promptly generate this heating impact. For the first time, our research identified a direct relation between the local temperature rise before transition was completed and the second-mode-induced dilatational heating.
By the joint use of TSP, PCB, Rayleigh-scattering and PIV techniques, the evolution of instability waves in a laminar hypersonic boundary layer and its relevance to aerodynamic heating were experimentally studied over a flared cone in the Peking University Mach 6 wind tunnel. Four unit Reynolds numbers, 5.4, 7.6, 9.7 and 11.7 × 10 6 m −1 , were investigated. PSE analysis and DNS were also conducted based on similar flow conditions. Several conclusions can be made as follows. with the second-mode (red) and low-frequency waves (blue) for (a) Re unit = 7.6 × 10 6 m −1 and (b) Re unit = 11.7 × 10 6 m −1 . Time-averaged heat generation rates (orange) for Re unit = 9.7 × 10 6 m −1 based on (c) PIV and (d) DNS. Each quantity is normalized by its maximum.
(i) Both the physical and numerical experiments showed that the evolution of a second mode brought about a local peak of surface-temperature HS that shifted slightly downstream. Increasing Re unit promoted the strength of HS but led to a stronger dissipation of the second mode. (ii) The heat generation rate induced by the dilatation and shear processes (denoted as, respectively, w θ and w ω ) were investigated. The second term included the pressure work w θ1 and dilatation viscous dissipation w θ2 . Different from the prevailing view for shear-induced heating, the aerodynamic heating in HS mainly arose from the high-frequency compression and expansion of fluid accompanying the second mode. The dilatation heating, especially w θ1 , was more than five times its shear counterpart. In some regions, even the underestimated w θ2 was larger than w ω . (iii) Based on a simplified two-dimensional, steady-state energy equation, heat transfer in the boundary layer and its effect on the surface temperature were analysed. Three factors were addressed, including streamwise heat convection, normal heat conduction and heat sources caused by viscous dissipation. The last factor heats the fluid passing through the source, resulting in a delay of surface-temperature 176 Y. Zhu, C. Lee, X. Chen, J. Wu, S. Chen and M. Gad-el-Hak rise in the streamwise direction. Such a mechanism was validated both by experiments and DNS, which explained why the peak of surface temperature shifted downstream compared to that of the second mode.
Two issues need to be further investigated in future work. First, a more accurate evaluation of the bulk viscosity coefficient in hypersonic boundary layers is necessary, which is an issue of crucial importance to the community. Second, since the phase angle between p and θ has played an important role in w θ1 , its detailed heating mechanism needs to be further clarified.