BODY FORCE MODELING OF THE AERODYNAMICS OF A LOW-SPEED FAN UNDER DISTORTED INFLOW

New propulsive concepts such as Boundary Layer Ingestion involve stronger interactions between the engine and its environment, and thus more complex ﬂows compared to classi-cal architectures. Usual turbomachinery design tools are inadequate, and new numerical methodologies are needed to accurately predict the engine performance with affordable CPU resources. The present paper examines the relevance of a reduced-order modeling approach, the body force modeling (BFM) method, for a low-speed cooling fan with inﬂow distortion. The formulation itself accounts for the blade metal blockage, compressible effects and it relies on a physics-based loss model, independent of CFD calibration. The BFM results obtained in the present work are assessed against full-annulus URANS results and experiments. The comparison shows that the BFM approach successfully quantiﬁes the fan stage performance. Furthermore, the distortion transfer across the stage is exam-ined and the ﬂow patterns observed are found to be the same as in the URANS results and in the measurements. Hence, this methodology, coming at a low CPU cost, is well-adapted to the early design phase of an innovative propulsion system.


INTRODUCTION
Further reducing the environmental impact of aviation requires the investigation of new aircraft architectures, among which embedded propulsion is a promising solution. Also known as Boundary Layer Ingestion (BLI), this approach increases the overall propulsive efficiency and potentially offers significant power savings, for instance up to 8.6 % according to Uranga et al. (2017).
With BLI, the fan constantly operates under severely distorted inflow conditions, with large fluctuations of the upstream total pressure and swirl angle, among other quantities. This can have a strong impact on the fan performance, and capturing the correct massflow redistribution and spatial distribution of work across the blade row currently represents a challenge. Distortion transfer predictions across the stage are also important with respect to the performance of downstream components. This emphasizes the need to estimate accurately the performance of the different engine components as early as possible in a pre-design cycle. Full-annulus 3D unsteady CFD simulations can be seen as a reference method to capture the aerodynamics of such flows, but the cost of this approach still remains prohibitive for early design phases. In the present paper, we examine the accuracy of a reduced-order modeling approach: the body force modeling (BFM) method. By using source terms to reproduce the passage-averaged effect of the blades on the flow, without actually meshing the blades, this method enables to simulate a full-annulus stage configuration under distorted inflow with a steady approach and reduced mesh size, therefore drastically reducing the cost of the simulation.
The BFM approach has already been applied to BLI propulsor configurations in the literature, but some open questions remain. Hall et al. (2017) used the BFM method to predict the flow in the low speed fan studied by Gunn & Hall (2014), with validation against experimental results. However, their approach is purely inviscid, and their body force formulation does not model losses but uses a diffusion factor instead as a surrogate measure to efficiency. Kim & Liou (2017) studied the N3-X hybrid wing configuration with a body force model that accounts for off-design losses, relying on a significant calibration process. They show good agreement with full annulus reference calculations, but validation against experimental data is not available. In this context, the main goal of the present paper is to assess a body force formulation that include loss modeling but with minimum calibration, and to validate the results against reference unsteady calculations as well as experimental measurements.
First, the test-case is presented: it is a low-speed cooling fan, for which experimental data are available with and without inflow distortion. Next, the numerical methods and the BFM model are described. Global performance are then validated against experimental data. Finally, inlet distortion simulations are analyzed. Downstream of the rotor, the agreement of the BFM and full-annulus unsteady RANS is very satisfactory. Discrepancies with the experimental data are found close to the hub, but at higher spanwise locations, experimental and numerical results match fairly well. Downstream of the stator, the flow redistribution is also found to be reproduced by the BFM simulations.

Test rig
The test case is a low-speed cooling fan stage, designed and produced by Safran Ventilation Systems, mainly used to regulate the temperature of heat exchangers when the aircraft is on the ground. Its main characteristics are given in Table 1.
Following the framework proposed by Binder et al. (2015) for windmilling studies, the flow and loading coefficients, expressed at a mean quadratic radius (Eq. (1)), are used to assess the performance of the stage. Furthermore, a reduced flow coefficient is defined in the following according to Eq. (4), where the flow coefficient is referred to the free-windmilling condition, The experimental measurements are carried out in the Département Aérodynamique Energétique et Propulsion (DAEP) of ISAE-SUPAERO. More details can be found in Ortolan's PhD (2017). Two different configurations are used: • a first one (Figure 1, left) where the test rig is equipped with an asynchronous electric engine, which is located far from the test sections; this configuration corresponds to uniform upstream flow conditions.
• a second one (Figure 1, right) where another electric engine is inserted in the fan body, so that there is no shaft outside of the body and distortion grids can be placed upstream of the fan stage in order to create non-uniform flow conditions.

Distortion grid
In order to keep the distortion pattern simple, a unique porosity grid with regular 1 mm sidewise cells is chosen, to ensure two uniform 180 zones ( Figure 2). This grid is inserted between the bellmouth and the instrumented section, close to Plane 1 location, around three  Figure 1: Schematic of the two experimental configurations, without / with distortion times the fan chord length upstream from Section 2 ( Figure 1, right). CFD results indicate that the distance between Section 1 and Section 2 is sufficient to insure the absence of potential effects in the grid neighboring. Detailed information about the sensors used in this study can be found in Ortolan et al. (2018). Although no global measurement is available with distortion, local measurements are available in Sections 2 and 3 ( Figure 1, right). They consist in 15 equally spaced azimuthal positions over the annulus (every 24 ) which are obtained by rotating the distortion grid relatively to the probe position. For each azimuthal position, 15 radial measurements are carried out using the five-hole probes previously described.
A repeatability study carried out on the global steady-state experimental data for the nominal operating point of the fan stage showed that the repeatability error was greater than the measurement uncertainties. The same conclusion applies to the local measurements coming from the five-hole probe located at rotor outlet. For this reason, the repeatability error is considered here instead.

NUMERICAL SETUP
Two different types of simulations are introduced : unsteady full-annulus RANS simulations and BFM simulations. A meridional view of the numerical domain is shown in Figure 3. This section gives more details about the numerical parameters used in each case.

URANS simulations
URANS simulations are carried out with the elsA software, a cell-centered finite volume solver (Cambier and Gazaix (2002)). A multi-block structured mesh is used, comprised of around 60 million cells for the full-annulus domain. The total quantities and the absolute velocity direction are imposed at the inlet. A throttle condition with radial equilibrium is used at the outlet. All the walls are considered to be adiabatic and the first cell size ensures y + ⇠ 1. The one-equation turbulence model of Spalart-Allmaras is adopted. Since no hot wire measurement is available upstream of the fan stage, a value of 1 is imposed for the µt µ ratio at the inlet, which is typical of internal flows. The convective fluxes are calculated using a second-order Roe scheme. An implicit second-order Gear scheme with the Dual Time Stepping (DTS) method is employed for time integration, with a maximum of 20 sub-iterations per timestep or a decrease in residuals of two orders of magnitude. After a timestep convergence study, the number of timesteps per revolution is chosen equal to 2N R N S , i.e. 782, which is in line with in-house best practice, with regard to the compromise between CPU resources and time precision. This number of timesteps, 46 / 34 per rotor / stator blade passing period, respectively, is enough to capture the first BPF harmonic (Gourdain (2011)).

Body Force Methodology
Initially introduced by Marble (1964), the concept of BFM has been successively developed in later works (Gong (1998) for example) and is currently used in order to address various problematics such as fan-intake interaction in short nacelles (Peters (2014)), vortex ingestion (Bissinger and Braun (1974)), BLI distortion or even acoustics (Defoe and Spakovszky (2012)). The basic principle of this method consists in modeling the forces exerted by the blades on the fluid instead of directly simulating the flow around the complete three-dimensional geometry. With this approach, the bladed areas are replaced by source terms and the effect of the blades on the flow are usually decomposed into two different contributions: • a force normal to the relative flow field, f n , which is responsible for the turning; • a force parallel to the relative flow field, f p , which accounts for the loss. The main advantages of the BFM method are the following: • the 3D mesh is simply obtained by extruding the 2D meridional mesh in the azimuthal direction. Not meshing the boundary layers around the blades yields a very low cell count.
• more importantly, this method enables to deal with multistage configurations while keeping a steady resolution, which considerably decreases the associated CPU cost.
• finally, the general formulation is very flexible since the expressions of f n and f p are userdefined and can easily be modified. The local definition of the two force components also lets some room for extra calibrations, which can rely either on preliminary BFM simulations or on more accurate RANS results.
The BFM formulation retained in the present is adapted from Hall's model (2017), which basically does not include any loss model. Metal blockage effects are added, as well as a compressibility correction in the normal force definition. According to Thollet's work (2017), a parallel force is introduced too. It accounts for the loss with a local friction coefficient C f , derived from an empirical turbulent flat plate correlation, and a local chordwise Reynolds number Re z . Of course, no specific phenomenon is expected to be captured, such as an endwall corner separation or a shock wave -boundary layer interaction.
This loss model can be completed by an additional calibration coming either from blade or preliminary BFM CFD calculations, which corrects the loss amount by taking into account the flow deviation. This "off-design" term involves to choose a reference operating point from which the distribution of deviation is extracted. The main effect of this term is to steepen the performance curve at high and low massflow rate regions (near choke and surge, respectively). As explained in the following, only the basic friction contribution (2C f term) was kept for the present study.
BFM simulations are performed with the unstructured solver of the Fine/Open TM package of Numeca. The mesh itself is created in IGG TM and then converted to an unstructured format in HEXPRESS TM . An azimuthal step of 1 is chosen to extrude the 2D meridional mesh, which leads to 200.000 cells mesh over a 5 angular sector for the clean case and around 21 million cells over 360 with distortion. The BFM mesh follows the same y + criterion than the URANS mesh at the endwalls. The same inlet conditions are imposed, as goes for the turbulence modeling. At the outlet, the massflow value is imposed by adapting the static pressure and a radial equilibrium is prescribed. A second-order Roe scheme is used to compute the convective fluxes and a Runge-Kutta scheme with implicit residual smoothing is used with a CFL number of 1000 for pseudo-time integration.
The OpenLabs TM feature of Fine/Open TM enables to implement the BFM source terms in the right-hand side of the Navier-Stokes equations. As the solver is unstructured, all geometric inputs (blade normal vector components n z , n r , n ✓ , blade metal blockage and its derivatives b, @b @z , @b @r ) are calculated in a pre-processing step and then interpolated as (z, r) crossed polynomials.

Numerical distortion
As explained previously, the objective of this work is to evaluate the BFM capacity to deal with non-uniform flows. The distortion pattern is voluntarily chosen as simple as possible : it consists in a uniform 180 zone of 1.7% total pressure deficit, i.e. P i min P imax = 0.983 ( Figure 2). As already mentioned, this value is measured experimentally in Section 2, which is very close to the rotor LE (around 5 mm) and thus, makes it impossible to impose numerically inlet boundary conditions at this location. Consequently, without information about the upstream flow conditions in Section 1, it is hypothesized that between Section 1 and Section 2, the total pressure variations can be neglected and a similar distortion pattern is imposed as a 2D (r, ✓) P i map at the inlet of the numerical domain of both URANS and BFM simulations.
Since the distortion is fixed in the absolute frame, the inlet blocks of URANS simulations are non-rotating and a sliding mesh is introduced shortly downstream of the inlet, in order to ensure the propagation of the distortion in the rotor blocks. On the other hand, in the BFM case, the additional source terms are expressed in the absolute reference frame, so that no particular interface is needed.

RESULTS AND ANALYSIS
Validation of BFM in the clean case Figure 4 shows the loading coefficient and the non-dimensional total-to-total isentropic efficiency of the fan stage at three rotational speeds (6500, 8000 and 9500 rpm). The design operating point, not measured here, has been added just for illustration (in magenta). The linear trends in the evolution of described in Binder et al. (2015) are illustrated here.The agreement between experimental measurements and BFM results is within 3% for and remains satisfactory for the efficiency, except at high massflow coefficients. Experimentally, local measurements were carried out separately for the rotor and the stator, keeping the same value forˆ ⇤ but at different rotational speeds. For this reason, Figures 5 and 6 represent the upstream and downstream radial profiles of total pressure and absolute flow angle for each blade row. Here, the results are shown forˆ ⇤ = 0.66. The results of a fullannulus URANS simulation are also plotted. The stagnation pressure is mass-averaged whereas the absolute flow angle is area-averaged. The time averaging consists in a simple arithmetic averaging.
The total pressure profiles match very well (maximal relative discrepancy of 0.3%) but there is a visible difference between experimental and numerical results regarding the flow angle, especially near the hub. This is due to the fact that the rotating shaft upstream of the fan stage ( Figure 1) is not simulated. The URANS profiles are closer to measurements at the endwalls, as could be expected since the body force simulations cannot capture secondary flows like corner separations or rotor tip flow. The maximum flow angle discrepancy compared to URANS values is always lower than 1 degree.
It should be emphasized that the body force formulation lets room for improving these trends by using user-defined calibrations (see Eq. (8)). However, on the whole, these results were considered satisfying enough to directly evaluate the method in the presence of an upstream flow distortion without resorting to any additional calibration.

Evaluation of BFM with distortion
Only one operating point is analyzed in the present study, characterized byˆ ⇤ = 0.79 without distortion. Table 2 gives the relative discrepancy between global performance obtained from BFM simulations, taking the URANS values as reference here, for the stage and also for each row. Due to the changes in the experimental setup, no global measurement is available in this case. This comparison shows that the stator is responsible for the discrepancies between the two simulations. Thus, if a calibration had to be applied, it would be relevant to increase first the loss level inside this row.
Local experimental measurements with distortion are not available downstream of the stator row. In the following, the axial and tangential velocity distributions in Section 3 are plotted in order to compare numerical and experimental evolutions across the rotor (Figures 7 to 10). Qualitatively, both numerical results stand in good agreement with the measurements. Near the hub, co-swirl and contra-swirl are visible where the rotor enters and leaves the axial velocity deficit area, respectively (✓ = 270 and ✓ = 90 ). The difference in relative flow incidence at the rotor LE between high and low axial velocity zones is around 10 . Since there is no flow separation on any blade, the relative flow angle at the rotor TE is the same in Section 3 for both zones, which explains higher values of the absolute azimuthal velocity, V ✓ , in the low V z area. A shift of the local extrema can be seen in the azimuthal direction for V ✓ but the amplitudes are correctly reproduced (Figures 9 and 10). A small offset in V z is also observed at mid-span and near the shroud (Figure 8).
The only major discrepancy concerns the V z distribution at the hub. It is attributed to the strong curvature upstream of the rotor root due to the front bulb (see Figure 3). As illustrated  by Figure 11, the axial velocity of the flow near the hub in Section 2 is not the same in the CFD and in the experiments. An iterative procedure to set the numerical inlet boundary conditions would be required to match the measurements at this location. In the future, a redesign of the fan stage is scheduled : it should include an aerodynamic shaped spinner upstream the rotor and the stage itself should be more loaded, in order to be more representative of propulsive fans. the extremal values during a rotor blade passing period. The main feature that distinguishes both simulations here is the presence of the stator wakes in the time-averaged URANS field. Finally, the flow pattern is also similar far downstream, behind the rear bulb ( Figure 15).
The satisfactory agreement between both numerical simulations proves that the BFM is well-adapted to treat this sort of configuration, provided that only large flow scales (i.e. larger than one inter-blade passage) are expected in the work or loss spatial redistributions. The re- maining difficulty in the present study comes from the fact that the upstream flow conditions are unknown. A fine characterization of the distortion grid would be needed, to ensure a proper comparison between measurements and CFD. Concerning CPU costs, a BFM simulation takes less than 24 hours to converge with 64 cores (around 1500 CPU hours). The full-annulus case took 5 revolutions to converge without distortion. Then 4 more revolutions were necessary to converge with distortion. On the whole, the URANS simulation presented here costed around 60000 CPU hours, which gives a ratio of 40 compared to the BFM simulation.

CONCLUSIONS
Body force simulations were carried out on a cooling fan stage. The model used includes compressibility effects, metal blockage effect and proposes a loss term accounting for viscous effects. The results are compared to full-annulus URANS on the one hand, and to experimental data on the other, in order to evaluate the capability of the BFM to predict global performance, as well as the flow physics.
With distortion, BFM and URANS results show a good agreement. Taking URANS as reference, the stage performance are quantitatively rather well predicted without resorting to any extra-calibration. Except for short length scale structures like stator wakes or secondary flows, the distortion transfer across both rows is the same in terms of velocity field, total pressure and total temperature. Some discrepancies exist when comparing with experimental results, which can be partly explained by different upstream flow conditions. More measurements would be needed just downstream of the distortion grid to reproduce exactly the same inlet conditions in the CFD. This limitation will be addressed in the coming year by a new experimental distortion campaign.
On the whole, with the simple distortion pattern studied, the body force simulation captures the main flow behavior, with a CPU cost lower than URANS by a factor of 40. This demonstrates that BFM can be a powerful tool to deal with performance prediction in innovative engine configurations during the first design loops, when the targeted accuracy is not too demanding. The next step in this research project is to evaluate this methodology with a more loaded fan stage, in order to be more representative of propulsive functions, and with a more realistic distortion pattern, closer to current BLI profiles.