Sensitivity Analysis of Transonic Flow over J-78 Wings
3D transonic flow over swept and unswept wings with an J-78 airfoil at spanwise sections is studied numerically at negative and vanishing angles of attack. Solutions of the unsteady Reynolds-averaged Navier-Stokes equations are obtained with a finite-volume solver on unstructured meshes. The numerical simulation shows that adverse Mach numbers, at which the lift coefficient is highly sensitive to small perturbations, are larger than those obtained earlier for 2D flow. Due to the larger Mach numbers, there is an onset of self-exciting oscillations of shock waves on the wings. The swept wing exhibits a higher sensitivity to variations of the Mach number than the unswept one.
In the 1960s through 1990s, numerical and experimental studies of transonic flow demonstrated formation of double supersonic regions on supercritical airfoils that comprise a flat or nearly flat arc [1, 2]. Such airfoils may arise in the course of a conventional profile redesign  or as solutions of an optimization problem . In particular, Jameson  designed four airfoils (J-75, J-78, K-80, and GAW-75-06-15) whose thorough studies demonstrated that they admit double supersonic regions and nonunique solutions of the Euler equations in narrow bands of the free-stream Mach number and angle of attack. Richter and Leyland  confirmed the nonuniqueness of transonic flow over the J-78 airfoil using an Euler solver with parallelization and adaptation of an unstructured mesh. Hafez and Guo  revealed multiple solutions for several symmetric flat-sided and wavy-sided airfoils. Caughey  analyzed the stability of inviscid transonic flow over J-78 and Hafez’s airfoils using time-dependent boundary conditions. A few more symmetric airfoils supporting nonunique transonic solutions were examined in . Meanwhile a nature of the nonuniqueness and the high sensitivity of the flow to variations of boundary conditions were not analyzed.
If a double supersonic region (S-region) forms on an airfoil, then it consists of a bow S-region, which expands as the free-stream Mach number increases, and a rear S-region which shrinks and either merges with the bow S-region or shifts downstream and vanishes. In 2002, Kuzmin  noticed that the coalescence of the bow and rear S-regions yields a pressure redistribution on the airfoil and a jump of the lift coefficient. Subsequent numerical studies of inviscid and turbulent flow confirmed the abrupt restructuring of 2D flow field due to the coalescence of S-regions at increasing or their rupture at decreasing . This phenomenon was scrutinized for a number of symmetric airfoils [10, 11] and also for symmetric ones, such as the J-78 airfoil, the Boeing 737 midspan profile , and a Whitcomb airfoil with aileron deployments . It was shown that both flow regimes (with the double and merged S-regions) can be either steady or exhibit self-exciting oscillations near the trailing edge due to the instability of the separated boundary layer. Figure 1 demonstrates a discontinuous dependence of the lift coefficient on the angle of attack α and , for example, for the steady 2D turbulent flow over the J-78 airfoil . The lift coefficient is calculated using the standard expression , where is the free-stream velocity module, is the free-stream density, is the wing area in planform, and is the normal force on the wing.
This paper addresses the 3D transonic flow over wings with J-78 airfoil at spanwise sections. The study is focused on the determination of Mach numbers and angles of attack that admit a high sensitivity of the flow field to small perturbations. Such adverse free-stream conditions may occur in practice at a vertical gust of wind in cruise flight of civil or transport aircraft.
2. Formulation of the Problem
We consider the fully turbulent flow past (semi) wings which extend spanwise from the root plane to the wing tip plane . Nondimensional coordinates , of the J-78 airfoil are given in . The trailing edge is normal to the root plane for both unswept and swept wings. In case of the swept wing, the leading edge makes an angle of 10 deg with the normal to the root plane, producing a taper ratio of 0.647346. The airfoil chord length in the root plane is set to 0.5 m.
The wing is inclined at an angle of attack α and enclosed in the box , ; see Figure 2. We prescribe a condition of symmetry on the root plane and the free-slip condition on the side, top, and bottom boundaries. On the inflow part of the boundary, , we impose the Mach number , the static temperature = 250 K, and condition that the velocity vector is parallel to the -axis. On the outflow boundary, , the static pressure = 50,000 N/m2 is given. The no-slip condition and vanishing heat flux are used on the wing, which is normally assumed to be smooth. In the end of Section 5 we mention also results obtained for a rough wing surface.
The air is treated as a perfect gas whose specific heat at constant pressure is 1004.4 J/(kg K) and the ratio of specific heats is 1.4. Then the temperature = 250 K determines the sound speed = 316.92 m/s. We adopt the value of 28.96 kg/kmol for the molar mass and use the Sutherland formula for the molecular dynamic viscosity. Initial data are parameters of the uniform free-stream, in which the turbulence level is set to 1%.
3. A Numerical Method
Solutions of the unsteady Reynolds-averaged Navier-Stokes equations were obtained with an ANSYS-15 CFX finite-volume solver of the second order accuracy, which is based on a high-resolution scheme for convective terms . An implicit backward Euler scheme is employed for the time-accurate computations. The code utilizes a linearization of the discretized equations and a multigrid accelerated factorization technique for solving the system of linear equations. We use a Shear Stress Transport k-ω turbulence model, which is known to reasonably predict aerodynamic flows with boundary layer separations .
Hybrid computational meshes were constituted by several millions prisms in 35 layers on the wing and over 107 tetrahedrons in the remaining region. The nondimensional thickness y+ of the first mesh layer on the wing was less than 1. Apart from the boundary layer region, mesh nodes were clustered in vicinities of shock waves and in the wake.
Test computations of transonic flow over the unswept J-78 wing on uniformly refined meshes of approximately 6 × 106, 14 × 106, and 20 × 106 cells showed that differences in the lift coefficient obtained on the second and third meshes do not exceed 5%. Global time steps of 2 × 10−5 s and 4 × 10−5 s yielded undistinguishable solutions. That is why we employed meshes of 14 × 106 cells and the time step of 2 × 10−5 s for the study of aerodynamic characteristics of the J-78 wings at various , α. The root-mean-square CFL number (over mesh cells) was about 2.
The solver was verified by computation of transonic flow for a commonly used test case of the ONERA M6 wing at , α = 3.06 deg. The Reynolds number based on the mean chord is 11.72 million. Figure 3 shows a good agreement of the surface pressure coefficients as computed with the ANSYS-15 CFX solver and a high-order discontinuous Galerkin solver by Bassi et al. .
Also verifications of the solver were performed by the simulation of an oscillatory transonic flow past a 18% thick circular-arc airfoil at , , and Re = 11 × 106. The obtained amplitude of lift coefficient oscillations was 0.35. This agrees well with the value of 0.37 obtained numerically in  using Spalart-Allmaras and Baldwin-Lomax turbulence models.
4. Unswept Wing: Lift Coefficient as a Function of at α = −0.6 deg
When the angle of attack is −0.6 deg, two-dimensional flow over the J-78 airfoil exhibits a lift coefficient jump at (see Figure 1). As mentioned in Introduction, the jump is caused by a coalescence of S-regions on the upper surface of the airfoil. Meanwhile, computations of the 3D flow over the unswept wing at = 0.84, = 50,000 N/m2, Re = 5.2 × 106 showed that S-regions reside at a considerable distance from each other. With the increase of from 0.84 to 0.8635 the distance gradually shrinks, while beneath the wing the flow velocity quickly rises. Figure 4 shows that the S-regions eventually coalesce on the upper surface at , and this only occurs in the root region.
Due to the larger Mach numbers, the terminating shock waves on the upper and lower surfaces of the wing are stronger than those in 2D flow. As a consequence, there is an onset of self-exciting oscillations of the shock waves and separated boundary layer at the rear of the wing. A dependence of the lift coefficient on time, for example, at , Re = 5.2 × 106, is depicted in Figure 5. The oscillations of are not periodic because of a superposition of different phases of the flow oscillations at different spanwise sections of the wing.
Figure 6 shows maxima and minima of the lift coefficient as functions of the free-stream Mach number. As seen, at Re = 10.4 × 106 the amplitude of lift coefficient oscillations is larger than that at Re = 5.2 × 106, though a jump of is nearly the same.
5. Swept Wing: Lift Coefficient versus at α = 0 and α = −0.3 deg
The leading-edge sweep angle of 10 deg provides a nearly uniform spacing between the bow and rear S-regions on the upper surface of the wing at various spanwise sections. Figure 7 displays the obtained sonic surfaces on the wing at = 0.845, α = 0. The height of the bow S-region is smaller than in Figure 4 due to the smaller Mach number, even though the angle of attack is larger than in Section 4.
It can be seen from Figure 7(b) that the height of the S-regions declines with increasing -coordinate. The decline sharpens near the wing tip where 3D effects are essential.
With an increase of from 0.845 to 0.849, the bow and rear S-regions approach each other and eventually coalesce in a vicinity of the wing root at . Contrary to expectations, the coalescence does not produce a jump of the lift coefficient (see the upper solid curves in Figure 8). A rise of is delayed to the Mach number at which the coalescence of S-regions accomplishes over a half of the wing span. Figure 9 shows plots of the pressure coefficient on the upper surface of the wing at and . We note that the difference of 0.003 in means a difference of only 0.95 m/s in the velocity . Though the lift coefficient oscillates in time, the plots of displayed in Figure 9 are stationary, as the oscillations are associated with static pressure fluctuations on the lower surface of the wing.
When the Mach number increases from 0.857 to 0.861, the lift coefficient declines due to the flow stabilization on the upper surface of the wing and an expansion of the supersonic region beneath the wing. The expansion yields a pressure drop on the lower surface and, hence, the drop of lift coefficient.
At a smaller angle of attack, α = −0.3 deg, a dependence of on is similar to the one discussed above for α = 0 though the rise of is shifted to higher values of (see the lower curves in Figure 8).
Dashed curves in Figure 8 demonstrate margins of the lift coefficient oscillations obtained at the vanishing angle of attack for a rough surface of the wing with a sand-grain roughness of 3 × 10−6 m. As seen, there is a good agreement of the results obtained for the smooth and rough surfaces.
For the swept and unswept wings with J-78 airfoil at spanwise sections and a span of two chord lengths, the numerical study revealed adverse bands of free-stream parameters , −0.6 deg , in which the lift coefficient exhibits a high sensitivity to small perturbations. In contrast to 2D flow over the airfoil, there are oscillations of the separated boundary layer and shock waves near the trailing edges of the wings. The sensitivity is expected to be even higher for larger wing spans.
Conflict of Interests
The author declares that there is no conflict of interests regarding the publication of this paper.
This research was performed using computational resources provided by the Computational Center of St. Petersburg State University (http://cc.spbu.ru/). The work was supported by the Russian Foundation for Basic Research under a Grant no. 13-08-00288.
R. T. Whitcomb and L. R. Clark, “An airfoil shape for efficient flight at supercritical Mach numbers,” NASA Technical Memorandum X-1109, Langley Research Center, Hampton, Va, USA, 1965.View at: Google Scholar
Ch. D. Harris, “NASA supercritical airfoils. A matrix of family-related airfoils,” NASA Technical Paper no. 2969, Langley Research Center, Hampton, Virginia, 1990.View at: Google Scholar
A. Jameson, “Airfoils admitting non-unique solutions of the Euler equations,” AIAA Paper 91-1625, American Institute of Aeronautics and Astronautics, 1991.View at: Google Scholar
R. Richter and P. Leyland, “Master-slave performance of unstructured flow solvers on the CRAY T3D,” in Parallel Computational Fluid Dynamics 1995: Implementation and Results Using Parallel Computers, pp. 297–306, Elsevier Science, 1996.View at: Google Scholar
M. Hafez and W. Guo, “Some anomalies of numerical simulation of shock waves—part II: effect of artificial and real viscosity,” Computers and Fluids, vol. 28, pp. 721–739, 1999.View at: Google Scholar
D. A. Caughey, “Stability of unsteady flow past airfoils exhibiting transonic non-uniqueness,” CFD Journal, vol. 13, pp. 427–438, 2004.View at: Google Scholar
A. Kuzmin and A. Ryabinin, “Transonic aerofoils admitting anomalous behaviour of lift coefficient,” Aeronautical Journal, vol. 118, no. 1202, pp. 425–433, 2014.View at: Google Scholar
F. Bassi, A. Crivellini, D. A. Di Pietro, and S. Rebay, “A high-order discontinuous Galerkin solver for 3D aerodynamic turbulent flows,” in Proceedings of the European Conference on ECCOMAS Computational Fluid Dynamics (ECCOMAS CFD '06), P. Wesseling, E. Oñate, and J. Périaux, Eds., Delft University of Technology, Delft, The Netherlands, 2006.View at: Google Scholar
T. J. Barth and D. C. Jespersen, “The design and application of upwind schemes on unstructured meshes,” AIAA Paper 89-0366, American Institute of Aeronautics and Astronautics, 1989.View at: Google Scholar
W. Geissler and L. P. Ruiz-Calavera, “Transition and turbulence modelling for dynamic stall and buffet,” in Proceedings of the 4th International Symposium on Engineering Turbulence Modelling and Measurements, Corsica, France, 1999.View at: Google Scholar