Prediction of Film Thickness of an Aero-Engine Bearing Chamber using Coupled VOF and Thin Film Model

In the present work, a coupled volume-of-fluid (VOF) model with Eulerian thin-film model (ETFM) approach is used to predict the film thickness in an aero-engine bearing chamber. Numerical studies are conducted for a wide range of shaft speeds with lubricant and air flow rates of 100 l/hr and 10 g/s respectively, at a scavenge ratio of 4 on a simplified bearing chamber test rig. Air-flow analysis inside the bearing chamber is also assessed. Primary and secondary airflow predictions are found to be in good agreement with the experimental results. The coupled ETFM+VOF approach is found to be sensitive enough to capture the qualitative trend of oil film formation and distribution over the chamber wall. Oil collection near the sump at a low shaft speed and a rotating oil film at a higher shaft speed are well captured.


INTRODUCTION
The bearing chamber of an aero-engine is a cylindrical enclosure of the turbine shaft and bearings. The lubricating oil supplied to the bearings is dispersed due to the angular velocity of the shaft and collected at the inner wall of the bearing chamber. Glahn et al. [1] identified bearings as a primary source of the formation of oil droplets where lubricating oil is atomized by the roller elements and the cage before ejecting into the bearing chamber. The oil shed from the bearings may appear in the form of droplets, mist, ligaments, thin-film, or combinations of these, Chandra et al. [2], as depicted in Fig. 1. Pressurized air is supplied to labyrinths seals to isolate bearing chamber from the gaseous environment of the engine. The interaction of air with this shed oil presents a complex two-phase flow environment. An oil film driven by gravity and shearing action of air is observed on the rotating shaft and the stationary walls of the chamber. Kurz and Bauer [3] conducted an experimental study to investigate the influence of shaft speed on the oil film formation. They identified two distinct oil flow regimes. At lower shaft speeds, the oil flow is mainly driven by gravity and a very thin oil film was observed on the chamber wall. At higher shaft speeds, a rotating oil film was observed with the shearing action of air influencing the oil film. The momentum flux of air and momentum transfer to the air by the rotating shaft was used as a parameter to determine the critical shaft speed for bearing chambers of varying diameter and length to demark the change in the regime. By using the later parameter, Kurz and Bauer [3] successfully scaled the critical shaft speed with oil flow rate for the bearing chambers of different length but of the same outer diameter. However, the correlation did not work for the bearing chambers of different outer diameters. The experimental studies of Gorse et al. [4,5] considering airflow without oil, revealed the influence of air flow rate and rotational speed on the secondary airflow. The airflow in the azimuthal direction (ϕ) is classified as primary airflow and airflow in axial-radial plane is considered as secondary airflow structure. A strong influence of rotational shaft speed and airflow rate was reported on the secondary air vortices. The flow behavior is classified either as rotational speed driven mode (RSDM) or sealing air driven mode (SADM). Above a velocity ratio (Ush/ Uax) of 4.2, rotational speed was found to drive the flow field and below 4.2 sealing air influenced the flow field.
Despite the experimental evidence of the influence of air on the development and spreading of the oil film, interaction of air with the oil is not well-understood inside a bearing chamber. It is difficult to quantify the interaction of air with the oil by experimental studies alone due to the complexity of the flow and limitations of measurement techniques. Hence, a reliable and efficient computational model based on computational fluid dynamic (CFD) is much desired to enhance the existing understanding of thin-film hydrodynamics inside the bearing chamber. Building on the previous work of Kakimpa et al. [6], a volume-of-fluid (VOF) model coupled with Eulerian thin-film model (ETFM) is implemented to predict the film thickness in an aero-engine bearing chamber. A multi-scale lubricant film is observed in the bearing chamber of an aero-engine. Utilizing the VOF approach alone to capture the formation of the thin-film and its interaction with the air-flows induced due to the rotating shaft is too computationally expensive as it requires relatively refined grids to capture the thin oil film. To cope with this difficulty, and to present a fast predictive tool, VOF coupled with ETFM was implemented by Kakimpa et al. [6] for rimming flows. A depth-averaged ETFM model is used in the region where oil film thickness is smaller than the affordable grid size required for VOF. A traditional VOF approach is retained in the domain where grid resolution is sufficient for the VOF studies.
In the study of Kakimpa et al. [6], a qualitative comparison of numerical results of coupled ETFM+VOF model and VOF were presented. Qualitative trends could not reveal much about the performance of the coupled model. Hence, the authors recommended to validate this model against experimental results. In the present study, the coupled ETFM+VOF approach is utilized to predict film thickness on a simplified bearing chamber. A quantitative comparison is required to understand the film-thickness prediction capability of the coupled ETFM+VOF model. Present work emphasized on accessing the coupled model capability to predict film thickness inside a simplified bearing chamber. In the coupled approach, ETFM and VOF are coupled through mass and momentum conserving source terms. In order to define the transition from ETFM to VOF or viceversa, the total liquid volume fraction in each cell is evaluated and depending on local grid resolution either the ETFM or VOF model is used. In this study, the commercial CFD code ANSYS/FLUENT, based on finite volume method (FVM) is used for numerical analysis. The coupling of ETFM and VOF is achieved by incorporating a user defined function (UDF).
The model validation for an aero-engine bearing chamber application is included in this study. The numerical results are compared with the experimental results of Kurz et al. [7] on a simplified bearing chamber for shaft speed ranging from 5,000 rpm to 15,000 rpm for the lubricant and air flow rates of 100 l/hr and 10 g/s respectively, at a scavenging ratio of 4. Shaft geometry of the bearing chamber is extracted from the experimental work of Kurz et al. [7]. A shaft of diameter, D=128mm rotates through the whole length of the chamber. The diameter and axial length of the bearing chamber are 1.73D and 0.52D, respectively.
A set of roller bearings are mounted on the one end of the shaft and labyrinth air seal on the other end of the shaft. The lubricating oil from the bearing enters into the chamber from the bearing side and air enters from the labyrinth air seal. Krug [8] reported a negligible effect of the bearing geometry on the oilfilm velocity distribution on the wall in case of roller bearings without air cross-flows. Hence, in the computational model, the oil inlet was modelled as a plain surface rather than modelling roller bearing, similar to the work of Peduto et al. [9], Krug et al. [10] and Bristot et al. [11]. Circular vent and scavenge pipes of diameter 10mm are located at the top and the bottom of the chamber along the width-wise centerline. In the case of an aeroengine bearing chamber, a known mass of the lubricant based on the pre-decided scavenging ratio is drawn from the sump. In order to model this condition numerically, a scavenge pump model approach based on the work of Robinson et al. [12] is adopted. Thus, the computational model with the key dimensions from the experimental rig of Kurz et al. [7] and with the above mentioned simplifications is shown in Fig. 1.

Coupling of Eulerian thin-film model with VOF
In the present study a volume-of-fluid (VOF) model coupled with Eulerian thin-film model (ETFM) developed by Kakimpa et al. [6] is implemented to predict the film thickness in an aeroengine bearing chamber. In this study, a finite volume method based commercial CFD code ANSYS/FLUENT is used. ANSYS/FLUENT has the capability to solve VOF and ETFM individually. Coupling of these two models is achieved by implementing a user defined function (UDF). The transition from one model to the other is triggered based on a local mesh resolution criteria. The grid resolution requirement of ETFM, VOF and coupled ETFM+ VOF method is depicted in Fig. 3. Fig. 3(c) also demonstrated the grid resolution criterion adopted in this study.
The local mesh resolution is assessed on the basis of total volume fraction (αl,total = αl,ETFM + αl,VOF) of the liquid resident in a given cell. Following the work of Kakimpa et al. [6], a cut-off volume criterion (αcut-off = 0.5) is adopted. When the volume fraction of a given cell reaches the cut-off volume computations are switched to a film-resolving VOF method and vice-versa. Film mass and momentum are transferred from one model to another according to the cut-off volume criterion. The flow governing equations viz. continuity and momentum for the individual model can be found in references [6] & [13]. Moreover, the continuity and momentum equation for the coupled model are given by equations (1-2) for VOF region and (3-4) for ETFM region.
where, α is volume fraction, u is velocity, m is equivalent density of all phases, p is pressure, g is acceleration due to gravity, h is film height, q is flow rate, σ is surface tension, k is interface normal curvature, αVOF is contribution of VOF in (1) αl,total > 0.5 volume fraction of oil, V cell, is cell volume, film is density of film material i.e. oil, Af is cell face area.

Turbulence Modelling and Free Surface Tracking
The flow inside a bearing chamber is turbulent. In the present study, Reynolds Averaged Navier-Stokes (RANS) is solved to account the effect of turbulence. Based on the prior experience of Nottingham group [11], Shear Stress Transport (SST) k-ω turbulence model is used to achieve closure of equation. Turbulent quantities such as k and ω are used to calculate velocity scale and length scale which are √ and √ respectively.
In the case of the VOF approach, fluids filling a cell are treated as a mixture and as a consequence the speed of all the phases is the same in a given cell. This situation leads to a lack in the definition of the interface which results in a non-physical turbulent momentum transfer between the liquid and gaseous phase. Egorov [14] implemented a source term in ω equation of k-ω model to introduce solid wall like damping of turbulence at the interface of both the phases. This source term is given by equation (9) . ∆ . ∆ where, A is interface area density, subscript i presents ith phase, the closure coefficient of the destruction term, ρ is density, B is a scaling coefficient known as turbulence damping factor, μ is viscosity and Δ is the near-wall cell height. Egorov [14] recommended a value of 10 for turbulence damping factor (B) for grids having Y + =1 at the interface and a higher value for the coarser mesh. Tkaczyk [15] and Adeniyi et al. [16] utilized this approach for the application of bearing chambers. Tkaczyk [15] recommended a value of 100 for the turbulence damping factor based on the numerical studies on stratified flows for a relatively coarser mesh at the interface. This value of 100 for the turbulence damping factor is therefore adopted in the present study.

BOUNDARY CONDITIONS AND MESH
The boundary conditions imposed on the computational domain are depicted in Fig. 4. Air enters in to the domain from an annular surface of height 0.7 mm which represents labyrinth air seal and oil enters from a representative surface for the bearing following the work of Peduto et al. [9], Krug et al. [10] and Bristot et al. [11]. A mass flow inlet boundary condition with a swirl angle, γ= 50° is specified at the inlet of these surfaces based on the guidelines of Bristot et al. [11] and Peduto et al. [9]. The shaft is considered as a rotating wall with the rotation direction of shaft clockwise when seen from the bearing side. A pressure outlet boundary is specified at the outlet of vent whereas a negative velocity boundary condition is imposed at the sump outlet. The value of velocity magnitude is set on the basis of scavenging ratio. Non-slip wall boundary conditions are imposed on all the other surfaces. Eulerian film wall conditions are specified at the chamber walls. Boundary conditions imposed on the computational domain are depicted in Fig. 4.
In the present study, non-uniform structured grids are used to discretize the computational domain. Grids with 1,817,260 cells were generated using ICEM/CFD, a pre-processor to ANSYS/FLUENT.

OPERATING CONDITIONS AND FLUID PROPERTIES
The operating conditions and fluid properties considered in the present numerical work are taken from the experimental study of Kurz et al. [7]. The shaft speed is varied from 5000 RPM to 15000 RPM. Air and oil are considered as working fluids. Fluid properties corresponding to the operating conditions are given in Table-1.
(a) (b) First order implicit scheme is used for temporal discretization and second-order upwind interpolation scheme is used to spatially discretize the governing equations. Coupled algorithm is used for pressure-velocity coupling.

Air Flow Analysis inside the Bearing Chamber
In order to assess the single phase airflow behavior of the investigated test rig, numerical simulations are carried out considering airflow only without oil. The experimental study which is considered as a benchmark [7] for the present numerical investigation, covers only scavenging efficiency and film thickness distribution, not airflow information. Hence, for the comparison of airflow inside the bearing chamber, a numerical study is also conducted on the experimental rig of Gorse et al. [4]. The test rig investigated by Kurz et al. [7] has identical dimensions other than the axial length of bearing chamber which was 1.4 times that of Gorse et al. [4]. Air flow rate and shaft revolution per minute considered in the simulations were 12.7 g/s and 9700 RPM for the test rig of Gorse et al. [4] and 10 g/s and 10,000 RPM for the test rig of Kurz et al. [7]. These operating conditions were selected so that qualitative comparisons could be made. The radial profile of nondimensional tangential velocity is plotted in Fig. 5 for the investigated test rigs, showing the experimental data from [4] as well as the two simulated cases. It can be seen from this figure that the radial variation of non-dimensional tangential velocity follows the trend as reported in the experimental study of Gorse et al. [4]. The numerical results predict the radial distribution of nondimensional tangential velocity with a maximum 10% difference in the tangential velocity close to the shaft compared to the experimental results. This difference decreases as the radial distance increases i.e. towards the chamber wall.
It can also be observed from Fig. 5 that the radial distribution of tangential velocity for the case simulating the test rig of Kurz et al. [7] follows a similar trend, which no airflow experimental data exists for this case. The difference in the magnitude can be attributed to the variation of the geometric and operating condition in these two studies. Gorse et al. [4] characterized the air-flow in tangential/ azimuthal direction as primary air-flow and the airflow in the radial-axial plane as the secondary air-flow.
The secondary air vortices are responsible for the spreading of the oil onto the walls of chamber. Gorse et al. [4] reported that the secondary vortices are either driven in rotational speed mode (RSDM) or by sealing air mode (SADM). In RSDM two counter rotating vortices are observed with a separating streamline close to the mid-span of the axial length of the chamber. The vortex towards the sealing air side keeps on increasing with the increasing flow rate of air. In SADM, the vortex on the bearing side disappears and a single vortex is observed. A transition between these two modes was reported for ma = 12.7 g/s and n = 9700 rpm i.e. the case considered for validation. The velocity vectors in the r-z plane are plotted in Fig. 7. The experimental results of Gorse et al. [4] are shown in Fig. 7 (a) and the present numerical results are shown in Fig. 7(b). It is evident from the figure that the present numerical study captured the transition mode accurately. The location and size of the  Despite the fact that the bearing chamber length is different in the test rig of Gorse et al. [4] and Kurz et al. [7] and the operating conditions are also slightly offset, the mode of secondary vortices can be analyzed based on the velocity ratio (R) criterion proposed by Gorse et al. [4] defined in Eq. (7). If R > 4.2, airflow is driven by the rotational speed of the shaft otherwise airflow is driven by the sealing air.

= ,
On the test rig of Kurz et al. [7] for ma = 10 g/s and n = 10,000 RPM, the velocity ratio R is 5.75 and hence two counter rotating vortices are expected corresponding to RSDM. The velocity vectors in the r-z plane obtained from the present numerical study on the test rig of Kurz et al. [7] are shown in Fig. 8 for the qualitative comparison. Fig. 8 Velocity vectors in r-z plane obtained from numerical analysis for ma = 10 g/s and n = 10,000 RPM on the test rig of Kurz et al. [7].
Two almost identical but counter rotating secondary air vortices can be identified in Fig. 8 which confirm the rotational speed driven mode for the given operating conditions. Apart from the qualitative comparison of the secondary air vortices, mean tangential velocity, of the core air-flow is assessed for the test rig of Kurz et al. [7]. The only available airflow data in the reference case of Kurz et al. [7] was the mean tangential velocity, of the core air-flow. Hence, the numerical results of mean tangential velocity are compared at 5,000 rpm, 10,000 rpm and 15,000 rpm. The deviation in the numerical results was 6% at 5,000 rpm, 2.2 % at 10,000 rpm and 1.3 % at 15,000 rpm. Hence, it can be concluded that the present numerical study is capable of predicting airflow inside a bearing chamber within an acceptable accuracy.

Film Thickness
Film thickness was measured at eight locations along the periphery of the bearing chamber in the experimental study of Kurz et al. [7] which is the benchmark case for validating film thickness. The diameter of the film thickness measuring sensors was 10mm in the experimental study. Therefore, in order to calculate a sensor area-averaged film thickness at a particular location, nine computational probes were placed in an equivalent probe area in the numerical model for each experimental sensor location. The experimental location of sensors and a scheme for a computational probe is shown in Fig. 9. A comparison of the film thickness obtained from the present numerical study carried out at shaft speed of 5,000 rpm, air-flow rate of 10 g/s and oil flow rate of 25.8 g/s is shown in Fig. 10, this also shows the result of a VOF analysis carried out alongside the coupled ETFM+VOF and results from a previous numerical study by Krug et al. [10]. The numerical investigation of Krug et al. [10] utilized VOF approach alone to predict film thickness; two grids with cell size 621,500 (M1) and 1,700,000 (M2) were investigated. The computational grids used in the present study is close to the second mesh (M2) of Krug et al. [10].
It can be observed from Fig. 9 that the present numerical results obtained from coupled ETFM+VOF approach is in better agreement with the experimental results compared to the numerical results of Krug et al. [10] obtained from VOF method.
It can be observed from Fig. 10 that the numerical results obtained from the coupled ETFM+VOF approach are underpredicting the film thickness as compared to the experimental results of Kurz et al. [7]. It can also be noted that the deviation in the numerical results at co-current chamber side (180° < ϕ < 360°) is lower as compared to the counter-current chamber side (0° < ϕ < 180°). At the lower shaft speed (5000 RPM) gravitational force dominates over the air shearing force. Consequently, the lubricating oil collects in the region close to the sump. Although, the shearing force exerted by the air is not sufficient to carry a bulk of oil with itself at counter-current chamber side where gravitational force and shear force acts in opposite direction, yet air displaces the oil. As a result a thick layer of oil is observed at location, ϕ = 135° in the experimental study. This phenomenon is not well captured in the computational studies. Fig. 10 Comparison of predicted film thickness with the experimental measurement of Kurz et al. [7] and numerical results of Krug et al. [10].
Apart from the coupled approach, a numerical study was also carried out with the VOF alone. For the case of VOF studies, the mesh was refined near the outer wall of the chamber to sufficiently resolve the oil film. The centroid of the first cell was placed 15 μm from the chamber wall. Based on the recommendation of Bristot [17] six nodes were kept in first 200 μm spacing from the chamber wall. The final mesh used for VOF studies had 3,206,572 cells. It can be observed from Fig. 9 that film thickness is over-predicted by VOF approach. This might be due to the false momentum transfer between the phases as reported by Krug et al. [10]. Turbulence damping factor is the controlling parameter to adjust the momentum transfer between the phases. In the present study, a value of 100 is taken for turbulence damping factor based on the numerical study of Tkaczyk [15] on stratified flows in a channel for coarser grids at interface. Bristot [17] reported a huge variation in the numerical results with a turbulence damping factor of 10 and 100 specifically at counter-current chamber side. Hence, further investigation is required to obtain an optimum value of turbulence damping factor for bearing chamber applications.
For the VOF only case computational time was increase by three times compared to the coupled ETFM+VOF approach. This indicates the importance of the coupled approach in developing an efficient method for predicting film-thickness inside a bearing chamber. Although, the results obtained from the coupled approach are not in complete agreement with the experimental results they do show a closer agreement than the VOF results reported in the literature for a similar mesh size. The coupled ETFM+VOF approach is still in development and further improvement in the modelling approach including a better understanding of the effects of switching criterion from ETFM to VOF or vice versa and local mesh resolution based on the switching criterion is recommended.
A qualitative comparison of oil distribution over bearing chamber wall is shown in Fig. 11 to Fig. 13 for shaft speeds of 5000rpm, 10,000rpm and 15,000 rpm respectively. It is evident from the iso-surface of oil volume fraction, αi,VOF = 0.5 that oil collected towards the sump in the region, 150° < ϕ < 180° for a shaft speed of 5,000 rpm, shown in Fig. 11 (a). This trend matches experimental observation [7] of an increase in oil film thickness towards the scavenge pipe. It can also be observed that in the region 150° < ϕ < 180°, the contribution of the Eulerian thin film model in the measurement of film-thickness is negligible. This indicates that the specified region is resolved by VOF approach. The contours of oil volume fraction shown in Fig. 11 (b) also confirms that the oil volume fraction in the region 150° < ϕ < 180° is higher than 0.5 which was the criterion adopted for switching from ETFM approach to VOF or vice-versa. It shows that the present coupling approach of ETFM and VOF is functioning as per the expectations. At the higher shaft speed i.e. 10000 rpm and 15000 rpm a rotating oil film is observed. The

Region resolved by VOF
Region resolved by ETFM contours of film thickness along with the iso-surface of oil volume fraction, α i,VOF = 0.5 (refer Fig. 12a and Fig.13a) reveals that the oil is not only near the sump but is distributed on the walls of the bearing chamber at the higher shaft speeds.

CONCLUSION AND RECOMMENDATIONS
In the present work, a coupled volume-of-fluid (VOF) model with Eulerian thin-film model (ETFM) approach is used to predict the film thickness in an aero-engine bearing chamber. Numerical studies are conducted for a wide range of shaft speed for the lubricant and air flow rates of 100 l/hr and 10 g/s respectively, at a scavenging ratio of 4 on a simplified bearing chamber test rig.
The numerical results revealed that the coupled ETFM+VOF approach is sensitive enough to capture the qualitative trend of oil film formation and distribution over the chamber wall. Oil collection near the sump at lower shaft speed and a rotating oil film at higher shaft speed is well captured. Under-prediction of the oil film-thickness over the chamber wall as compared to the experimental results, was noticeable. The coupled ETFM+VOF approach is still in development and further improvement in the modelling approach include a better understanding of the effects of switching criterion from ETFM to VOF or vice versa and local mesh resolution based on the switching criterion is recommended. The turbulence damping source term is also an important factor which governs the momentum exchange between the liquid and gas phases. A detailed study to understand this source term and scaling of turbulence damping factor for bearing chamber application is also required.