Evaluation of vortex-blade interaction utilizing flow feature detection techniques

The extraordinary growth in computational capabilities over the last few decades has enabled the numerical simulation of massive and complex flow problems with high accuracy. Not surprisingly, CFD has become a crucial tool in the design of pioneering aircraft engine architectures, such as CROR engines. Noise and performance requirements lead the design process from a very early stage, thus requiring deep investigation of the acoustic and aerodynamic behaviour. Monitoring the trajectory of the vortices generated at the tip of the front rotating blades is of critical importance to understand and prevent vortex-blade interaction with subsequent stages, as this non-linear flow topology strongly influences the aerodynamic performance and acoustic footprints of the engine. This manuscript follows a flow feature-based approach to visualize and track these coherent structures, for the particular case of CROR. The paper starts by evaluating the suitability and performance of four typical region-based (RB) vortex detection criteria, which delimit rotating flow regions, and one line-based (LB) method that allows the reconstruction of the imaginary centre lines of the vortices. Then, two novel methodologies are introduced that improve the original assortment of seeds required by the tested LB method, as they increase the probability of these seeds to grow into a tip vortex line, providing faster and more accurate answers during the design-to-noise iterative process.


Introduction
In recent years Counter Rotating Open Rotor (CROR) engines have become of prime interest in the aeronautical industry, in search for more efficient aircraft configurations. Amongst all potential future propulsion systems, the open rotor concept offers the most significant fuel consumption and CO2 emissions improvement. The biggest drawbacks in these particular engine types are the high levels of noise generated, both broadband and tonal, due to the absence of a cold stream duct, thus design focus is placed on trying to reduce it in order to meet tighter noise regulations within the airspaces. As it is well known a major contribution to tonal noise is caused by the first stage rotor blade tip vortices interacting with the second stage rotor. This impact condition is also undesirable from the aerodynamic and structural point of view, as it penalizes drag, and increases significantly vibrations and fatigue loads nearby the impact regions. On the other hand, broadband noise is directly associated with the breakdown of trailing vortices. An effective and efficient tracking of these two vortex topologies is therefore indispensable to the designers that look for advantageous configurations. Thanks to the developments in Computational Fluid Dynamics (CFD) techniques and the exponential growth in computational power, today it is possible to obtain with high accuracy detailed flow behaviour predictions under normal operating conditions. Not surprisingly, over the last few decades the reliance on CFD in the design process of pioneering aircraft engine architectures has increased significantly. Noise and performance requirements lead this design process from a very early stage, thus requiring deep investigation and complete understanding of the acoustic and aerodynamic behaviour within numerical solutions. A new problem arises from current industrial simulations, which is the large amount of CFD data that needs to be processed to be able to derive valuable conclusions. In particular, the simulations focusing on noise prediction require very small time steps and very large meshes which increases the burden on data reduction prior to analysis on the part of the designer. Most methodologies in noise prediction move from the CFD analysis directly to noise propagation models based on pressure data around known noise sources, which gives quantitative noise information in the form of sources and footprints at high computational costs [1]. Alternatively, flow feature detection techniques focus the analysis on regions of relevance within large CFD data domains. These techniques aim to extract flow features that are of interest for a certain engineering problem in a numerically efficient manner. Tip vortices and vortex breakdown are the features of interest for the particular case of CROR engines in the current design-to-noise mindset, where ideally original data would only be accessed for the reconstruction of such structures. Without adding excessive computational burdens to the CFD, this methodology offers an important size and complexity reduction to the original data and a higher level of abstraction [2,3], providing faster and more accurate answers to the designers. In this manuscript, one first evaluates the suitability of traditional region and line-based vortex detection techniques, available in the literature, for the particular case of CROR. Four region-based (RB) methods based on the velocity gradient tensor and one line-based (LB) algorithm were selected for the visualization and tracking of the tip vortices produced at the first engine stage and interacting with the second. The strengths and weaknesses of each individual algorithm and its computational performance are also discussed in this manuscript. Following this comparative analysis, two new initialization methods are presented to improve the original selection of seeds required by the tested LB method. The first methodology uses the information given by the four RB detection methods previously tested as an input for the aforementioned LB method. For the second method, the selection of seeds is driven by high friction drag and high pressure gradient criteria.

Numerical methodology
The CROR test case used for this work comes from a Z08 open rotor configuration, based on the geometry of the generic 14ft diameter open rotor AI-PX7 from Airbus Operations at scale 1/7. It features eleven blades in the first rotating row, and nine blades in the following counter-rotating part. A multi-block structured grid was generated, comprising 825 blocks and around 95 million grid points, and modelling the entire annular section. The first rotating row was connected to the second counter-rotating row through a sliding plane. Fig.  1 shows the geometry, the computational domain and the structured mesh of the CROR test case used in the present work. For additional details, the reader is directed to [4].
The elsA solver [5,6], developed by ONERA, was used to perform time-dependent, compressible Reynolds-Averaged Navier-Stokes (URANS) simulations, adopting the oneequation Spalart-Allmaras turbulence model. The computations allow the simulation of the actual CROR engine conditions in terms of Reynolds and Mach numbers (Re = 4.5x10 6 and M = 0.2). The working fluid was assumed to be a thermally and calorically perfect gas. Spatial discretization was carried out using the cell-centred finite volume method with second-order upwind. The Roe scheme was used to compute the viscous fluxes. For the temporal discretization, the implicit dual time stepping method was selected with backward Euler integration with implicit LU schemes. The CFL number was adjusted to 1. Four RB vortex detection methods relying on the velocity gradients (Q-criterion [7], Kinematic vorticity number (Nk) [8], Δ-criterion [9] and λ2-criterion [10]), and one LB scheme (Predictor-Corrector (PC) method by Banks & Singer [11]) were selected and implemented within an internal Airbus object-oriented post-processing library. This tool follows a pipeline programming philosophy written in Python and making use of scientific and visualization libraries, such as NumPy, SciPy and VTK. Detailed description of each of the selected vortex detection methods is given in [12].

Search of vortex core regions
The iso-surfaces of different threshold values for the parameters Q, Nk, Δ and λ2 applied to one snapshot of the CROR solution are shown in Fig. 2. For the lower absolute threshold values selected, all these four methods successfully predict the trajectory and the shape of the vortices developed at the tip of the blades comprising the front rotating rotor, and a tipvortex impact condition with the subsequent rotor stage can be observed. Downstream of the second stage significant vortex interactions occur, and none of the four RB methods is capable to provide a clear representation of all the individual coherent vortex regions. As one increases the absolute value of these four parameters the low energy structures can be filtered, allowing a clearer representation of the tip vortices. However, if the selected thresholds are too high, the trajectory of the tip vortices cannot be visualized correctly. For the highest threshold values of Q and Δ, tip vortices are not extracted. However, these two RB methods extract always regions representing boundary layers. This is not the case for the other two RB methods. Applying the kinematic vorticity method, boundary layer regions are removed almost automatically. Using the λ2-criterion, one is able to vanish the boundary layer in the shaft imposing strong negative values. All the four RB methods are applied locally to each grid cell, making them relatively easy to parallelize across multiple processors. When passing from 2 to 12 processors the running time decreases about 85%. As expected, the most time consuming method, the λ2-criterion, requires about 10 minutes of computation time using 12 processors. These figures would drop considerably if one takes the same number of processors used for the CFD simulations, making the post-processing time required by these RB methods irrelevant in comparison to the simulation.

Search of vortex lines
In order to investigate the suitability of the PC method for the particular case of CROR engines, three representative candidate seeds (P1, P2 and P3) were selected relying on both low static pressure and high vorticity magnitude thresholds. The original PC method was then initiated individually on each of these three candidate seeds. Fig. 3-left shows (in blue) the vortex axis developed from seed P1 using the PC method, together with iso-surfaces of λ2 = -5 (in red). This line clearly defines the skeleton of the corresponding vortex core region extracted with the λ2-criterion, both in upstream and downstream directions. In the backwards direction to the flow, the line grows from its starting point P1, to the onset location of the matching tip vortex TV1. Initially, in the forward direction to the flow, the skeleton line follows correctly this tip vortex, until there is an impact with a second stage blade. At this intersection point the line extracted for TV1 stops its growth, and then a different skeleton line belonging now to the tip vortex TV2 starts to develop in the backwards direction, stopping at the exact location where TV2 is originated. From Fig. 3left one can also observe some vortex line oscillations near the sliding plane region, that the PC method cannot avoid. Similar results can be obtained starting the PC algorithm from seed points P2 and P3, as shown in Fig. 3-center and Fig. 3-right respectively. Note that P2 is located near the midregion between CROR Stage 1 and Stage 2, while P3 is close to the tip of a front rotating blade of Stage 1. The vortex lines obtained from these two candidate seeds have a very similar shape and are almost superimposed in space, as they belong to the same tip vortex TV3. These results are shown here only as an illustrative example about the suitability of the PC method to the present CROR test case, once the original algorithm already expects that a pre-selected candidate point belonging to an already extracted vortex is eliminated from the list of possible candidate seeds, avoiding the growth of the same vortex (as seen in Fig. 3-center and right) and any possible redundancies, and thus less examinations. The biggest problem encountered with the original PC method is related with its high computational demands. At each pressure correction step, the method needs to interpolate over the computational domain in order to reach the local pressure minimum, and for a big test case such as the full annulus CROR this process can take long computation times (vortex core line growing from P1 -16 min, from P2 -27 min and from P3 -23 min, using one single processor). Besides that, the method is extremely sensitive to the quality of the prior initialization. For the present CROR test case, the original selection of candidate seeds relying on pressure and vorticity thresholds provides a significant amount of points that do not evolve in a vortex line, increasing the already high computation costs of the method.

Improving the initialization of the PC method
The original PC method makes the selection of candidate seeds according to high vorticity and low static pressure threshold criteria. Normally this methodology is enough for simple academic test cases, with a low number of grid points. However, when dealing with a large scale industrial case with high flow complexity, such as the present CROR test case with 95 million points, the original initialization may return a huge amount of candidate points. This will directly penalize the subsequent PC step, once it has to be started from each one of those candidate points, resulting in excessive and prohibitive computational costs. Furthermore, by relying exclusively on those two threshold criteria, there is not a physical guarantee that the selected thresholds contain the most relevant features. For the present CROR test case, by setting a threshold for cells whose vorticity magnitude is higher or equal to only 1% of the maximum vorticity of the computational domain, this criterion fails to extract cells related to the tip vortex emerging from the first rotating row. The output of the aforementioned vorticity threshold is only associated with the boundary layer of the rotating surfaces, where vorticity can also be high. Alternatively, for the lowest static pressure values of the domain, the corresponding threshold filter retrieves no more than points or cells located in low pressure zones of the blades. Furthermore, significant trial and error tests were required to correctly tune a static pressure threshold that allowed the extraction of cells or points related to tip vortices. To address these limitations, two novel initialization methodologies are proposed.

Initialization based on RB vortex detection criteria
This initialization uses the previous RB vortex detection methods to build a set of candidate points for the PC method, increasing its probability to find candidate seeds that will develop into a tip vortex core line. With the inclusion of these four gradient-based methodologies, the extraction of candidate seeds being related with a tip vortex becomes also less sensitive to the selection of thresholds, contrary to what happens with the thresholds of vorticity magnitude and static pressure. For additional details concerning its algorithm, the reader is directed to [12,13]. For very high threshold values of Q, Nk and Δ tip vortices from the first rotating row are not extracted, as previously shown in Fig. 2. However, for all the selected thresholds of λ2 shown in this figure, cells belonging to tip vortices are always extracted, demonstrating that the suggested initialization is much less sensitive to the selection of those values, when confronted to pressure and vorticity magnitude thresholds. For the latter, note that even with a threshold of 1% the maximum vorticity magnitude of the computational domain, no cell associated with the tip vortex is extracted. And for the static pressure, only with a threshold from 6 times its minimum value inside the computational domain, the extraction of cells related with the tip vortex is possible. For higher thresholds of Q, Nk and Δ and lower values of λ2, the ratio of volumes extracted over the total volume of the computational domain decreases, lowering the number of candidate seeds extracted and thus reducing the computational burden of the PC method, once less evaluations are required. Moreover any threshold value shown in Fig. 2 of the aforementioned gradient-based parameters combined is able to extract cells belonging to tip vortices, enhancing the probability of the proposed initialization to extract candidate cells or points that have the potential to grow into a tip vortex core line.

Initialization based on high pressure gradients and high friction drag
Similarly to the previous method, an initialization relying on high pressure gradients and high friction drag allows to improve the quality and to reduce the uncertainty of the original initialization of the PC method, as it increases as well the probability of the aforementioned LB method to find candidate cells or points that will develop into a tip vortex core line. On the one hand, adverse pressure gradients in the direction of the flow induce the formation of the "separated" tip vortex from the first blade row. The extraction of cells or points where this parameter is higher than zero enhances thus the probability to find good candidate seeds. On the other hand, the proposed methodology searches also for surface cells or points where the drag friction coefficient is high, suggesting high local aerodynamic losses. In CROR engines, the formation of tip-vortices and possible tip-vortex-blade impact conditions have a significant contribution to the overall aerodynamic losses of the system. With the combination of these two superficial aerodynamic parameters, a more intelligent selection of candidate seeds is granted, but also the quality of the initialization process becomes less dependent upon arbitrary threshold selection, in comparison to the thresholds of vorticity magnitude and static pressure from the original initialization. The pseudo-code used for the presented methodology can be found in [12,14]. The proposed initialization demonstrates to be much less sensitive to the selection of thresholds, when confronted to pressure and vorticity magnitude thresholds. It also enables a significant drop of the number of candidate seeds extracted, allowing thus an important reduction of the computational burden of the original PC method for the particular case of CROR. As an illustrative example, Fig. 4 depicts the seeds extracted for an initialization based on high pressure gradients (a-b) and high friction drag coefficient (c-d). Fig. 5 shows the vortex core lines obtained with the PC method from the seeds illustrated in Fig. 4(d) using this new initialization scheme. Note that prior to the beginning of the growth of each vortex line, the original position of each seed is also refined in order to match the local point of maximum Nk. As it can be observed from Fig. 5, all candidate seeds correctly end up in lines that accurately describe the trajectory of the corresponding tip vortices. This illustration also confirms that all the extracted vortex lines are impacting the second stage counter rotating blades of the tested CROR configuration. The computation time required for the growth of one of these lines was, in average and using one processor, 8 minutes and 35 seconds.

Conclusions
For large and complex industrial simulations, flow feature detection techniques can provide fast and very accurate information about existing non-linear flow features. Using a typical industrial CROR test case, this work shows that both RB and LB vortex detection methods are very convenient tools to visualize and to track the vortices generated at the tip of the blades, providing feedback about possible vortex-blade impacts. All the four implemented RB methods successfully depict the trajectory of the tip vortices, nonetheless they strongly depend on fully arbitrary thresholds. Boundary layer removal is almost automatic using the Nk method, a valuable feature considering the boundary layer is of no interest to the current study. For the particular case of CROR, the tested PC method provides a precise tip vortex tracking and a near automatic impact evaluation, once its algorithm breaks when a vortex core line intersects a blade. However, the method cannot guarantee the correct development of a tip vortex line when it interacts with other rotating structures. Furthermore its original seeding selection (or initialization) is somehow limited. Two initialization methodologies are proposed to address the solution of this last problem, by increasing the probability of finding good candidate points that develop into a tip vortex core line. The first initialization scheme relies on the aforementioned RB methods, and can be applied to any generic dataset containing rotating structures. The second initialization method is geared towards particular flow problems where fluid-structure interaction exist, and where the formation of vortices is linked to high pressure gradients and friction drag on the structures, as in this CROR case.