A numerical insight into contra-rotating open rotor in-plane loads

– This paper presents the comparison of results obtained from two codes with diﬀerent degrees of accuracy and computational time on the case of a contra-rotating open rotor (AI-PX7). The main goal is to assess the capacity of a code based on singularity methods to predict global open rotor performance, the in-plane loads (also called 1P loads), and the blade load distribution during a full blade cycle. Moreover, installation and unsteady correction models have been implemented in this code in order to evaluate the possible improvements of initial predictions. Finally, a numerical method to quantify some of the mechanisms governing open rotor thrust and 1P loads is proposed.


Introduction
Due to the increase in fuel prices, airframe and engine manufacturers are looking for step-changing technologies that might enable them to produce more efficient and cleaner aircraft by 2020 and later. Among a number of emerging concepts, the contra-rotating open rotor engine (CROR) is likely to be an interesting option for the significant reduction of the aircraft fuel consumption and polluting emissions (Fig. 1). Nevertheless, a number of key aspects are not still sufficiently mastered from the early development steps of the engine. Therefore, manufacturers are investing in the development of methodologies adapted for preliminary design phases in order to assess some of these key parameters.
As a step to provide these reliable methodologies for preliminary design of CROR, a code based on the unsteady lifting-line theory (HOST [1]) has been assessed by comparing its results with the ones obtained with a more complex and accurate CFD code solving the unsteady Reynolds-averaged Navier-Stokes (URANS) equations. The test case used for comparison is AI-PX7, an Airbus' generic contra-rotating open rotor geometry, at high-speed conditions and at 1 • of incidence. Rotor performance, blade loading, and induced velocity fields over a rotation are compared, in order to show the strengths and weaknesses of the lifting-line code with respect to a Corresponding author: i.gonzalez.martino@gmail.com the CFD solution. Finally, a method to quantify the impact of different mechanisms on thrust and 1P loads by induced velocity decomposition is proposed and applied to the same test case. The force developed by a propeller comes from the integration of the forces developed by each of its blades. If these forces are not the same for all blades along a rotation, the resultant force vector may not be aligned with the rotation axis. This means that, in addition to the thrust which is aligned with the rotation axis, another component of the force appears on the propeller plane. This in-plane force is commonly called the 1P load and can be decomposed in a norm and a phase lag. The phase lag is defined as the angle between the 1P load vector and the plane formed by the freestream velocity and the rotation axis, as shown in Figure 1.
As the magnitude of these 1P loads can be of the order of the thrust, they constitute a dimensioning factor on the engine installation systems, on aircraft handling qualities, on structural vibrations and fatigue limits, etc. Thus, it is critical to understand the major aerodynamic mechanisms behind these 1P loads and to predict them accurately from the pre-design simulations.
HOST was conceived by Eurocopter as a comprehensive code for the aeromechanical simulation of helicopters. It uses a modular structure in order to be able to simulate from an isolated rotor to a "complete" helicopter. This modular structure makes it easier to extend the use of HOST to other applications such as single propellers and Axial freestream velocity (m.s −1 ) BPF Blade passing frequency (s −1 ) URANS Unsteady Reynolds-averaged Navier-Stokes equations contra-rotating open rotors. The aerodynamic computations in HOST use a method based on the lifting-line theory, a singularity method for which the blade is reduced to its quarter-chord line. In the results hereafter, induced velocities from the wake and the surrounding blades are calculated using MINT, a high-order free-wake module developed by Onera [2]. As it has been detailed in a previous communication [3], an unsteady airfoil model for curved and swept blades has been implemented in this lifting-line code. This airfoil model is based on the unified lifting-line theory developed by Guermond and Sellier [4] and following the adaptation to numerical applications proposed by Devinant [5]. The global problem is solved by matched asymptotic expansions between the outer domain, a classical lifting-line problem, and the inner domain, the airfoil and its wake at the considered blade station. The inner domain is solved using a linearized unsteady thin-airfoil theory [6] which adds, to quasi-steady airfoil polar data, the extra lift produced by the pitching and plunging motion of a infinite-span flat plate.
Finally, the effect of the nacelle on the blade aerodynamics can be taken into account thanks to a velocity perturbation field, obtained from a blade-off nacelle CFD simulation.
Two main contributions are therefore obtained from this paper. On one hand, as public experimental data are not yet available, they provide a first step in the detailed assessment of a simplified lifting-line code by comparison with CFD computations. On the other hand, a method for quantifying the impact of a number of aerodynamic mechanisms on thrust and 1P loads has been developed, providing a useful methodology for blade designers. Figure 1 shows the AI-PX7 open rotor geometry 1 used in the present study. Table 1 presents the main geometrical characteristics and the flight conditions of the considered case. Four different simulations are presented in this paper: three different lifting-line simulations using the HOST code are compared to a URANS CFD computation using the elsA solver performed by François et al. [7].

Test case description
The CFD computation uses a mesh with 53 million nodes, divided in a fixed far field domain and two rotating domains, one for each rotor, connected by the sliding mesh technique. This technique allows relative motions between the non-coincident interfaces and interpolates the fluxes between them.
Different lifting-line simulations are performed to assess the effects of installation effects and an unsteady airfoil model. First, the "isolated" lifting-line simulation considers only the propeller blades, and although no nacelle is considered, the freestream Mach number is slightly increased (from 0.73 to 0.75) to simulate the mean acceleration due to the nacelle. Second, the effects of the hub, also called "installation effects", can be taken into account by introducing a steady velocity perturbation field for both rotors. This field is obtained from an initial CFD simulation of the nacelle in blades-off configuration. Finally, "unsteady corrections" stand for the addition of new terms in the airfoil lift calculation so as to consider its unsteady motion and three-dimensional blade forms.
While the CFD computation has been performed in 8 days on 128 cores and with a time step equivalent to a rotation of 0.5 • , lifting-line simulations have been performed in 11 h on 8 cores and with an equivalent 2 • time step, for a matter of reasonable computational costs in a design process.
In fact, the computational cost of lifting-line simulations are of the order of O(N 3 step ) due to the simulation time, O(N step ), and the free convection of the wake, O(N 2 step ), where N step is the number of time steps. Moreover, as in the lifting-line approach there is no interpolation between rotors and no viscous diffusion, the wake propagation is less sensitive to time step in lifting-line simulations than in CFD computations. Thus, a 2 • time step simulation has proven to be a good compromise to capture the main part of the interaction between rotors, i.e. the effect of incidence. Indeed, as the main frequency   of the incidence effect is a full revolution, 180 time steps are computed per period, which is largely sufficient to capture this mechanism. Moreover, blade-wake interactions might also be captured with a certain accuracy, as 8 time steps are computed per blade passage period.

Code assessment: comparison between elsA and HOST/MINT
This section presents the main results obtained in the lifting-line simulations and compared with the reference CFD simulation. First, open rotor performance is presented, then global blade thrust is compared, then blade thrust distribution is presented, i.e. the contribution of each section to the thrust, and finally induced velocities in the field are compared. Table 2 shows the propeller performance comparison between CFD simulation and three lifting-line simulations: "isolated", with installation effects, and with both installation and unsteady corrections. For HOST simulations, pitch angles have been adapted in order to obtain the same global thrust levels T = T FR + T RR and power ratio P FR /P RR than in the CFD simulations. An acceptable modification of the blade pitch setting, i.e. around 1 • , is required in order to achieve the thrust and the power ratio target levels.

Rotor performance and 1P loads
Three non-dimensional parameters are analyzed in this section (see Eq. (1)): the thrust coefficient, the 1P load coefficient, and the 1P load phase lag for front and rear propellers, respectively.
(1) Figure 2 plots the comparison between the different simulations. The total thrust level is obtained with a slight overestimation of the front rotor thrust (+3.5%) and a slight underestimation of the rear one (−3.5%). These differences can be explained by a lack of slipstream from the rear rotor on the front rotor. Indeed, as the lifting-line method models the blade as a vortex filament at the blade quarter-chord line, the effects of the chord and the blade thickness cannot be captured. As the distance between front and rear blade quarter-chord lines is more important than between the front blade trailing edge and the rear blade leading edge, the potential interaction effects between rotors may be underestimated. In this open rotor configuration, even if results are still satisfactory, this deficiency starts to be noticeable. When comparing the 1P loads of the front rotor, it can be noticed how installation effects play an important role in predicting better its norm: the offset is reduced   from −17% to +10% for the worst case. Rear rotor 1P loads norm is slightly modified by installation effects and the unsteady airfoil model.
The 1P load phase lag of the front rotor is reduced by the effect of the hub as it induces mainly a vertical velocity component. This leads to an increase in the vertical component without modifying the side component, thus increasing the 1P load norm and reducing its phase lag. The offset of front rotor phase lag is increased due to installation effects and the unsteady model, but offsets are still reasonable, i.e. under the azimuth step of 2 • . Mean values of the rear rotor 1P load phase lag are closer to CFD results when using the unsteady airfoil model, i.e. from +12 • to +7 • . Indeed, the airfoil model adds some lift components that tend to accelerate the lift evolution under a local incidence variation, thus they reduce the propeller 1P phase lag.
Oscillation amplitude due to the passage of the front rotor wake, however, is slightly overestimated. Three factors may explain these differences in the oscillations due to the blade-wake interactions: (a) differents time steps, (b) neglecting viscous effects, and (c) reducing the blade to its quarter-chord line.
First, with a time step corresponding to 8 steps per blade-wake interaction, HOST simulations can capture the correct frequency but not necessarily the oscillation amplitude. Second, HOST uses a potential wake represented by a series of vorticity panels without thickness, which means that no viscous effects are considered. Finally, as in the lifting-line approach the blade is reduced to its quarter-chord line, each blade section is represented by a point. The induced velocities from the wake on this point represent the impact on all the airfoil. This generates a more impulsive response to the wake passage than in the CFD simulation, where the influence of the wake varies all along the blade airfoil.  rear blade loads in all lifting-line and CFD simulations. Notice however that for the front rotor the thrust coefficient is slightly overestimated in the downward moving blade (between 60 • and 150 • ) when installation effects are added. This offset can also be attributed to the lack of slipstream from the rear rotor.

Blade thrust
Notice also the lack of oscillations in lifting-line results. As explained before, this can be explained by the fact that lifting-line methods reduce the blade to its quarter-chord line, losing all chordwise information. As the front quarter-chord line is further from the rear rotor than its trailing edge, effects from rear rotor are not well captured.
When observing rear blade loads, a slight underestimation appears in the upper part of the blade cycle for all HOST simulations. Unsteady corrections tend to increase blade loads around 90 • and around 270 • , where the incidence variation is the most important. On rear blades, oscillations are fairly captured in lifting-line simulations due to the passage of the front wake across the rear rotor. Notice however that a certain phase lag appears between the peaks, which is due to the point of emission of the wake: while in elsA it is convected from the trailing edge of the blade, HOST convects the wake panels from the quarter-chord line. Figure 4 shows the thrust coefficient distribution along the blade ∂τ/∂ξ for the same four simulations. Thrust coefficient distribution is defined as follows:

Blade thrust distribution
Left-hand figures show the mean value of the thrust distribution for front and rear rotors. Fair predictions are obtained in lifting-line approach compared to CFD results for both rotors. However, the wall effect near blade root is not captured by the lifting-line code as no hub model is implemented in it, i.e. the loads are imposed to zero at the root section. In the central part of the blade, ∂τ/∂ξ is overestimated. Finally, the blade maximum loading is underestimated (−8%) and slightly closer to the blade tip. No significant difference is observed between the different lifting-line simulations. Right-hand side figures show the first mode of ∂τ/∂ξ, i.e. the incidence effect. Installation effects increase these oscillations on both rotors. The unsteady model corrects slightly the modulus on the rear blade and reduces the phase lag, getting significantly closer to CFD results over an important part of the blade.
Higher modes corresponding to the BPF are not shown as their amplitude is not well captured in HOST and it is negligible with respect to the incidence mode.

Induced velocity fields
This paragraph compares the axial and circumferential induced velocities as predicted by CFD and lifting-line simulations for two planes normal to the rotating axis: one upstream the front rotor and the other one between the rotors (see Fig. 5). The circumferential component for the first and the second plane are calculated respectively in the direction of rotation of the front and the rear rotor.
For the plane upstream the front rotor (Figs. 5a and 5c), important mismatches are observed. Nevertheless, as it has been shown in the previous paragraph, the comparison between blade loads predicted by both simulations are very similar. To explain this it must be considered that the upstream plane is very close to the leading edge of the blade, and so the passage of the blade is very well captured in the CFD simulations. It is not the case in the lifting-line ones, as the blade is reduced to its quarter-chord line. Therefore, this comparison shows that volume effects are not negligible in the field near the blade.
In Figures 5b and 5d, even if some differences can be observed, the general pattern is well captured by the lifting-line approach, specially the circumferential component. The axial component calculated by the lifting-line code presents higher values in the zone of the downward moving blade, which is less evident in CFD results. On the contrary, this zone in the circumferential component is noticed both in lifting-line and CFD results. As the lifting-line theory considers a potential wake, no velocity deficit is observed due to the viscous effects in the wake. Moreover, similar to what was observed for rear rotor blade loads, a certain phase lag in the tip vortex position can be observed due to a difference in the emission points of the wake.

Aerodynamic mechanisms of 1P loads
The last section of this study is devoted to the analysis of the aerodynamic mechanisms behind open rotor 1P loads by using a number of lifting-line simulations.
From a purely geometrical approach, 1P loads can be attributed to a difference in the relative velocities of the downward and the upward moving blades, which generates a difference in their loads and thus a net force in the propeller plane. However, the impact of induced velocities, v ind , on 1P loads cannot be neglected as it has been shown in Section 3. To analyze their contribution, the v ind are going to be decomposed in several components. Under linearity hypothesis, a quantification of the impact of one single v ind component on thrust and 1P loads can be estimated by removing that component in a new lifting-line simulation and comparing the results with respect to the original lifting-line simulation. The induced velocities represent the influence of potential vorticity wakes, installation effects and airfoil motion on the airflow on a particular point of the field. To obtain the aerodynamic behavior of a propeller, lifting-line methods calculate the v ind on each blade section. Thus, a decomposition of v ind can be done as follows: (a) unsteady corrections due to the airfoil motion; (b) installation effects due to the hub; (c) auto-induced velocities, v ind by the wake of a propeller on the propeller itself; and (d) mutually-induced velocities, v ind by the wake of the other propeller. Both, auto-and mutually-induced velocities are decomposed in mean value (mode 0), 1/rev or incidence mode (mode 1), and the Blade Passing Frequency mode (mode BPF). Figure 6 plots the contribution of each v ind component on thrust, 1P load norm, and 1P load phase. Notice that the percentage value accounts for the importance of the impact of each component, whereas the sign accounts for its positive or negative contribution. Moreover, percentages in red should be considered carefully, as they are a consequence of non-linear effects, as it will be explained hereafter. These results have been obtained by storing the v ind field for a complete CROR simulation and then removing one v ind component at each simulation, in order to quantify its impact.

Aerodynamic mechanisms of thrust
The mean value of v ind is composed mainly of an axial and a circumferential component, i.e. the swirl. In the present method both components are considered together. Notice this is the main mechanism impacting the propeller thrust, but in a different way between autoinduced and mutually-induced velocities. On one side, in auto-induced velocities both axial and twist decrease the local angle of attack and thus the thrust. On the other side, in mutually-induced ones twist increases the angle of attack whereas axial component decreases it, but the overall effect tends to increase the angle of attack.
The front rotor thrust predicted by the lifting-line simulations is almost not affected by rear rotor wake. Thus we can consider that its response will be similar to the case of a single propeller in incidence. On the contrary, the rear rotor thrust is more impacted by the front rotor mean v ind (bar 5) than by its auto-induced velocity (bar 2). Furthermore, the contribution of the swirl of the front rotor wake increases the blade incidence and thus the rear rotor thrust.
Finally, as hub tends to accelerate the airflow in the axial and upward directions, it has a negative impact on the thrust of both rotors (bar 1), though it is less important than the wake components (−9%).

Aerodynamic mechanisms of 1P load norm
The main v ind mechanism impacting the 1P load norm is the first mode of the wake (bars 3 and 6), due to the propeller incidence. Indeed, the vorticity shed in the wake is proportional to the evolution of the angle of attack of the blade and acts in the opposite direction. For example, an increase in the angle of attack generates a wake with more vorticity and more v ind which tend to reduce that angle of attack. As the 1P load norm is directly related with the oscillation of the angle of attack of the blades, both auto-and mutually-induced velocities tend to reduce this 1P load norm. As for the thrust, the most important contribution on the 1P load norm of the rear rotor is the first mode of the front rotor wake. This helps to explain what was noticed in Section 3: that front rotor wake reduces the rear rotor incidence and thus the 1P load norm.
Notice the important contribution of mean v ind (bar 2). As plotted in the scheme of Figure 7, this is due to non-linear airfoil data, which is not desired in the present linear analysis of 1P loads. Indeed, when removing mean v ind from the input perturbation files, local incidence increases up to angles after the airfoil stall. As it is shown in the figure, for the same incidence oscillation, lift coefficient almost does not oscillate in the used airfoil data. Consequently, the contribution of mean v ind on the 1P load norm and phase lag is supposed to be very important, but it would not be the case if the angles of attack remained around the working point, in a linear part of the airfoil data. Thus, in the 1P load norm and phase plots the percentages of mean auto-and mutually-induced velocities (bars 2 and 5) should be considered carefully.
Finally, the hub in incidence has a positive contribution on the 1P load norm as it accelerates the airflow upward (see the scheme in Fig. 8), increasing the difference in relative velocity seen by the downward moving blade and the upward moving blade, and thus increasing the rotor vertical force and consequently increases the 1P load norm (bar 1).

Aerodynamic mechanisms of 1P load phase
One of the most important contributions comes from the installation effects (bar 1). Indeed, as it has been explained in the previous paragraph, the hub increases the vertical force generated by the propeller without modifying its side force. This is why installation effects tend to reduce the 1P load phase lag.
Notice also that the unsteady airfoil model has an important positive contribution on the phase lag (bar 8). In addition to the v ind from the airfoil motion, this model includes the effect of the near wake. The first one tends to accelerate the lift evolution when the incidence is modified, whereas the second one tends to lag it. The overall contribution is an increase in the lag between incidence and lift, which generates an increase in 1P load phase lag.
As explained before, mean v ind has a non-linear impact on 1P load norm and phase lag, and should be considered carefully (bars 2 and 5).
The first mode of v ind (bars 3 and 6) tends to decrease the amplitude and lag the blade loading oscillation due to the propeller incidence. This explains why the 1P load phase lag of the rear rotor is around two times the one of the front rotor. Figure 9 tries to illustrate how the front rotor induced velocities contribute to lag the angle of attack of a blade section of the rear rotor. First, the dashed curve represents the angle of attack of a single propeller induced by its wake; its minimum value is found in the downward moving blade zone, because there the wake has more vorticity. Second, the dash-dotted curve represents the angle of attack induced by the front rotor wake; its minimum value is placed around 270 • , where the front rotor wake has more vorticity. Finally, the solid curve shows the addition of both contributions, as in the case of a CROR. It presents a reduction of the amplitude of the Δα and an increase in its phase with respect to the single propeller case.
Finally, the contribution of BPF modes in lifting-line simulations is negligible (bars 4 and 7). However, an analysis of the sensitivity of the 1P phase lag to time step presented by François et al. [7], it was noticed that BPF modes had an impact of around 3.5 • for the front rotor (−17%) and 3 • for the rear rotor (−7%). These offsets can be a consequence of the three reasons exposed in Section 3: using different time steps, neglecting viscous effects in the lifting-line approach, and reducing the blade to its quarter-chord line in the lifting-line approach.

Conclusions
A detailed code-to-code assessment between elsA CFD solver and HOST lifting-line code, has been performed on a generic contra-rotating open rotor geometry (AI-PX7) at high-speed conditions and at 1 • of incidence. Three different types of HOST simulations have been performed in order to assess the impact of nacelle effects and unsteady corrections on aerodynamic performance and 1P loads.
Total thrust and power ratio target values have been achieved for all lifting-line simulations within reasonable pitch angles modifications, i.e. ∼1.36 • . On the one hand, offsets obtained in the 1P load norm have been reduced when adding the effect of the nacelle. On the other hand, the unsteady airfoil model implemented in the lifting-line code allows reducing mismatches on the rear rotor 1P loads phase lag.
Fair results have been obtained in global blade thrust. Nevertheless, the comparison has also highlighted some of the shortcomings of the present method in the prediction of rear potential effect on the front rotor. The mismatches with respect to the CFD results have been explained by the fact that chord effects are not considered in lifting-line methods, and they become important for a rotor close to another.
In the comparisons of thrust distribution along the blade, the shape predicted by the CFD simulations was also obtained with the lifting-line ones, but with some minor offsets. They were noticed near the blade root due to the wall effect. Moreover, acceptable differences in maximum mean loads were also observed (−8%) and their position was slightly closer to the blade tip than in CFD simulations. The first harmonic was well captured by all lifting-line simulations, though installation effects overestimate its amplitude on the front blade. Unsteady corrections reduced the phase lag of the first mode, and hence results get closer to the CFD predictions.
Induced velocity fields have shown, on one side, the impact of neglecting chord effects in the lifting-line approach and, on the other side, the capacity to capture correctly the potential effect of the tip vortices and the global effect of the front wake on the rear rotor.
A method using the HOST code to quantify the impact of a number of aerodynamic mechanisms on the open rotor thrust and 1P loads has been proposed. Though some non-linear effects appear in some of the results, the method is still valid and enables to determine the quantitative impact of each mechanism on global performance of open rotors. Front rotor behavior predicted by the lifting-line simulations is very close to a single-rotating propeller, as the rear rotor has a very small contribution on its performance and 1P loads (<10%). These results have to be analyzed taking into account what has been exposed in the Section 3, i.e. that the lifting-line approach neglects all blade chord and thickness effects, leading to an underestimation of the potential effect of the rear rotor on the front rotor.
On the contrary, the most important mechanism on rear rotor response is the front rotor wake effect. Indeed, more than 50% of its thrust is due to the combination of positive axial and negative circumferential v ind from the front rotor wake. Moreover, this wake decreases the 1P load norm (−41%) and increases its phase lag significantly (+18%).
Finally, while installation effects play an important role in the prediction of 1P loads for both rotors, the unsteady corrections have shown to be important for the 1P phase lag, mainly in the case of the front rotor, where more important incidence variations occur.