Grid resolution assessment method for hybrid RANS-LES in turbomachinery

ABSTRACT Hybrid RANS-LES is a promising method for analyzing the complex flow structure in turbomachinery due to its excellent compromise between accuracy and computational cost. However, it is challenging to employ the existing URANS and LES grid resolution evaluation approaches in a hybrid RANS-LES simulation, thus commonly the required grid resolution is unknown, increasing the risk that LES relies on overly coarse grids and leading to unreliable predictions. Hence, spurred by this deficiency and driven by the aspects of coupling grid resolution evaluation with the interest quantity during turbomachinery design and mechanism analysis, this work proposes a novel grid resolution evaluation method suitable for both LES and URANS. The suggested technique considers three grid resolution criteria: no effects on the time-averaged flow field, major unsteadiness, and minor unsteadiness. The developed method is verified employing a classic scenario, i.e. Pitz-Daily backward-facing step. Furthermore, we demonstrate our method’s applicability through a T106A turbine cascade example, confirming the feasibility and reliability of the proposed scheme. When applied to hybrid RANS-LES turbomachinery simulations, the research results highlight our method’s capability to reduce uncertainty and improve reliability during gridding.


Introduction
Hybrid RANS-LES is a promising approach for predicting the flow field in turbomachinery owing to its appealing computational accuracy and resources balance, especially for high Reynold's problems (Iim et al., 2012;Riéra et al., 2016;Ruiyu et al., 2019;Shi & Fu, 2013;Su et al., 2019;Yamada et al., 2017;Yan et al., 2019). For the nonzonal hybrid RANS-LES method, such as the detached eddy simulation (DES) (P. R. Spalart et al., 1997) and its improved variants, i.e. delayed detached eddy simulation (DDES) (P. R. Spalart et al., 2006) and zonal detached eddy simulation (ZDES) (Chauvet et al., 2007), the LES region highly depends on the local grid size (grid resolution) (Xia et al., 2018). While the grid resolution CONTACT Limin Gao gaolm@nwpu.edu.cn should be fine enough to capture the physical interest phenomenon accurately, it should also be as coarse as possible to reduce the computational complexity, which is much concerned from an engineering point of view. Thus, considering a suitable grid resolution is crucial for the hybrid RANS-LES approach, where suitable refers to minimizing the computational cost while still affording reliable predictions. Nevertheless, assessing grid resolution is quite challenging in a hybrid RANS-LES simulation due to the different grid resolution requirements in the URANS and LES simulations. Indeed, for a URANS simulation, managing a time-averaged flow field that does not change with grid refinement is sufficient, i.e. achieving grid independence. Accordingly, the grid resolution during LES should be able to resolve 80% of the turbulence energy. However, the region where URANS or LES is executed in the non-zonal hybrid RANS-LES approach is blurry and might change over time. Thus, evaluating if grid resolution meets the requirement in a hybrid simulation domain is still an open research case.
One way to solve the above problem in a hybrid RANS-LES simulation is to adopt the method utilized in URANS. The goal is to verify that the simulated time-averaged flow field presents a minor change as the grid density increases (commonly, the tolerance ranges from 1% to 5%). Benchmarking includes both global parameters (Jiang et al., 2021;Lei et al., 2017;Lizarose Samion et al., 2019;Shi & Fu, 2013;Tao et al., 2016;Zhao et al., 2017), i.e. pressure ratio and efficiency, and detailed flow field parameters (Barnes & Visbal, 2013;Sa et al., 2015;Xia et al., 2018;Xiao & Luo, 2015), i.e. vortex structures and time-averaged surface pressure distribution.
However, a sufficient grid resolution during a URANS simulation does not guarantee that LES is appropriately performed, as several flow-separation cases such as cylinders and triangular bluff bodies based on hybrid RANS-LES simulation, highlight that although the gridscale converges for time-averaged results, the impact on dynamics characterized by root-mean-square (RMS) is still remarkable (Dong et al., 2018;West et al., 2017;Xiao & Luo, 2015). In the context of industrial studies, an alternative solution for the above problem is performing simulation on the grid, which is designed using the finest resolution possible given the available computing resources (Gant, 2010). However, this pragmatic approach does not ensure that the maximum grid resolution can resolve 80% of the turbulent energy. Hence, it is unknown if the grid resolution sufficiently captures the dynamic physical quantities we are interested in.
To solve the difficulty of evaluating grid resolution in the hybrid RANS-LES simulation and improve the reliability of the dynamic turbomachinery results, this work presents a novel grid resolution evaluation technique suitable for both LES and URANS. Our method considers flow unsteadiness and is also driven by coupling grid resolution evaluation with the interest quantity during the turbomachinery design and mechanism analysis. In the following sections of the paper, the principle and verification of the grid resolution evaluation technique will be presented, and an example to demonstrate the method's applicability in turbomachinery will be given.

Grid resolution criteria
The present paper proposes three levels of grid resolution criteria based on different demands: • Criterion 1: No effects on the time-averaged flow field.
When the grid-scale is reduced according to a proportion, the tolerance of the time-averaged results (referred to as Res space ) should be within e s . The Res space is defined as Res space = |(ḡ 1 −ḡ 2 )/ḡ 2 |, whereḡ 1 andḡ 2 represent the time-averaged physical quantity of the flow field before and after grid refinement. Moreover, the choice of e s depends on your accuracy requirement, usually between 1% and 5%. The quantitative description of the grid-scale is discussed in section 2.1.4. • Criterion 2: No effects on major unsteadiness. Based on satisfying Criterion 1, when the grid-scale is reduced according to a proportion, the unsteady frequency tolerance (referred to as Res time ) of major fluctuation modes between two grid cases should be within a specific error range, i.e. (Res time ) maj ≤ e t . The 'major fluctuation modes' means the modes where the fluctuation energy ratios are relatively higher than others after the POD decomposition. The definition of (Res time ) maj and the value of e t is discussed below. • Criterion 3: No effects on minor unsteadiness. Upon satisfying Criterion 1 and 2, when the grid-scale is reduced according to a proportion, the tolerance of the minor unsteady frequency (referred to as (Res time ) min ) should be within e t . Similarly, the 'minor unsteady frequency' corresponds to the frequency in the POD modes containing lower fluctuation energy, ensuring the grid scale has negligible effect on detailed flow structure capturing.
Choosing the evaluation criterion should be made according to the practical requirements. Hence, if accurate time-averaged results are required, 'Criterion 1' is sufficient. If main unsteadiness is considered, 'Criterion 2' is more appropriate, while if detailed flow structures are considered, 'Criterion 3' is needed. It can be speculated that the grid number gradually increases from satisfying Criterion 1 to Criterion 3, with commonly, Criterion 2 is sufficient. A detailed description of these criteria is discussed below.

Indicator choice
The grid resolution is assessed utilizing the Res space or Res time indicator, coupled with the physical interest quality in the turbomachinery design and mechanism analysis.
In Criterion 1, the time-averaged physical interest quantity is selected to construct the indicator, involving overall performance parameters, i.e. pressure ratio, efficiency and stability margin, or local aerodynamic parameters, i.e. static pressure distribution, total pressure loss coefficient, and wall shear stress. This criterion considers the same grid-independence analysis as in the RANS simulation.
Criteria 2 and 3 describe the grid resolution effect on the dynamic flow field, with indicator Res time employing the unsteady frequency of the aerodynamic parameters Rotating stall (Gao et al., 2015) Frequencies of rotating stall cell.
Frequencies related to vortex separation or shedding. Vortex induced blade vibration (Lau et al., 2007) The frequencies of pressure fluctuation on the blade surface Frequencies not directly related to vibration Flow-induced noise (Sharma et al., 2019) Dynamic pressure fluctuation frequency with high power spectrum density (dB) Dynamic pressure fluctuation frequency with low power spectrum density (dB) Flow loss (Ruiyu et al., 2019) Large-scale fluctuation frequency in high loss region. Such as frequencies of tip leakage vortex in tip clearance region and frequencies of passage vortex, corner vortex in corner separation region.
Fluctuation frequencies in low loss region Unsteady flow control (Li et al., 2010) Frequencies of the vortices induced by the flow control technology and the frequencies of the vortices mixed with it Frequencies not directly related to unsteady flow control Rotor-stator interference (Gourdain, 2015) Frequencies of wakes, Frequencies of large eddy frequency close to tip region.
Fluctuation frequencies related to the vortexes inside the wakes such as pressure, velocity, or temperature fluctuations. The reasoning for that choice is that unsteady frequency is one of the most crucial parameters in the highperformance turbomachinery design or analysis of physical phenomena (Ez Abadi et al., 2020;Gao et al., 2015;Ghalandari et al., 2019). Indeed, the unsteady frequency always acts as a surge warning signal of a compressor rotating stall (Mao & Liu, 2016), as an essential indicator in an unsteady flow control method (Werder et al., 2020), it plays a vital role in improving the fatigue life of a vortexinduced vibration blade, and reduces the flow-induced noise.
The main frequency, obtained by fast Fourier transform (FFT) of the unsteady physical parameters (velocity, pressure, or temperature), refers to the dynamic characteristics that play a crucial role in the research problem. Compared with the minor frequency, it is often characterized by higher fluctuation energy and a larger vortex scale, which is why we can extract it through POD analysis. Table 1 presents a detailed description of the major and minor frequencies in typical turbomachinery problems.

The definition of Res time as well as the value of the e t :
Assume that the normalized frequency in the k-th mode based on low grid density is f k and is f k std based on larger grid density. Then, the Res time is defined as: The error range e t is defined as: where N min is the minimal number of samples in a period and t physical denotes the physical time step. The paper takes the theoretical minimum of periodic fluctuation captured for the sin type, i.e. N min = 5. The meaning of the definition is that the difference in physical time steps between flow fields of two cases cannot capture a new sin-type fluctuation.

Proper orthogonal decomposition
The proper orthogonal decomposition (POD) method is used to analyze the influence of the grid density on the unsteadiness of the flow field (Guan et al., 2021;Shen et al., 2020). It decomposes the flow field into different modes according to unsteady fluctuation energy to provide the most powerful way to capture the dominant unsteadiness in the flow field. This section briefly introduces the POD method.
The POD method originated in the statistical analysis of vector data and was first proposed by Lumley (Lumley, 1967). The input data for the POD consist of N snapshots sampled in time, and the flow field at each snapshot is described by K spatial grid cells. Considering that timeaveraged flow fields are compared separately, the POD method is only applied for the fluctuation term. The fluctuation flow field at the n-th moment is defined as: Thus, the input matrix is constructed as G K N = {g 1 , g 2 , · · · , g n , · · · , g N }.
The basic idea is to decompose the unsteady flow field (G K N ) into an algebraic sum of a time coefficient function related only to time (B N N ) and a standard basis function (referred to as 'mode') related only to space ( K N ), as given in Eq. (4).
The key to the POD method is to find a set of standard orthogonal basis functions K N = {ψ 1 , ψ 2 , · · · , ψ N } that satisfy: (5) where operator | · | denotes the norm and (·, ·) is the inner product. I represent a unit matrix.
The correlation matrix C N N , signifying the time correlation of the flow field at any two moments, is constructed as: Then, solving the eigenvalue problem: where λ n represents n-th eigenvalue. A N N is the eigenvector matrix. The eigenvalues and eigenvectors are arranged in descending order of the eigenvalues, which is λ 1 ≥ λ 2 ≥ · · · ≥ λ N . The orthogonal basis function The standard orthogonal basis function is given by Eq. (9) and ψ n / √ λ n represents the n-th mode to satisfy the demand in Eq. (5).
The physical meanings of the POD analysis are listed below: (1) Fluctuation energy The eigenvalues characterize the fluctuation energy of modes; thus, the energy ratio of the i-th mode to the total energy is: (2) Unsteady characteristic of fluctuation The time coefficient, representing the unsteady characteristics of each mode, is defined by: Each column of B N N symbolizes the evolution of the fluctuation over time, and we can extract the frequency in each mode by FFT analysis of the corresponding mode. Each row of B N N exemplifies the contribution of each mode in the flow field at the moment.
A detailed description of the POD method can be found in Chen et al., 2012.

Quantitative description of the grid-scale
The grid-scale quantitative description intuitively reflects the grid resolution size and improves grid generation portability, i.e. once the grid-scale is determined, it can be employed as a guideline while generating a grid of similar geometry without repeating the grid resolution evaluation when the same numerical solver is used. This is a crucial issue to reduce computational cost.
To quantificationally describe the grid-scale in a certain region, we define a parameter named the averagenormalized grid-scale ( l + ), which represents the average of the smallest turbulent length scales that can be resolved when the LES model is applied. The length scale l definition relies on the filter width employed in the hybrid RANS-LES method. Current literature suggests three commonly used definitions where the length scale of grid cell k can be expressed as: a local maximum cell size (Eq.(12)) (P. R. Spalart et al., 1997), a local cell volume's cubic root (Eq.(13)) (Breuer et al., 2003), and a vorticity-based length scale (Eq. (14)) (Riéra et al., 2016).
The normalization of the grid scale l at cell k is defined as: Here, C denotes the characteristic length, which is the chord length in the turbomachinery cascade; Re τ is the Reynolds number and is defined by the wall shear stress. Local Re τ , k (grid cell k) is given as: where μ l,k and τ w, k represent the laminar turbulence coefficient and the magnitude of the wall shear stress at grid cell k, respectively. Re τ is defined as the maximum of Re τ , k . Considering that the grid cells close to the leading edge and trailing edge are generally refined, the maximum value is selected after removing the values within the range of 5%C near the two edges, that is: Since the grid is nonuniform, giving the locally normalized grid-scale in the entire computational domain is not convenient for summary and analysis. Therefore, we defined the volume-averaged l + of a specific region as the criterion parameter, given in Eq. (18). The determination of the region is related to the gridding strategy, which is described next.

Gridding partition strategy of turbomachinery
The grid resolution requirements in a turbomachinery flow field vary with the adverse pressure gradient magnitude and the flow complexity extent. Therefore, it is necessary to partition the entire computational domain and evaluate the grid resolution separately. In 2001, Spalart (R. Spalart, 2001) provided grid guidelines for external flow field simulations using the DES method, such as aerofoils. Similarly, according to the flow characteristics of multiple walls, complex physical phenomenon like shockwave and vortex in the turbomachinery flow field, the computational domain can be divided into four parts: the focus region (FR), RANS-LES region (RL), inlet region (IR), and RANS region (RR). Take a three-dimensional (3D) compressor rotor blade, ROTOR37, as an example, as presented in Figure 1. The characteristics of each region are listed below: (1) FR region: The most critical region is colored red in Figure 1, which contains significant complex flow phenomena in the region, such as separation, shock wave, and vortex interaction. In the turbomachinery flow field, complexity always appears near the blade tip, and hub corner due to tip leakage flow, corner separation, wake, and they have a significant effect on its performance, which deserves detailed analysis. The turbulence in the region is characterized by anisotropy and nonequilibrium. Thus, it cannot be adequately described by the RANS or URANS method with the eddy viscosity turbulence model. The grid-scale in this region is relatively small compared with others.
(2) RL region: The region is located in the blade passage without the complex flow phenomena mentioned above and is colored yellow in Figure 1. Controlling the airflow through the blade passage is the basic principle of turbomachinery; thus, the passage flow deserves our attention. While the flow fields in some regions are relatively simple, such as the region near the mid-span or with an accelerating pressure gradient, we define it as the RL region. (3) IR region: The region contains inlet extension from the inlet to the position near the leading edge and is colored blue in Figure 1. The flow field in the region is relatively uniform and stable. The grid number in this region is relatively small due to the low-pressure gradient and unresolved inlet flow unsteadiness. (4) RR region: The region close to the solid walls and outlet is colored green in Figure 1. The RANS method is expected to be employed in the region. The grid number in the boundary layer close to the wall is approximately 20. The region one chord away from the trailing edge is also defined as the RR region because the intensity of flow mixing has been reduced, and the large-scale vortex has been considerably dissipated.

Grid resolution evaluation process
The flowchart of grid resolution evaluation process for hybrid RANS-LES is displayed in Figure 2. First, the demand for accuracy should be clarified. 'Demand 1' means the time-averaged flow field does not change with grid refinement; 'Demand 2' means based on 'Demand 1', the major unsteadiness caused by large-scale fluctuation does not change with grid refinement; similarly, 'Demand 3' means the grid scale has little effect on not only the major unsteadiness but also the minor unsteadiness. Second, we separate the entire calculation domain into four regions based on the Gridding partition strategy. Third, the grid in different regions, especially the FR and RL regions, should be refined with a particular proportion, and several cases with different l + are obtained. After DES simulation, we can analyze the time-averaged results and fluctuation flow field separately. For 'Demand 1', satisfying 'Criterion 1' is sufficient. For 'Demand 2', Criterion 2 should be achieved; otherwise, we need to continue to refine the grid. Similarly, only by achieving Criterion 3 can we meet Demand 3.
As a result, l + can be quantified to meet different demands, which can be used as a basis for gridding for simulations with similar geometries.

Destination of verification
In turbomachinery, the unsteadiness can be decomposed in three parts. Take pressure as an example: p = p + p + p . The term p represents the time-averaged part, the term p represents the periodic part of the unsteadiness (e.g. rotor wakes) and the term p is related to turbulent fluctuations (with p = 0) (Gourdain, 2015).
According to the grid resolution criteria definitions of Section 2.1, the influence of the grid resolution on a timeaveraged flow field and the periodic frequency (the major unsteadiness) can be evaluated. Therefore, the grid resolution assessment method is suitable for URANS, as only these two terms can be resolved in this grid resolution evaluation approach.
Given that most of p is resolved in LES, except for a small part provided by a sub-grid model, verification shall focus on whether the grid size that meets the grid resolution criteria proposed in this paper is sufficient for LES simulation when predicting large scale flow separation.

Verification strategy
We consider the classic Pitz Daily turbulence case in the OpenFOAM tutorial for the verification strategy presented in this sub-section. The corresponding model geometry is illustrated in Figure 3, comprising a short inlet, a backward-facing step, and a converging nozzle as the outlet. The simulation is applied on a two-dimensional grid, involving a single 1 mm length cell in the span direction. Considering that we investigate whether this paper's grid resolution evaluation method is suitable for LES simulation, rather than flow mechanism analysis, we employ the two-dimensional geometry model used in LES, despite the results having no physical meaning.
A well-resolved LES should aim to resolve 80% of the turbulence energy, i.e. IQ = k res /k total ≥ 0.8. The total turbulent kinetic energy is considered to comprise the sum of the modelled (k mod ), resolved (k res ), and numerical (k num ) components, where commonly, the numerical component is neglected. The resolved component can be calculated by k res = 0.5(u 1 u 1 + u 2 u 2 + u 3 u 3 ), while the modelled component (k mod ) is obtained through the subgrid turbulence model. Ten cases are obtained by altering the grid size, with the grid's quantitative description for these cases presented in Table 2. We focus on whether the grid resolution that satisfies the suggested criteria 2 and 3 satisfies IQ ≥ 0.8.
The LES simulations are performed on the Open-FOAM platform utilizing pisoFOAM solver, and a second-order backward time integration as well as a second-order central differencing in space. The dynamic k-equation is set as the subgrid-scale model, and the physical time step is 2 × 10 −5 s. Table 3 describes the boundary conditions.

Verification results discuss
The corresponding IQ results for the ten cases examined are illustrated in Figure 4. The latter figure highlights that case g1 to g5 the IQ ≥ 0.8 is satisfied in the large  separation area, indicating that the grid resolution in these cases meets the LES requirement. In contrast, case-g6 to g10 presents an insufficient grid resolution for LES due to regions with IQ < 0.8. Specifically, for cases g6 to g8, the poor grid resolution appears mainly close to the outlet, while for cases g9 and g10, the grid resolution is poor both in the region close to the outlet and where separation is generated, i.e. the region after the steps. The grid resolution of the ten cases is evaluated based on the grid resolution criteria proposed in this paper. The stream-and span-wise fluctuating velocity field is decomposed by POD, containing 500 snapshots, where ten-time steps separate adjacent snapshots.
The energy percentage of each mode after POD decomposition is illustrated in Figure 5. The abscissa is the mode number, and the ordinate is the proportion of fluctuation energy per mode, calculated according to Eq.(10). Figure 5 highlights that the energy ratio of the first two modes is significantly higher than the remaining ones. Therefore, modes 1 and 2 are regarded as the major fluctuation modes.
Considering mode 1 as an example, Figure 6 presents the vorticity magnitude based on the fluctuation velocity. For cases g1 ∼ g5 (Figure 6(a-e)), the vorticity contours are similar, showing the flow separation after the step, the flow separation on the upper-end wall, and the mixture of these two. However, for cases g6 ∼ g10, the vorticity contours distribution is different, as the flow separation generated at the upper-end wall is different from dissipating compared with the first five cases. This is because the sub-grid model cannot correctly dissipate large-scale vortices caused by an insufficient grid resolution.
The mode difference is also manifested through frequency. Figure 7 shows the relationship between Δl + and the frequencies obtained after POD decomposition for all cases examined in this work. In this verification case, the frequencies contained in the main modes, i.e. modes 1 and 2, are not unique. Therefore, to facilitate analysis, we extract the frequency corresponding to the highest amplitude, and considering the minor modes, we extract the frequency corresponding to the highest amplitude other than the major one. Diamonds and circles in Figure 7 are the major and minor frequencies of the above ten cases, respectively, while the cases within the yellow area satisfy IQ ≥ 0.8. Figure 7 indicates that when criterion 3 is met, i.e. the major and minor fluctuation frequencies are both within the error range, the grid resolution meets the LES requirement. When criterion 2 is met, the reliability of the major fluctuation frequencies can be guaranteed.   Although this criterion cannot guarantee that the grid resolution of all LES simulation domains meets the requirement, the grid size in the major area where the separation is generated is guaranteed to meet the requirements, which is also the area we are interested in. Hence, in conclusion, the grid resolution assessment method proposed in this paper is suitable for URANS and LES simulation.

Application in a turbine cascade
In this section, we take the turbine cascade T106A as an example to demonstrate a detailed calculation process for the criterion and quantitatively present the grid generation guidelines. The blade-to-blade (B-B) section is the basic unit of turbomachinery, and it represents the principle of energy transformation. Thus, here we focus on the B-B section.

T106A cascade
T106A is a two-dimensional turbine cascade with a relatively low Reynolds number. It is a classic standard example of computational fluid dynamics that does not consider the problem of heat exchange and cooling. The geometry is shown in Figure 8. The geometric and aerodynamic parameters are listed in Table 4 and Table 5, respectively. The relevant data can be found in Refs. (Luo et al., 2019).  According to Refs. (Garai et al., 2015;Stieger & Hodson, 2002), the spanwise extent of the computational domain was set to 15% of the chord length to ensure that vortices can develop in the spanwise direction. The computational domain was extended by 1.0 upstream and by 2.0 axial chords (C ax ) downstream of the leading and trailing edges, respectively.
The original blade was defined by 162 points (Luo et al., 2019), and it was parameterized and redefined by a larger number of points (563) to avoid geometry discontinuity.

Gridding strategy
Based on the gridding partition strategy promoted in section 2.2, the calculational domain is separated into IR, RL, RR, and FR regions. For the flow field in T106A, complex separation occurs near the region close to the trailing edge and the end of the suction surface. We regard it as the FR region, which is exhibited in Figure 9 and colored red. Although multiwall surround the blade passage, the flow is relatively stable, and no separation occurs, which can be regarded as the RL region, as shown in Figure 9 and colored yellow. The IR region contains inlet extension, and the RR region contains the region close to the wall in addition to the region downstream of the outlet extension.
For the RR region, the grid-scale only needs to meet the requirements of RANS, which means: y + wall ≤ 1, maximum grid aspect ratio < 1000, and maximum grid extension ratio 5. Here, we only focus on the grid density in the IR, RL, and FR regions. The average-normalized  Change in the grid-scale in FR region Change in the grid-scale in RL region Change in the grid-scale in IR region l + (F 5 )| IR < l + (I 2 )| IR < l + (I 1 )| IR grid-scale ( l + ), defined in section 2.1.4, is used to describe the grid-scale quantitatively. It can be calculated according to Eq. (15). The local Re τ , k is computed using the RANS method and according to Eq. (16). The distribution around the blade is illustrated in Figure 10. We can obtain that Re τ = max (Re τ , k ) = 4400 according to Eq. (17). Eleven cases are set by separately changing the l + in those three regions, as described in Table 6. The cases in which the grid scales of the FR region changed are labeled as F-series (including F1, F2, F3, F4, and F5). Similarly, the cases where grid scales of RL and IR regions transformed are labeled as R-series (including R1, R2, and R3) and I-series (including I1 and I2), respectively. Case B is a basic grid; the grid number in the FR and RL regions is based on the LES case in Michelassi et al., 2002. Since the inlet and outlet extension were longer than those indicated in the literature, the number of grids was correspondingly higher.
The average-normalized grid scales ( l + ) of different regions are presented in Figure 11. In each series, the grid-scale of only one region (FR, RL, or IR) is changed. The decreasing proportion of l + in the Fseries is approximately 1.3 and is approximately 1.5 in the R-series. In both the F-series and R-series, the l + of the other two unchanged regions are the same as that of case B. Similarly, the I-series (including I1 and I2) are obtained by changing the l + of the IR region based on case F5, and the grid-scale reducing proportion is approximately 2. The grid number of each case is given in Table 7. The total grid cell number of B was 4 × 10 6 .

Calculation setting
The DES method, most sensitive to the grid, is used for numerical calculation. The flow solver EURANUS, integrated into the FINE/Turbo TM software package, was used in the present study. A finite volume method with the second-order central difference scheme was employed in this solver for space discretization. In the DES simulation, the second-order and fourth-order artificial viscosity coefficients were set to 0.2 and 0.05, respectively, to reduce the influence of numerical viscosity. Meanwhile, in the RANS simulation, these parameters were set to 1.0 and 0.1 to ensure stability. The dual time step method is applied for time marching and the physical time step t physical is set as l 0 /U max , where U max is typically more than 1.5 times the inlet velocity, and l 0 is the averaged grid-scale, which is 2 × 10 -6 s in the current work. The maximum inner literation is 40. The pitch-wise boundary condition is set to periodic. The boundary condition at span-wise is set as symmetric. All walls are set as no-slip adiabatic. The inlet total pressure, total temperature and flow incidence are imposed at the inlet, whereas the static pressure is specified at the outlet according to the experiment results given in Table 5. The outlet Mach number is 0.404 and the Reynolds number based on axial chord is 51800.
For the unsteady convergence criterion, we set pressure monitors near the blade trailing edge. Convergence was achieved if the main frequency was the same over 500 consequent iterations. The primary frequency was obtained using the fast Fourier transform algorithm (FFT).

Demand 1
'Criterion 1' should be satisfied to meet 'Demand 1' discussed in section 2.1, meaning the changes in the timeaveraged flow field should be within e s when refining the grid. We set e s as 1% here, which is the same as the grid-independent analysis in the URANS simulation. The time-averaged static pressure coefficient C p of the blade surface predicted for different cases is compared with the corresponding experimental values, as shown in For the FR region (Figure 12(a)), case F1 fails to satisfy 'Criterion 1' because the difference between F1 and F2, i.e. Res space (F1, F2), is 3.3% which exceeds the e s . All the cases except F1, R1 and I1 satisfy 'Criterion 1'. Thus, l + should be at least 70 in the FR region, which is the value of l + in case F2. Similarly, l + should be at least 59 in the RL region and 210 in the IR region.

Demand 2
'Criterion 2' should be satisfied to meet 'Demand 2' discussed in section 3.4, meaning the frequency changes of the main flow unsteadiness should be within e t (see Eq. (2)) when refining the grid. The POD method was applied to extract the main flow unsteadiness, strongly affecting the entire flow field. Five hundred instantaneous snapshots available from the DES approach were adopted for computing the POD. The time between adjacent snapshots, i.e. t physical was 2 × 10 -6 s. Due to the large number of grids and time steps, the computational resources used for solving the entire domain using the POD method are not feasible. Therefore, a specific region of interest (ROI) was selected to reduce resource consumption (Lei et al., 2017). The selected ROI is presented in Figure 13 in blue. All the areas with pronounced unsteady fluctuations were already included.
Given that F1, R1, and I1 did not meet Criterion 1, they were not considered in the unsteadiness analysis. Figure 14 demonstrates the energy of the first ten modes after the POD decomposition of eight test cases. Although the energy ratio of each mode to total energy varies, the ratio of the first two modes is significantly higher than the other modes. The energy ratio of the first two modes E 1+2 (E 1+2 = E 1 + E 2 ), shown in Figure 15 as a line graph, accounts for more than 60% of the total energy. That is, the first two modes are the major fluctuation modes.
We extract the unsteady frequency in the first two modes and present it in Figure 15 as a bar chart. We can see that only case F2 is different from the others. According to Eq (1), the frequency difference comparing cases F2 and F3 is: which is larger than the e t , i.e. N min t physical = 10 −5 . All the cases mentioned in Figure 11 except case F2 satisfy 'Criterion 2'. Therefore, considering the major unsteadiness, the ranges of l + that meet 'Demand 2' are summarized as follows: l + | FR ≤ 54, l + | RL ≤ 59, and l + | IR ≤ 210.

Demand 3
'Criterion 3' should be satisfied to meet 'Demand 3' discussed in section 2.3, meaning the frequency changes of main flow unsteadiness and minor flow unsteadiness should be within e t when refining the grid.
The key to this analysis is to find secondary modes that contain minor unsteadiness. From Figure 14, we can see that the mode energies of Mode 3 to Mode 5 are very close to each other. Although the energies of these three are  significantly smaller than those of Mode 1 and 2, they are relatively higher than other modal energies beyond Mode 5. Therefore, they can be regarded as modes containing minor unsteadiness, and we call them secondary modes. The frequency of a coupled flow field of Mode 3 + Mode 4 + Mode 5 is analyzed in this section.
(1), Res time (F3, F4) is 1.3 × 10 −4 , which is greater than e t (10 −5 ); that is, case F3 does not meet the demand. Thus, l + for the FR region should be no more than 44; (2) For the R-series ( l + | R2 = 59, l + | R3 = 40), the high-frequency differences between R2 and R3, as well as R3 and B, are less than the e t . For low frequency, case R2 shows a large difference compared with B. The Res time of that is 1.3 × 10 4 , which means the R2 does not meet the demand. Thus, the l + for the RL region should no more than 40; (3) For the I-series ( l + | B = 100, l + | I1 = 210), the high frequency can satisfy the demand, while none can capture the low frequency; thus, l + in the IR region should no more than 100.
Based on the analysis above, we can quantitatively summarize the average-normalized grid scale l + that meets different demands, given in Table 8. DES simulations with similar flow field structures can use these guidelines, and commonly Demand 2 is sufficient. Compared with the basic grid case B, in which the grid density in the RL and FR regions satisfies LES demand, the grid number referring to Demand 2 can be lowered by 53%.

Grid resolution analysis of multiple regions
The above grid generation bases were proposed based on varying the grid-scale of a single region. To analysis  whether the grid scale that meets specific requirements obtained based on a single region can also meet the above requirements when multiple regional grid scales change simultaneously. The grid cases in which the grid scales of the IR, RL, and FR regions satisfy the Demand 2 and Demand 3 are named H2 and H3, respectively. The detailed description of grid in these two cases is given in Table 9.
Here, we use the same numerical simulation and data analysis strategy as discussed in section 4.3 and 4.4. The time-averaged results are shown in Figure 17, and case H2 and H3 meet the requirement of 'Criterion 1'.
The normalized frequencies of the main POD modes or the first two modes are given in Figure 18. The primary frequencies of H2 and H3 are the same as that of B, meaning they meet the requirement of 'Criterion 2'.
For frequencies contained in secondary modes, the results are presented in Figure 19. The high frequencies of H2 and H3 are the same as that of B, while the low frequencies show a difference. The Res time (B, H2) of low frequency is 1.4 × 10 −4 and is greater than et, which means that only H3 meets the requirement of 'Criterion 3', as we expected.   Overall, it can be concluded that the grid generation guidelines obtained by the grid resolution assessment method proposed in the paper are feasible and reliable.

Conclusion
In this paper, an innovative grid resolution evaluation method suitable for hybrid RANS-LES simulation was proposed, and couples the grid resolution evaluation with the interest quantity during the turbomachinery design and mechanism analysis. Three grid resolution criteria considering different accuracy demands are proposed below, making flexible choices according to the research problem, and it has been verified employing a classic scenario. Commonly, Criterion 2 is sufficient.
Moreover, a parameter named the average-normalized grid-scale parameter ( l + ) was proposed to quantificationally describe the grid. By quantifying the grid sizes that meet different criteria, the blindness of gridding can be reduced to a certain extent. A gridded partitioning strategy was also given to lower the consumption of computing resources in turbomachinery simulations. According to the flow characteristics of multi-walls, adverse pressure gradient, and complex flow phenomenon, e.g. flow separation, shock-wake, vortex interaction, and shock-wake vortex interference, the computational domain is divided into four parts: the inlet region (IR), the RANS/LES region (RL), the RANS region (RR), and the focus region (FR), which can further lower the consumption of computational resources.
In addition, the T106A turbine cascade was used as an example to demonstrate the application of the grid resolution evaluation method, confirming the reliability of the method. Three grid schemes were finalized, which lowered the grid number demand up to 53% compared to LES standard in the example. The quantitative grid generation guidelines can be extended to simulations with similar geometries and the same numerical solver, which can reduce the blindness of grid generation.
The limitation of the proposed method is that for a new problem, the computational cost for quantitatively determine the grid scale is still expensive. While once the grid-scale is determined, it can be employed as a guideline while generating a grid of similar geometry without repeating the grid resolution evaluation. Although the grid resolution evaluation method proposed in this paper has been successfully applied to a turbine cascade with wake and small-scale flow separation, the grid resolution for complex flow phenomenon in turbomachinery with three-dimensional separation and high adverse pressure gradient should be analyzed further.