Assessment of Steady and Unsteady Full Annulus Simulations Predictivity for a Low-Speed Axial Fan at Load-Controlled Windmill

A steady mixing plane approach is compared with the time-averaged solution of an unsteady full annulus calculation for a conventional fan operating at load-controlledwindmill.The objective is to assess the added value of a complete unsteady calculation compared with a more classical approach, especially concerning the effect of the spatial and temporal periodicity release in such an unusual operation as windmill. Experiment with global steady measurements and rotor radial characterizations was conducted. Numerical analysis demonstrates that windmilling global performances obtained with the time-averaged solution of the unsteady simulation are not far different from the steady case, especially in the rotor. Somedifferences arise in the stator, particularly regarding the velocity field. Temporal periodicity release in this row has clearly a significant effect on the flow unsteady response. A detailed analysis highlights that generic patterns of windmilling flows recorded on a steady approach are also reported on the unsteady case.


Introduction
In a near future, the need for onboard electrical energy will become a subject of great interest for aircraft manufacturers.From this perspective, several studies are currently carried out to investigate any potential energy sources among existing onboard equipment.This is the context of the present study which is focused on the analysis of a low-speed lowdimension fan.
This fan aims at cooling heat exchangers when the aircraft is on the ground.It is located downstream of an air intake.In flight, the mass flow created through the ram air ensures this function; the fan is of no use and windmills freely.In this context, the windmilling operation is no longer considered as a marginal functioning mode only occurring in accidental situations (engine flame-out) but is regarded as a potential energy recovery solution.The windmilling configuration studied in this paper is called "load-controlled windmill" and is located between the freewheeling condition, well documented in the literature, and the locked rotor configuration.
In the literature, several studies on turbofans demonstrated that windmilling operation is characterized by highly negative incidences on both rotor and stator rows which leads to massive separation [1][2][3][4].This phenomenon is more important at load-controlled windmill than at freewindmill.For a better understanding, an illustration of the velocity diagrams obtained for a conventional fan in compressor mode, at freewindmill and at load-controlled windmill, is given in Figure 1.
Windmilling flows complex topology is problematic regarding the numerical convergence and predictivity.Different methods and turbulence models were tested in unpublished previous works either to assess the windmilling performance or to better understand the unsteadiness of the flow in this regime.
Regarding the predictivity, on turbofans at freewindmill, it is generally reported that numerical radial distributions at rotor outlet are in good agreement with experimental data, which is not the case at stator outlet (steady mixing plane simulation [5] and unsteady simulation on 1/12 and 1/10 sector of the full annulus [6,7]).This comes from the fact that losses mechanisms are not well captured.The overall performance is therefore not well predicted.At loaded windmill on small axial fans, numerical radial distributions at rotor outlet are no longer in accordance with the experimental ones, which worsen the discrepancy regarding the global performance (steady mixing plane simulation [8]).This comes from the fact that, in loaded cases, the flow on the rotor is more separated than at freewindmill.In this paper, a steady mixing plane calculation is compared with a time-averaged full annulus calculation at load-controlled windmill, using experimental data as a reference.The impact on performances of the mixing plane approach as well as the periodicity release is investigated.It is important to point out that the windmilling operating point selected for this study is the one for which many convergence difficulties were present in the steady state mixing plane simulation.This is typical of an operating range located between freewindmill and highly loaded windmill.In these cases, the difficult convergence is suspected to come from intrinsic aerodynamic unsteadiness.
To check this assumption, the unsteady calculation presented here was carried out.
Regarding the unsteadiness of windmilling flows, a 2D NLH approach [9] from Dufour et al. [5,10] and an unsteady calculation on 1/12 sector of full annulus from Goto et al. [6] underline the existence of a vortex-shedding phenomenon occurring in the stator.According to Dufour, the rotor wakes that periodically impact the stator blades are responsible for a continuous modification of the stator angle of attack which leads to a more or less intense separation.Goto et al. [6] try to analyze this vortex-shedding on a 3D point of view but because of the complexity of the pattern more investigation is needed.One of the objectives in this paper is to analyze the temporal evolution of the 3D stator separation, in particular by taking into account its radial unsteady behaviour.The unsteady computation presented in this paper is therefore used to analyze the unsteadiness of the flow in order to improve the understanding of windmilling flow topologies and draw generic patterns.A last objective is also to check the trends concerning the rotor topology that were underlined on a steady calculation presented in a previous paper of the present authors [11], by the use of an unsteady simulation.Consequently, the scope of the present paper is to bring light on aerodynamic properties of conventional fans operating in unusual situations, namely, windmilling conditions.The following topics are discussed here: (i) Global and local predictivity of a steady mixing plane simulation and a time-averaged full annulus calculation compared with experimental data (ii) Comparison of secondary flow topologies between the full annulus calculation and the steady mixing plane approach (iii) Analysis of the 3D structure of windmilling separated zones on both rotor and stator on the full annulus calculation

Machine and Methodology
The machine studied is a conventional axial fan of rotor-stator type.As the geometry of the fan is confidential, only some basic properties are presented in Table 1.
A specific representation is used in this paper to account for the global performance of the conventional compressor: the ( ψ, φ) formalism.These parameters, suitable to off-design configurations according to Binder et al. [12] and Turner and Sparkes [13], are defined at the mean quadratic radius: (1)

Experimental Setup
In the Department of Aerodynamics, Energetics and Propulsion (DAEP) of ISAE-SUPAERO, a test facility, called WIL-LOW for WIndmilL LOW-Reynolds, has been implemented to analyze the flow behaviour on low-dimension fans (diameter ≤ 200 mm) in windmilling conditions [14].The test rig is composed of an air vacuum system and an asynchronous electrical engine, the combination of which enables us to reach every possible operating point from compressor surge to locked rotor configuration.This facility allows us to control independently the flow generation by suction and shaft loading.The rotational speed and torque are measured and controlled by a differential transformer torque meter (Magtrol TM 306) coupled with the electrical engine.The test facility is visible in Figure 2. As illustrated, the engine was deported far from the test section to protect the measurements from electromagnetic disturbance.If required, the instrumented section is reversible to invert the flow direction.Fan geometries consist of fixed removable rotors and stators.Measurements between inlet (Plane 1) and outlet (Plane 5) were carried out to assess the global performance.It consists of 3 thermocouples (Type K, basic accuracy 1 ∘ ) and a 6-hole Kent chamber connected to a differential sensor (Rosemount, 0/100 mbars, basic accuracy 0.15% Full Scale (FS)).This enables us to measure the mean total temperature and the mean static pressure.In addition, the outlet (plane 5) is equipped with a fixed five-hole probe used to get the total pressure (Keller, −150/ + 100 mbars, basic accuracy 0.047% FS).In addition, three wall static pressure ports are localized at fixed azimuthal position at rotor inlet, rotor outlet, and stator outlet.The uncertainties were assessed on both global and local performance based on that of the sensors.In addition, a repeatability study was performed for three operating points in compressor, freewindmill, and load-controlled windmill by using Student's -distribution applied over a sequence of several tests for a 95% confidence interval.The result is   presented for the global measurements in Table 2 and for the local measurements of the five-hole probe at rotor outlet in Figure 3. On the whole, the repeatability error was found greater than the uncertainties.The latter are not indicated in the figures of this paper for clarity reasons.
Within the scope of this study, one particular operating point in load-controlled windmilling operation was investigated.It corresponds to a situation where the flow is not fully separated spanwise on both rotor and stator rows.This point is characterized by a flow coefficient of about φ = 1.66, a rotational speed of 1500 rpm, and a mass flow of 0.455 kg⋅s −1 .The presented results were corrected in terms of atmospheric temperature and pressure.

Numerical Setup
Numerical results were carried out with the EURANUS solver of FINE6/Turbo V9.0.2 developed by NUMECA.All computations were performed on an in-house supercomputer, composed of 24 nodes of 16 processors and 8 nodes of 24 processors.The numerical domain for steady mixing plane and full annulus simulations are illustrated in Figure 2. In both cases, the shaft, bell mouth, and rear support were not numerically represented.
Regarding the steady mixing plane approach, a mesh with an overall size of 3.3 M pts was used.The number of flow paths is 109 in the rotor including 29 points in the tip gap and 89 in the stator.The number of cells in the pitch-wise direction is 51 in the rotor and 43 in the stator.Each blade surface (suction and pressure) is discretized with 105 grid points in the rotor and 85 in the stator.The mesh quality parameters are given in Table 3.
The unsteady computation is characterized by a mesh of 60 M pts (rotor: 32 M, stator: 24 M), which corresponds to the repetition of the steady case.According to the aforementioned in-house calculations means, for such a big mesh size the numerical domain had to be reduced (Figure 2).Despite this measure, 240000 CPU hours were still necessary to reach the convergence.A more refined full mesh was not possible for the same reasons.One revolution of the rotor was discretized in 1564 angular positions corresponding to 92 time steps per rotor period and 68 per stator period.The convergence was reached for 10 rotations.
The numerical parameters were set for both computations as follows.The domains are discretized with a multiblock structured mesh generated with AutoGrid 56.The mesh consists of O-grid blocks located around the blades and Hgrid blocks in the passage (O4H mesh).The first cell at wall was set to 5 m.The mean value of the  + parameter at walls is about 0.3 ( + max ≈ 1.5).Global tests data were used to set the total uniform inlet boundary conditions (pressure, temperature) and to impose the mass flow at outlet.The twoequation SST model was used for turbulence closure.The turbulence intensity was calculated based on the Reynolds number (about 3%).The advance in pseudo-time is ensured by a four-stage Runge-Kutta scheme for the steady simulation.Regarding the full annulus case, time discretization is achieved through a dual-time stepping method proposed by Jameson.Iterations in the inner loop are set to 10.The discretization in space is based on a cell-centered finitevolume approach.Convective fluxes are determined by a Jameson second-order centered scheme with added artificial dissipation.

Comparative Analysis
In this section, the predictivity of steady mixing plane and full annulus unsteady RANS simulations are compared in terms of global performances, radial distributions, and local flow topology.Global performances are extracted from the same axial planes as the experiments and mass-averaged.All radial numerical profiles are extracted on experimental axial planes and obtained after an area-weighted azimuthal average.It can be assumed that numerical pitch-wise averaging in the relative frame is equivalent to experimental temporal averaging at a fixed azimuthal position since the computed flow in the rotor is steady [5].In the stator, several radial distributions in the pitch-wise direction are experimentally  measured and averaged over a blade passage.The unsteady calculation has been time-averaged over 6 rotor periods,   .This last parameter is defined as a function of the number of stator blades   as indicated below: 5.1.Global Performance.The global flow coefficient φ is well estimated by both numerical approaches as shown in Table 4.This good agreement comes from the fact that both rotational speed and mass flow are numerically imposed.
On the contrary, a worse prediction is observed on the loading coefficient (6% error) for both simulations.This can come from the estimation of the experimental aerodynamic power by considering the mechanical losses and the use of a mass flow-weighted average for numerical results.The two predicted total-to-total isentropic efficiency is also similar to the measured one.The only difference between the calculations is reported on the stator static pressure drop.Indeed, stator separations for the unsteady case are less extended and more intense than for the mixing plane simulation.With the effective flow area being more important, the axial velocity at stator outlet is smaller and the static pressure is greater than the one recorded on the steady case.In addition, both simulations show the same discrepancy from the experiment in terms of rotor static pressure drop.simulations are equivalent, as expected, and well match the tests data from midspan to shroud.However, a significant discrepancy is recorded in the blade lower half.The lack of predictivity near the hub can be explained by the insufficient development of the numerical boundary layer.At rotor outlet, no differences are reported between the mixing plane approach and the time-averaged unsteady case.On the whole, the axial velocity is rather well predicted along the span.This means that the vein blockage (the ineffective flow area out of the geometric area [15]) is well estimated by the simulations.At stator outlet, greater discrepancies are recorded between both methods.The time-averaged solution of the full annulus case is reported to be in better agreement with the experimental data.This result is not surprising since this method takes into account the rotor wakes effect on the stator topology and is therefore expected to be more representative of the actual flow reality in this row.

Radial
The radial distributions of the rotor loading coefficient and pressure ratio are considered in Figure 5.The work exchange does not evolve monotonously from hub to shroud at load-controlled windmill contrary to classical fan operation.The general shape is caught by both simulations, in particular the slope inversion.However, a similar overestimation is recorded from hub to 30% of the span as well as an underestimation from this height to shroud.The latter indicates a bad behaviour of the methods in the stalled areas.The bad prediction of the tangential velocity is believed to be responsible for the bad agreement of the work exchange with the tests data.Regarding the losses, the numerical prediction is good from hub to shroud.Surprisingly, the steady simulation seems slightly in better agreement with the tests data.More accurate results would have been expected from the full annulus case since more flow physics is taken into account.
The predictivity of the numerical methods is now checked in the stator.Radial distributions of incidence and pressure ratio are given in Figure 6.At inlet, the flow is well predicted by both methods from hub to 30% of the blade span.From this height, great discrepancies are recorded.The mixing plane approach appears to give better results than the unsteady case.This can be explained by the location of the axial plane from which the radial distributions are extracted: the rotor/stator axial distance is too small to put more than one five-hole probe.The extraction plane for stator inlet is the same as that of rotor outlet and is numerically located before the mixing plane.Regarding the stator losses, both simulations far underestimate the pressure ratio.

International Journal of Rotating Machinery
On Figure 7, an illustration of the vein blockage, provoked by the stator massive separation, is given for all methods.Significant differences arise between the numerical approaches when the azimuthal variation is considered.This result clearly discredits the mixing plane approach and underlines the relevance of the unsteady case to predict the flow physics in the stator, despite a slight overestimation of the maximum velocities at the hub.
Finally, in the rotor, the flow topology was rather well captured by both methods, except the work distribution in the stalled regions and near the hub.This point leads to suspect the turbulence closure.In-house studies were carried out to investigate the influence of a large number of turbulence models.However, no significant improvement was witnessed, although the SST and Spalart-Allmaras models were reported to be the best ones.On the whole, the permeability was relatively well predicted throughout the stage for both methods.In the stator where most of the differences arise, the unsteady approach seemed more relevant than the steady one, in particular when the azimuthal direction is considered.However, in both cases, total pressure losses in the stator are badly assessed, underlining the potential irrelevance of RANS methods to capture the physics of separated flows at windmill.A first enhancement of the simulations would be to use a laminar-turbulent transition model to account for Reynolds number effects.This solution was tested by using the model proposed by Langtry et al. [16,17] on a more refined mesh (14 M pts) in order to respect the method requirements.No differences were witnessed in the rotor.In the stator, no improvement was reached: a slight difference in pressure loss is the only point worth being reported.The next step would be to try LES calculations to capture more flow physics.

Local Flow
Topology.The objective of this part is to check if the generic trends, typical of windmilling flows, that were recorded with a mixing plane approach on a previous paper of the present authors ( [11]) are observed on the time-averaged unsteady computation as well.
Figure 8 shows static pressure contours with streamlines near the rotor hub.Both numerical results are quite similar.Curvature effects inside the hub boundary layer generate a crosstream directed inversely compared with the rotor rotation.
Likewise, the radial velocity contours given at rotor outlet on Figure 9 emphasizes that the same results are obtained on both simulations.In addition, each passage shows exactly the same pattern.
The aforementioned results lead finally to the same 3D structure for the rotor separation, as visible on Figure 10.The flow is separated on the entire leading edge but propagates towards the shroud due to positives values of radial velocities and decrease in solidity.On the whole, all the results regarding the rotor topology are confirmed by the unsteady calculation.
Figure 11 presents the contour of axial velocity at stator outlet for both calculations.Extreme values are more important for the steady mixing plane simulation than for the timeaveraged solution.The flow topology is smoothed in this latter case due to a homogenization enabled by temporal and spatial constraints release.
What emerges from the previous numerical results is the fact that generally speaking no significant differences exist between the steady mixing plane calculation and the timeaveraged unsteady case in the rotor.This result comes from the fact that stators only affect the system permeability and do not modify the rotor topology [8].On the stator, discrepancies are reported on extreme values and flow topologies.This last point will be investigated in the next section.

Analysis of Stator Separation at Windmill
In this section, the objective is to analyze the unsteady 3D structure of the stator separation, which is presented on Figure 12 at a given time step and for 8 successive passages.Conveniently, for a better understanding, the numbering of these passages is indicated on all the figures presented in this section.The flow separation is visualized thanks to a zero axial velocity isosurface.The main differences are found at the trailing edge from one passage to another, where the separation structure highly evolves in the three directions.In the tip region, the separation expands axially and radially from passages 1 to 4. In this last passage, the separation is truncated which means that the azimuthal width of the separation is spread on at least two passages.Radially, a vortex      International Journal of Rotating Machinery 9 seems to emerge from the shroud and then migrates towards the hub while advected to the outlet.Figure 13 describes the evolution in time of entropy contours on blade-to-blade views taken at different time steps during one rotor period.One rotor and one stator blades are filled to follow the rotor rotation.The corresponding thetaaveraged radial distributions of the stator angle of attack at the given times are presented with black thick lines (other time steps are in grey).The red line corresponds to the time-averaged profile during the rotor period.This figure confirms the phase locking of the stator separation on the rotor passage which was first reported by Dufour et al. [5,10].A structure, identified at stator outlet, seems to evolve in the opposite direction of the rotation.The figure illustrates that the vein blockage is typically more important in the stator than in the rotor.In windmilling operation, secondary flows are responsible for a tipward motion of the rotor separation which widens the wakes impacting the stator [11] near the shroud.This may explain why most of the absolute flow angle variations are located in the tip region.
For a complete understanding of the 3D structure of the stator separation, its radial evolution is now investigated.Entropy contours and zero axial velocity isolines are given at stator outlet in Figure 14 over a rotor period.In addition, the radial velocity distributions are presented for the same time steps and were extracted for a fixed azimuthal angle (black thick lines).The red line corresponds to the time-averaged profile during the rotor period.As can be seen, there is a great heterogeneity of entropy although the highest values are located near the shroud where separation is the most massive.The sign of the radial velocity seen by the passage number 3 changes several times during a rotor period.This underlines that the fluid in this passage is alternatively moving towards the hub and towards the shroud, attesting of a radial unsteadiness of the stator separation.This behaviour is a consequence of the periodic confinement induced by the displacement of the equally spaced high entropy areas (spatial mode 3 turning in the opposite direction of the rotation).To the authors' knowledge, so far this unsteady radial behaviour of the separation has never been reported in the literature concerning windmilling operation.
Numerical results revealed a highly unsteady 3D separation occurring in the stator.Thorough understanding of this topology remains to be achieved but separated zones are clearly locked on the phase of the rotor.This phenomenon is particularly important in the shroud region where the rotor wake is the widest due to the effects of secondary flows.The URANS calculation revealed an interesting topology as far as the windmilling flow understanding is concerned.However, such a costly approach has no added value concerning global performance estimation or understanding of rotor local topology.Experimental preliminary results also show some indication of complex unsteady behaviour but the complete experimental campaign is still ongoing.

Conclusion
The first objective of the present paper is to compare global and local performance between a mixing plane steady calculation and the time-averaged solution of a full annulus unsteady calculation, for an axial fan operating at a selected operating point at load-controlled windmill.The second objective is to analyze the stator separation during one rotor period to better understand its 3D structure.The following findings emerged from the previous discussions: (i) Little discrepancies are visible in the rotor in terms of global and radial performance between the mixing plane calculation and the time-averaged solution of the unsteady case.(ii) Significant discrepancies are reported between the two approaches in the stator where the full annulus approach seems to be more representative of the flow, especially in the pitch-wise direction.(iii) Regarding the rotor topology, little discrepancies are recorded between numerical and experimental results from 30% of the blade span regarding the work prediction.This is imputed to the 3D triangular shape of the rotor separation which barely affects the region near the hub.The bad predictions of the intensity of the stalled region as well as the position of the reattachment line are suspected to be involved.In addition, in the stator, the numerical prediction of the losses is far from the experimental data.All in all, the relevance of RANS methods for fans operating at load-controlled windmill can certainly be questioned.(iv) At windmill, the structure of the stator separation is highly tridimensional and locked to the phase of the rotor.

Figure 1 :
Figure 1: Velocity diagrams from compressor to load-controlled windmill.

Figure 2 :
Figure 2: Illustration of the deep windmilling test facility (willow) and the numerical domains.
Local radial distributions were recorded to analyze the rotor flow topology.Two directional five-hole probes are, respectively, located around 7 mm upstream the leading edge (Plane 2) and 5 mm downstream the trailing edge (Plane 3) of the rotor (Keller, −150/ + 100 mbars, basic accuracy 0.047% FS, Rosemount −50/50 mbars, basic accuracy 0.15% FS).Radial probing of velocity and pressure fields as well as flow angles is conducted from hub to shroud through around 35 positions.A third directional five-hole probe is used around 27.5 mm downstream the trailing edge of the stator (Plane 4) (Keller, −150/ + 100 mbars, basic accuracy 0.047% FS, Rosemount −50/50 mbars, basic accuracy 0.15% FS).In this case, both azimuthal and radial probings are enabled over 45 ∘ in the pitch-wise direction with a discretization of 1 ∘ .This leads to around 14 measured distributions within a stator blade passage for about 15 different radii.The azimuthal movement of the five-hole probe is motorized and ensured by a ring gear mounted on the casing.A Newport high performance rotation stages device ensures the rotation of this gear over an angular sector of 45 ∘ , corresponding to about 3 blade passages, depending on the considered rotor.The angular resolution of the engine is 0.0025 ∘ with a mechanical clearance of about 0.01 ∘ .Data acquisition is achieved by a composite National Instruments (NI) cDAQ-9172 device connected to 3 kinds of NI modules that are used for current (NI 9215), voltage (NI 9203), and temperature (NI 9211) recordings.The measurements are acquired on 5000 samples with a sampling frequency of 20 kHz.Rotation and translation for the five-hole probes are carried out by a 2-axis NI motion device (NI PCI-7344).

Figure 3 :
Figure 3: 95% confidence intervals for local measurements at rotor outlet in compressor mode, freewindmill, and load-controlled windmill.

Figure 4 :
Figure 4: Experimental and numerical rotor and stator axial velocities.

Figure 5 :
Figure 5: Experimental and numerical rotor work and losses.

Figure 6 :
Figure 6: Experimental and numerical stator incidence and losses.

Figure 7 :
Figure 7: Experimental and numerical reduced axial velocity at stator outlet.

Figure 8 :
Figure 8: Numerical rotor static pressure contours and streamlines near the hub.

Figure 12 :
Figure 12: Numerical zero axial velocity isosurfaces on stator blades at a given time step for 8 successive passages.

Figure 14 :
Figure 14: Numerical entropy contours and zero axial velocity isoline at stator outlet.Radial velocity distributions at a given time step (mean profile in red).

Table 1 :
Basic properties of the conventional compressor.

Table 2 :
95% confidence intervals for global measurements in compressor mode, free windmill, and load-controlled windmill.

Table 4 :
Synthesis of experimental and numerical global performance.