Direct numerical simulation of TS-waves over suction panel steps from manufacturing tolerances

To determine allowable tolerances between successive suction panels at hybrid laminar wings with suction surfaces, direct numerical simulations of Tollmien–Schlichting waves over different steps are carried out for realistic suction rates on a wind tunnel configuration. Simulations at given suction panel positions over forward and backward facing steps are carried out by the use of a high-order method for the direct simulation of Tollmien–Schlichting wave growth. Comparisons between high-fidelity direct numerical simulations and quick linear stability calculations have shown capabilities and limits of the well-validated linear stability theory design approach.


Introduction
The aerodynamic design of laminar profiles with suction is generally carried out by Linear Stability Theory (LST) as a quick and reliable prediction tool for the transition position [1]. Recently, alternatives such as direct numerical simulation (DNS) of transitional modes have shown promising results for the understanding of the physical background under special flow conditions such as geometrically singular cases, where physical fundamentals of classic LST are not valid. A typical geometric singularity of this class is the forward or backward facing step, violating the parallel flow assumption at corners that require a special singularity treatment, as demonstrated by Zahn, Edelmann and Rist [2,3]. For Tollmien-Schlichting modes (TS-modes), a two-dimensional approach of the perturbed boundary layer flow-field is not very time-consuming and allows the direct investigation of TS-modes over even complicated geometries, including non-linear growth effects. For the presented study, Crossflow instabilities (CF-modes) and their influences are not considered, to keep comparisons simple in the investigated region. This simplification by TS-mode treatment, including suction, can be found in Zahn [2] as well.
The assessment for integration of HLFC systems into a long-range aircraft has to answer the following two questions: how does this integration into a given long-range a/c configuration influence aircraft performance? What is the optimal HLFC aerodynamic and system configuration to obtain maximum performance benefit?
Principal feasibility of Hybrid Laminar Flow Control (HLFC) by suction sheets for large transport aircraft was shown for example by Airbus with the flight tests of an HLFC system on the vertical-tail plane of the A320 aircraft [5], where it was mainly used for the control of Attachment line and crossflow (CF) transition. The suction systems for these tests were designed to explore the limits of HLFC and were rather complex to determine the envelope of applicability of this suction technology. After having shown that HLFC fulfills the aerodynamic requirements, simpler and lighter systems have been developed to obtain the overall benefit for the aircraft.
A major step was the simplified suction system developed within the European ALTTA project [7,8]. This system works without complex internal structures of classical approaches [9,10] and was significantly refined within the German national project HIGHER-LE (Fig. 1). A numerical investigation of the suction on the boundary layer (BL) is carried out in [10] and [11]. Aside of CF-modes near the leading edge, under off-design conditions, TS-modes can be amplified as well. For this reason, all wing design studies at DLR take care of CF amplification as well as TS influences separately, as shown for example in [7].
Within different EU-projects, such as Clean Sky I and II, suction through a porous leading edge of a vertical-tail plane (VTP) is investigated by wind tunnel tests and numerical analysis [2,12]. As mentioned, the 2D TS-wave character allows simplified DNS of the flow field using a well-validated high-order numerical approach [13,14] which is used for a detailed simulation of transitional modes in the region of interest.
For experiments in the DNW-LLF wind tunnel, a verticaltail model with nearly no scaling and laminar leading-edge design will be investigated with respect to recently extended suction surface concepts. This so-called TSSD concept (Tailored Skin Single Duct) is a technical realization of different layers without the use of single suction chambers (Fig. 2). The outer skin is especially manufactured as a thin metallic layer with accurately etched suction holes, in contrast to former approaches of costly laser drilling [6]. For structural stabilization, different layers of metallic fabric are added underneath this thin surface, where this fabric is arranged by refining mesh diameters towards the suction surface. Especially under off-design conditions, TS-modes are no longer damped at the end of the suction surface. Future experiments may use even longer suction panels to control TS-transition as well, for this reason, 2D TS growth rates are considered at steps in the following.
The manufacturing of the etched top layer is limited by material dimensions of maximum 600 mm in streamwise direction and consequently needs welding with the neighbor metallic layer. This technically advanced connection process still generates small forward-facing steps on the surface of Fig. 1 Overview of suction panel geometry in flight experiments. The sketch of the A320 tail is taken from [4] approx. 40 microns (Fig. 3), influencing the stability of the BL by TS-modes.
Stabilization of the laminar flow by suction is assumed only in front and downstream of the welded region, where the connected foils prevent air flow through the surface. Consequently, an implementation of suction boundary conditions in the DNS-code is necessary for numerical modeling.
To demonstrate sufficient transition suppression by the suction panel with or without welding of the neighbor panel, direct numerical simulations of critical TS-modes are carried out here [14], including post-processing and comparison with growth rates on unperturbed surfaces and stabilization by uniform suction. The strong stabilizing effect of such suction techniques was already calculated by Rist and Zahn [2,15], using DNS and LST approaches.
Aside of inaccuracies due to the parallel flow assumption, non-linearities at the step may violate LST requirements and consequently additional empirical assumptions are necessary for prediction (see [16][17][18]). At the moment, very few physically reasonable corrections of this kind are available. In comparison with the LLF experiment at the fin, the investigated step at Re x ≈ 6 · 10 6 would appear at about 25% chord, which is in fact located downstream of the suction panel, but still in a realistic region for a HLFC design. In this region, the N-Factor is found by LST in a range between 3 and 4. Therefore, the investigation of amplified TS-modes over steps, generated by suction panel intersections, is of vital interest for calculations of resulting ∆N-Factors and the shift in transition positions. As shown in Fig. 4, the pressure gradient is nearly zero, for this reason a Blasius BL was chosen for the generic step geometry. Since the stability diagrams in Fig. 4 are generated by an incompressible LST approach, an overprediction of the amplification factors can be expected under compressible flow conditions (see [19] for in-flight experiments).  Since the free stream Mach number of the LLF experiment is only 0.295, these effects will have small impact. Nevertheless, local Mach numbers of M = 0.4 near the step may influence some LST results.
Within the critical area of the step, time-accurate simulations of the instabilities up-and downstream of the separated region are carried out, including grid refinement studies by DNS. For the extracted BL, unsteady perturbations were added at the inflow plane, which resemble TS-waves approximately in frequency and wall-normal perturbation profile. Determining growth rates of these TS-waves allows prediction of the accuracy of LST approaches and the prediction of the laminar-turbulent transition.
Local DNS on high-quality meshes of different resolution at different perturbation frequencies and amplitudes over millions of iterations were carried out. The process chain, reaching from the steady inflow BL profile to the implementation of temporal perturbations at the inflow plane, shows the capabilities of the presented approach.
For these calculations, a laminar incoming base flow is assumed for FFS (Forward Facing Step) as well as BFS (Backward Facing Step) geometry along a generic flat plate configuration with Blasius BL. For all cases, laminar separation appears upstream of the FFS, while for pronounced step heights another separation will appear downstream of the edge. The resulting pressure singularity at the step leads to additional amplifications of Tollmien-Schlichting waves behind the step, which will be shown by LST as well as DNS studies and comparisons of both techniques.

High-order DNS solver
All calculations in this paper are carried out with the DLR FLOWer-code [20] by solving the compressible full Navier-Stokes equations on block-structured two-dimensional grids, enabling treatment of complex aerodynamic configurations with any mesh topology. Dummy layers around each block are chosen to maintain the chosen accuracy in space at block intersections. While the standard version of this code provides second-order finite volume formulations for accurate and robust steady RANS simulations [20], unsteady direct simulations of BL-modes and their growth can be predicted by a high-order unsteady approach [13,21].
The basic equations and numerical schemes of the compressible high-order FLOWer-code will be described in the following. Generally, a fourth-order central differencing scheme based on a standard compact approach is used for the DNS. By this approach, two-and three-dimensional high-order simulations can be carried out for arbitrary geometries on block-structured grids.
High-order compact filters, which are applied at the end of each time step and sponge-zone boundary conditions, are optional to reduce reflections. For the present work, a sixth-order filter and the standard conservative form of the Euler terms is chosen. Time advancement is applied by a five-step second-order Runge-Kutta method or a standard fourth-order Runge-Kutta scheme [22].
The solved equations for a perfect gas with density , velocity components u i , pressure p and internal energy e , are written in conservation law form as Especially for a temporal TS-wave simulation or in cases, where the base flow must be fixed, forcing terms f i and g are included on the right hand side of the equations, such that a specified prescribed base flow ̄(y),ū i (y),Ē(y) is time independent for comparisons of DNS results with the temporal LST data.
In practice, these terms are evaluated numerically within the code by computing and storing the initial residual. For spatial simulations of TS-waves, no source terms are necessary; in these cases, f i and g are always zero since the BL will grow naturally. The equations are closed with the perfect gas law and the constitutive relations for q i and τ ij [13]. This approach is used without modifications as a basis for all calculations on different structured grids.

Linear stability theory
For the following investigations, the local linear approach is applied which is a subset of the non-local stability equations. The LILO-code [4,23], which is a spatial linear stability solver, is a development of Airbus and can be used for linear stability analysis as well as parabolized stability equations.
The stability equations are derived from the conservation of mass, momentum and energy, governing the flow of a viscous, compressible, ideal gas, formulated in primitive variables. All flow and material quantities are decomposed into a steady laminar base flow q and an unsteady disturbance flow q: The perturbation is represented by a harmonic wave: with the complex valued amplitude function q . Since LILO is a temporal code, frequencies are complex values, while wavenumbers and are real quantities. The growth rate is defined by i , which is the imaginary part of the complex frequency. For validations by TS-waves, the chosen DNS approach is inherently temporal, so a direct comparison between DNS and LST is possible. LILO is validated by several test cases against published results, including DNS, PSE (parabolized stability equations), multiple scales methods and LST.
The boundary layers for the LST calculations are generated by the BL-solver COCO, using pressure distributions from the CFD results. The use of original BL profiles from CFD is a possible, but not preferable approach. While results from a BL solver generate well-defined second derivatives of the near-wall velocity profiles, which is crucial for solutions of stability equations, CFD output usually generates nearwall oscillations of these second derivatives.
Near steps, the LST solver LILO uses profiles even in the separated flow, by omitting the step singularity itself. For this reason, a very small range around the step has to be omitted for the growth rate calculation. This is only useful for small surface imperfections, where the near-step range does not play an important role. To define sufficiently small imperfections, this study compares the described approach with complementary DNS output.

DNS step geometry and flow conditions
The investigated geometries are derived from a flat plate configuration with different forward and backward facing steps. Additional calculations without steps are carried out for comparison, especially between DNS and LST. Blasius boundary layers as generic velocity profiles are prescribed at the inflow, perturbed by TS-waves which develop in the growing BL. As well known, this Re is very large for a Blasius BL, but it was chosen from the wind tunnel model at near-zero pressure gradient.
The geometrical data of the step position and height as well as the flow conditions are taken from the experiments within the scope of the HIGHER-LE campaign, where a realistic laminar tail fin was investigated in the DNW-LLF [7,24]. In the following, results for the five chosen step geometries and three frequencies (see Table 1) are presented.
The flow conditions represent a laminar inflow BL from transonic wind tunnel experiments. For this reason, a local Mach number of M = 0.4 in front of the step and a displacement thickness Reynolds number Re 1 = 3735 is chosen. All references are taken from the displacement thickness of the inflow profile at x in = 0.5156 m, Re x,in = 4,715,380 (see Fig. 5).
Forward and backward facing steps of 57 and 114 μm are investigated, equivalent to 0.14 δ 1 and 0.28 δ 1 with a suction coefficient of C q = 1 ⋅ 10 −4 in front or downstream the step location at x step = 0.64636 m and Re x,step = 5,911,172 (Fig. 5). Grids of different densities near the step are generated to resolve TS-modes in the recirculation region. The resulting step Reynolds numbers are consequently Re h = 522 and Re h = 1042 , respectively.
Three TS-modes are at different frequencies and wavelengths as shown in Table 1, including wavenumber and phase velocity. All of them are strongly amplified at the step position for a Blasius BL, so they allow comparison of wavelength and frequency influences on growth rates at different steps. For these amplified modes, LST as well as DNS with two different step dimensions and suction panel positions are carried out.

CFD grids
Using the formally defined simplifications of the simulation region, a rectangular CFD grid is generated around the step position (Fig. 6). Refinement is carried out by considering strong velocity gradients at the suction panel in front of the step and the separation region.
The coarse grids contain 1440 × 240 cells while fine grids with 4320 × 320 are generated for refinement studies. The large step of 114 μm is resolved by 32 cells on the coarse grid in wall-normal direction.
The step dimensions of all FFS, generated by surface welding of suction sheets, result in step Reynolds numbers well in the non-critical regime, which is in good agreement with former empirical studies. For BFS, the larger step Reynolds number Re h is about 10% above the critical one defined by Nenni [27]. Since these criteria are well known as conservative estimates, we expect even the higher BFS as still sub-critical.
In the following, all coordinates are given in DNS grid scale, starting with x∕ 1 = 0 at inflow (x in = 0.5156 m, Re x,in = 4,715,380 from Fig. 5) and referenced by the

Unsteady inflow conditions
In accordance with the DNW-LLF experiments, the inflow Mach number is set to M ∞ = 0.4 in the far field at an inflow Reynolds number of Re 1,inflow = 3735 , where the inflow displacement thickness δ 1,inflow will be denoted as δ 1 in the following.
The steady inflow profile was generated by a Blasius solver, since wind tunnel conditions have shown good agreement with the BL profile of a flat plate near the investigated step position. Former detailed simulations with Falkner-Skan profiles could demonstrate the applicability of the Blasius approach for this investigation.
Approximations of TS-modes at small amplitudes A 0 of 1 ⋅ 10 −5 U ∞ in the wall-normal velocity perturbation maximum are added to the steady inflow velocity profile in each time step, using a harmonic amplitude variation at a given frequency. A 0 is the maximum v-perturbation of the amplitude function near 1.7 δ 1 wall-normal distance. The analytical Approximations in the v-component will eventually generate the amplified mode after some travel through the BL, which is filtering out all non-amplified parts. The small amplitude is chosen to stay safely in the linear regime of the perturbation, everywhere in the configuration's BL. From the prescribed perturbation at inflow, the resolved boundary layer upstream of the step will amplify the final TS-wave in this region. Its additional amplification by the step will be investigated further downstream in comparison with the smooth surface.
This is a similar approach as chosen for the linear stability solver, using velocity profiles from the BL solver COCO, adding perturbation functions locally, and solving the resulting perturbation equations. The resulting growth rates are consequently comparable in regions, where both approaches are valid. This comparison is carried out for different step geometries in the following validation studies.

Suction boundary condition
Relations between suction velocity and pressure difference in the suction chambers are usually taken from ALTTA studies (Fig. 7) where the following empirical formula was found  where the suction pressure is Δp sc , and w s is the suction velocity at the wall. To demonstrate these linear and quadratic influences of the above relation, the ALTTA diagram in Fig. 7 is added.
For the following studies, only a prescribed suction velocity, perpendicular to the wall is implemented by specification of the suction coefficient: For comparisons with wind tunnel data, this coefficient is taken from DNW-LLF studies by C q = 1 ⋅ 10 −4 .

Overview of the post-processing
The post-processing of TS-wave results from DNS needs to determine perturbation amplitudes from the respective modes. Visible temporal and spatial perturbations appear in the surface output of C f . Generally, C f is an appropriate quantity for the post-processing of transitional simulations and fits with shear-stress sensors in different experimental campaigns [20].
Suction is prescribed either in front or downstream of the welding region. To validate the quick local LST transition prediction at steps, comparisons with DNS data are carried out. Parts of the step region are omitted for the LST simulations to allow convergence of the LILO-code for all BL profiles. It is found that the LST with this simplified design approach over-predicts the influence of forward-facing steps as shown in Sect. 4.4.

Grid refinement studies
Simulations of TS growth with or without step at different grid densities are carried out for different steps with or without suction. For most cases, only coarse grid simulations of TS-amplification are chosen to find tendencies of step influences, especially for the small step Reynolds number. For all cases, where the small recirculation in the corner of the step geometry and especially the narrow recirculation on top of the step show pronounced influences, fine-grid simulations are mandatory (Fig. 8).
Nevertheless, the coarse grid results provide a good approximation of ∆N at the step for all cases. Figure 9 shows pressure distributions on the surface of all five investigated FFS, BFS, or flat plate geometries. The pressure rise in front of the forward-facing steps is well visible, followed by the discontinuity and another pressure rise  Fig. 9 Surface pressure distribution of all steps in comparison with a flat plate for FFS and BFS geometries on the coarse grid behind the step. These distributions are in good agreement with studies from Edelmann [3]. No such singular behavior appears on backward facing steps, but a maximum is reached behind the step. Finally, the surface pressure converges towards the flat plate pressure without gradients. These pressure gradients are important for the LST predictions of TS-wave growth and finally the transition line, since insufficient pressure distributions would generate inaccurate velocity profiles in the BL-solver COCO and finally wrong n-factors for the transition prediction.

DNS/LST comparison
For application of linear stability theory, at first, a surface pressure distribution is extracted from steady high-order calculations. For a step geometry, parts of the BL up-and downstream of the corner must be omitted (Fig. 10) to allow convergent LST results and solutions of the BL-solver COCO.
Since arbitrary suction distributions can be prescribed directly in COCO without any additional treatment, at first, comparisons of LST and DNS are carried out for smooth surfaces with suction but without steps. The suction region is placed in front of the position where steps are located, as well as behind this position, as shown in Fig. 5.
For comparison, three modes with different growth rates are added at inflow between 1850 and 2780 Hz (see Table 1). For these modes, the inflow amplitudes of DNS calculations are adjusted for each frequency to be comparable with the respective n-factor from the linear LST approach. This study is carried out without any step but includes the suction panels in the respective regions (Fig. 5). The amplitude growth from DNS and LST of these modes is in good agreement, as shown in Fig. 11 for suction up-or downstream of the step position.
Another study of DNS and LST growth rates is carried out on different steps without suction using the unstable TSmode at 2230 Hz. Comparisons of the linear n-factors from LST and the logarithmic wall shear-stress amplitudes are shown in Figs. 12 and 13 over FFS and BFS without suction. For these comparisons, only the slope and the relative deviations from the flat plate results are significant. The global offset between n-factor and logarithmic amplitude is just a result from the chosen inflow amplitudes, which is meaningless for transition calculation. The FFS simulations are carried out on a refined DNS grid, which was necessary to resolve the recirculation vortex properly.
In front of the steps, all predicted amplifications are in good agreement, as expected from the last comparison (Fig. 11).   Comparisons with existing empirical ∆N-correlations from Crouch [25] and Perraud [16] at small step Reynolds numbers are shown in Table 2 for different FFS and BFS. The empirical correlations in [25] and [16] were formally validated in comparison with experimental data by Methel in [26]. Comparisons were carried out for DNS results at 2230 Hz, since the other TS-modes show comparable n-factor shifting and consequently ∆N is expected in the same range as ∆n. Furthermore, the maximum ∆n behind the step and the nearly constant shift downstream at 600δ 1 is added to the table. The larger FFS would need more sophisticated empirical ∆N-factor calculations for ∆N Perraud [16], since Re h is above 700. For simplicity, this was not carried out, to get an approximate, but constant ∆N comparison. Generally, better agreement between Crouch and DNS-Max prediction is found.
The differences between LST and DNS results are still in the expected rage, which was already shown by various studies at different flow conditions. The approach of Rist and Zahn [2] for LST at steps has shown much better agreement at larger Re h , but requires the identification of nearwall modes behind the step and a mode-following technique, which is not implemented in the LILO stability solver for simplicity. The following comparisons between LST with the extended approach and DNS are intended as crosscheck for design predictions of the stability solver near steps. Deviations in N-Factors in first order will be not unexpected.
Small differences appear for the steps with h = 0.14 1 , especially for BFS. For FFS, where the necessary omitted region is larger, LST predicts earlier transition from the steps, while DNS predicts an influence of the separation in front of steps, but nearly none in the downstream area.
Substantial differences are visible for the higher steps of h = 0.28 1 . These differences could be reduced by implementing empirical ∆N-corrections into LST at the step singularity, which is still a difficult task under general assumptions. Further numerical enhancements near steps, as shown by Zahn and Rist in [2], will definitely improve the LST prediction near the step position, but the implementation for the utilized LST design tool LILO is still complex on general wing configurations.
Nevertheless, for the investigated steps and flow conditions, a successful cross-comparison of DNS/LST using quick simplified approaches near step-regions could be demonstrated for the welding region of interest. Tendencies of quick LST predictions may be extracted from these results; especially n-factors behind forward-facing steps are globally overpredicted by this very simplified approach and need additional confirmation by additional DNS output.
Nonetheless, at a step height of h = 0.14 1 , the deviation is not more than Δn ≈ 0.5 , which is generally acceptable for transition prediction. For BFS, even the larger height would be reasonable for this kind of LST design predictions, but in contrast to FFS, the transition location is overpredicted in comparison with DNS. Though substantial deviations between DNS and LST results are visible in Figs. 12 and 13, only the largest FFS will have a strong impact on the predicted transition position by LST.

Influence of suction on TS-mode growth
In the following, the influence of steps on the suction efficiency with respect to TS-mode damping will be examined. All steps are calculated with suction boundary condition in front or behind the step position respectively. For these test cases, the simplified approach of LST with extended BL calculation near the step will be discussed as well.
The amplitude growth of the C f perturbation is shown in Figs. 14 and 15 for the TS-mode at 2230 Hz in comparison with LST results.
For all test cases, the slope of the curves downstream of the steps is similar and the largest deviation of the n-factor, represented by the v amplitudes, is found for the BFS with the highest step Reynolds number.
In fact, the amplitudes behind steps with downstream suction are slightly smaller far behind the step position. An opposite behavior is visible directly downstream of the steps,  where suction in the upstream region reduces TS-mode amplitudes and consequently the start level of ln(A∕A 0 ). Comparison with LST in these diagrams show similar deviations for different steps as already discussed for Fig. 12 in the last section.
As a result, suction in front of steps can be recommended if the perturbation level is already large enough for nearly transitional behavior when transition is assumed directly behind the step. Otherwise, suction in the downstream region is slightly more efficient.

Conclusion
Within the presented study, direct simulation of TS-wave growth was carried out on generic steps with suction sheets at realistic step dimensions between 57 and 114 micron, realistic suction rates and suction panel size. For this configuration, grid refinement studies have shown good resolution of the separation regions around the steps. All step geometries are well in the range of non-critical step Reynolds number in agreement with former empirical studies.
Suction boundary conditions for distributed suction velocities were implemented in the high-order CFD-code for HLFC simulations. In addition, for the DNS with suction, grid refinement studies were carried out successfully for TS growth rates. Furthermore, non-reflective outflow boundary conditions are validated by long-term simulations, investigating the process of reflecting perturbations at outflow.
For the DNS calculations of TS-wave growth, unsteady mono-frequent velocity perturbations were added at inflow boundaries and a realistic behavior of TS instabilities was found behind the step using the high-order DNS solver. In agreement with LST predictions, amplified TS-waves were investigated downstream of the step. As a result, the development of growing unstable instability modes was demonstrated successfully downstream of the step region.
Within the numerical study, comparisons of direct TSwave simulations by a high-order code have shown reasonable agreement with LST results from BL data with step, especially for backward facing steps. In the cases with forward-facing steps, LST calculations show higher wave amplitudes compared to DNS results. In a first step, test cases without suction were carried out.
Considering these results, growth rates of TS-waves can be compared with different LST predictions including empirical extensions in the step region.
Further parameter studies have to be carried out by varying inflow perturbations of the boundary layer at different step geometries. Furthermore, the suction boundary condition of the high-order code needs ongoing testing.
Nevertheless, the applicability of this technique just for design purpose is successfully demonstrated at hand of the NACOR wind tunnel configuration with suction panels. As a result, further DNS of TS-transition for flow conditions of future wind tunnel configurations are applicable with reasonable computational effort.
Even validations of modern stability solvers for linearized Navier-Stokes equations, recently used for calculations of steps and gaps, are well in the scope of this approach.
As a conclusion, an improved prediction of TS-transition demonstrated at different surface imperfections with suction, is well in the scope of numerical DNS studies. Once again, the applicability of the chosen local DNS of TS instabilities as a supporting technique for LST design is well confirmed.
Funding Open Access funding enabled and organized by Projekt DEAL. Funding was provided by Deutsches Zentrum für Luft-und Raumfahrt (DE) and the EU CleanSky2 project NACOR.