Tutorial - Thermo-Hydro-Dynamic (THD) Bearings#
This tutorial covers the thermo-hydro-dynamic (THD) bearing classes available in ROSS. These classes compute dynamic coefficients by solving the Reynolds equation (pressure field) and the energy equation (temperature field) for the lubricant film. They model the physics of the oil film in detail and are the preferred choice for industrial-grade analysis.
Class |
Bearing type |
Model |
|---|---|---|
|
Fixed-geometry journal bearing |
THD, multiple geometries |
|
Tilting-pad journal bearing |
THD, adiabatic or full thermal |
|
Tilting-pad thrust bearing |
THD, axial load capacity |
|
Squeeze film damper |
Analytical short-bearing, 3 geometries |
All four classes inherit from BearingElement and can be placed in a rs.Rotor assembly exactly like the simpler classes described in the general bearings tutorial.
Section 1: PlainJournal#
This tutorial demonstrates how to use the PlainJournal class from ROSS to simulate the thermohydrodynamic behavior of fixed-geometry journal bearings (cylindrical, lemon-bore, elliptical, and lobed) and to determine their dynamic coefficients for different operation conditions.
1.1: Example of application#
First, to run the test case provided as an example in the documentation of the PlainJournal class, use:
from ross.bearings.plain_journal import PlainJournal
from ross.units import Q_
import plotly.io as pio
pio.renderers.default = "notebook"
plainjournal = PlainJournal(
n=3, # Nodal location of the bearing in the FE model of the rotor
axial_length=0.263144, # Axial length of the bearing
journal_radius=0.2, # Journal radius
radial_clearance=1.95e-4, # Bearing radial clearance
elements_circumferential=21, # Number of elements in the circumferential direction per pad
elements_axial=11, # Number of elements in the axial direction per pad
n_pad=2, # Total number of pads
pad_arc_length=176, # Angular length of the pads
preload=0, # Preload factor
geometry="circular", # Bearing geometry type
reference_temperature=50, # Oil supply temperature
frequency=Q_([900, 1200], "RPM"), # Journal rotational speed
fxs_load=0, # External load applied to the shaft in the horizontal direction
fys_load=-112814.91, # External load applied to the shaft in the vertical direction
lubricant="ISOVG32", # Lubricant type
sommerfeld_type=2, # Equation used to calculate the Sommerfeld number
initial_guess=[
0.1,
-0.1,
], # Initial guess to find the equilibrium position of the shaft using an optimization process
method="perturbation", # Method used to calculate the bearing dynamic coefficients
operating_type="flooded", # Lubrication condition in the bearing
groove_factor=[
0.52,
0.48,
], # Mixing factor between the inlet oil and recirculating oil
oil_supply_pressure=0, # Oil supply pressure
)
/home/gsabino/Documents/dev/rosstest/fernandarossi/venvross/lib/python3.13/site-packages/tqdm/auto.py:21: TqdmWarning: IProgress not found. Please update jupyter and ipywidgets. See https://ipywidgets.readthedocs.io/en/stable/user_install.html
from .autonotebook import tqdm as notebook_tqdm
/home/gsabino/Documents/dev/rosstest/fernandarossi/ross/ross/bearings/plain_journal.py:739: OptimizeWarning: Covariance of the parameters could not be estimated
popt, pcov = curve_fit(viscosity, xdata, ydata, p0=(6.0, 1.0))
As shown in the example above, when method="perturbation" is selected, small perturbations are applied to the rotor around its equilibrium position, with the stiffness and damping coefficients calculated using finite-difference approximations of the derivatives of the bearing forces with respect to rotor displacements and velocities.
On the other hand, when method="lund" is selected, the solution of the Reynolds equation is perturbed, with the stiffness and damping coefficients determined by integrating the dynamic pressure fields over the bearing surface.
In addition, to help the user set up a simulation case using the PlainJournal class, the table below specifies the required parameter for each operating_type option:
Option |
Required parameter |
|---|---|
|
|
|
|
1.2: Bearing representation#
To display a schematic representation of the bearing geometry (number of pads, groove locations), including the direction of the external load applied on the shaft, use:
fig = plainjournal.plot_bearing_representation()
fig.show()
Note: The way the PlainJournal class was built, the first groove is always located at the positive horizontal axis.
1.3: Optimization convergence#
To display in table and graph formats the magnitude of the residual force vector acting on the shaft throughout the iterations of the optimization process for each analyzed speed, use:
plainjournal.show_optimization_convergence(show_plots=True)
=========================================================
OPTIMIZATION CONVERGENCE - 900.0 RPM
=========================================================
+---------------------------+---------------------------+
| Iteration | Residual [N] |
+---------------------------+---------------------------+
| 0 | 1.0381e+05 |
| 1 | 1.0284e+05 |
| 2 | 1.0189e+05 |
| 3 | 9.9408e+04 |
| 4 | 9.6102e+04 |
| 5 | 8.8939e+04 |
| 6 | 7.6378e+04 |
| 7 | 5.3391e+04 |
| 8 | 5.3391e+04 |
| 9 | 5.3391e+04 |
| 10 | 5.2609e+04 |
| 11 | 5.2609e+04 |
| 12 | 5.2609e+04 |
| 13 | 5.1748e+04 |
| 14 | 5.1748e+04 |
| 15 | 5.1675e+04 |
| 16 | 5.0488e+04 |
| 17 | 5.0338e+04 |
| 18 | 4.6676e+04 |
| 19 | 4.6676e+04 |
| 20 | 3.8705e+04 |
| 21 | 3.8705e+04 |
| 22 | 2.7423e+04 |
| 23 | 9.9313e+03 |
| 24 | 9.9313e+03 |
| 25 | 9.9313e+03 |
| 26 | 9.9313e+03 |
| 27 | 9.9287e+03 |
| 28 | 4.5320e+03 |
| 29 | 3.7950e+03 |
| 30 | 3.7950e+03 |
| 31 | 1.1300e+03 |
| 32 | 1.1300e+03 |
| 33 | 1.1300e+03 |
| 34 | 8.1374e+02 |
| 35 | 7.4933e+02 |
| 36 | 3.2156e+02 |
| 37 | 3.2156e+02 |
| 38 | 3.2156e+02 |
| 39 | 1.2052e+02 |
| 40 | 1.2052e+02 |
| 41 | 1.2052e+02 |
| 42 | 5.8065e+01 |
| 43 | 5.8065e+01 |
| 44 | 5.8065e+01 |
| 45 | 5.4652e+01 |
| 46 | 2.7175e+01 |
| 47 | 1.0612e+01 |
| 48 | 1.0612e+01 |
| 49 | 9.6573e+00 |
| 50 | 6.3791e+00 |
| 51 | 2.5117e+00 |
| 52 | 2.5117e+00 |
| 53 | 2.5117e+00 |
| 54 | 1.8404e+00 |
| 55 | 1.4085e+00 |
| 56 | 1.0492e+00 |
+---------------------------+---------------------------+
=========================================================
=========================================================
OPTIMIZATION CONVERGENCE - 1200.0 RPM
=========================================================
+---------------------------+---------------------------+
| Iteration | Residual [N] |
+---------------------------+---------------------------+
| 0 | 9.7770e+03 |
| 1 | 3.5396e+03 |
| 2 | 3.5396e+03 |
| 3 | 3.5396e+03 |
| 4 | 3.5396e+03 |
| 5 | 3.5396e+03 |
| 6 | 2.7468e+02 |
| 7 | 2.7468e+02 |
| 8 | 2.7468e+02 |
| 9 | 2.7468e+02 |
| 10 | 2.7468e+02 |
| 11 | 2.7468e+02 |
| 12 | 2.7468e+02 |
| 13 | 2.7468e+02 |
| 14 | 2.7439e+02 |
| 15 | 9.0084e+01 |
| 16 | 9.0084e+01 |
| 17 | 9.0084e+01 |
| 18 | 6.2551e+01 |
| 19 | 4.2383e+01 |
| 20 | 2.7320e+01 |
| 21 | 1.9681e+01 |
| 22 | 1.3802e+01 |
| 23 | 1.0657e+01 |
| 24 | 8.5002e+00 |
| 25 | 5.2910e+00 |
| 26 | 3.2309e+00 |
| 27 | 2.2083e+00 |
| 28 | 1.8081e+00 |
| 29 | 1.1152e+00 |
| 30 | 1.0161e+00 |
+---------------------------+---------------------------+
=========================================================
1.4: Execution time#
To display the total time spent during the execution of the optimization process to find the equilibrium position of the shaft for all speeds specified as input in the PlainJournal class, use:
plainjournal.show_execution_time()
Execution time: 6.27 seconds
1.5: Pressure and temperature distributions#
To display the pressure and temperature distributions (2D contours and 3D profiles) for the last analyzed speed, use:
plainjournal.plot_results(show_plots=True);
1.6: Circumferential pressure profile#
To display the circumferential pressure profile at the middle plane of the bearing (or in the plane closest to it, depending on the number of elements in the axial direction) for the last analyzed speed, use:
plainjournal.plot_pressure_distribution()
Optionally, to display the circumferential pressure profile at a specific axial plane of the bearing, pass the axial element index as an argument to the plot_pressure_distribution method:
plainjournal.plot_pressure_distribution(2)
1.7: Coefficients comparison#
To display in table format the stiffness and damping coefficients of the bearing for all analyzed speeds, use:
plainjournal.show_coefficients_comparison()
=============================================================================================================================================================================================================================================================
DYNAMIC COEFFICIENTS COMPARISON TABLE
=============================================================================================================================================================================================================================================================
+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+
| Frequency [RPM] | kxx [N/m] | kxy [N/m] | kyx [N/m] | kyy [N/m] | cxx [N*s/m] | cxy [N*s/m] | cyx [N*s/m] | cyy [N*s/m] |
+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+
| 900.0 | 2.3587e+09 | 7.2577e+08 | -3.0873e+09 | 2.5047e+09 | 3.3502e+07 | -3.1875e+07 | -3.8124e+07 | 9.5236e+07 |
| 1200.0 | 2.2310e+09 | 7.5598e+08 | -2.9689e+09 | 2.3502e+09 | 2.3617e+07 | -2.0366e+07 | -2.4895e+07 | 6.7264e+07 |
+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+
=============================================================================================================================================================================================================================================================
1.8: Show results#
To display in table format the optimization results and the bearing dynamic coefficients for each analyzed speed, use:
plainjournal.show_results()
=====================================================================================
PLAIN JOURNAL RESULTS - 900.0 RPM
=====================================================================================
+---------------------------+---------------------------+---------------------------+
| Parameter | Value | Unit |
+---------------------------+---------------------------+---------------------------+
| Operating Speed | 900.0 | RPM |
| Eccentricity Ratio | 6.5259e-01 | - |
| Attitude Angle | -4.3884e+01 | deg |
| Load Fx | 0.0000e+00 | N |
| Load Fy | -1.1281e+05 | N |
| kxx (Stiffness) | 2.3587e+09 | N/m |
| kxy (Stiffness) | 7.2577e+08 | N/m |
| kyx (Stiffness) | -3.0873e+09 | N/m |
| kyy (Stiffness) | 2.5047e+09 | N/m |
| cxx (Damping) | 3.3502e+07 | N*s/m |
| cxy (Damping) | -3.1875e+07 | N*s/m |
| cyx (Damping) | -3.8124e+07 | N*s/m |
| cyy (Damping) | 9.5236e+07 | N*s/m |
| Optimization Success | True | - |
| Function Value | 1.0492e+00 | - |
| Iterations | 58 | - |
| Evaluations | 109 | - |
| Execution Time | 1.4276e+01 | s |
+---------------------------+---------------------------+---------------------------+
=====================================================================================
=====================================================================================
PLAIN JOURNAL RESULTS - 1200.0 RPM
=====================================================================================
+---------------------------+---------------------------+---------------------------+
| Parameter | Value | Unit |
+---------------------------+---------------------------+---------------------------+
| Operating Speed | 1200.0 | RPM |
| Eccentricity Ratio | 6.0528e-01 | - |
| Attitude Angle | -4.1085e+01 | deg |
| Load Fx | 0.0000e+00 | N |
| Load Fy | -1.1281e+05 | N |
| kxx (Stiffness) | 2.2310e+09 | N/m |
| kxy (Stiffness) | 7.5598e+08 | N/m |
| kyx (Stiffness) | -2.9689e+09 | N/m |
| kyy (Stiffness) | 2.3502e+09 | N/m |
| cxx (Damping) | 2.3617e+07 | N*s/m |
| cxy (Damping) | -2.0366e+07 | N*s/m |
| cyx (Damping) | -2.4895e+07 | N*s/m |
| cyy (Damping) | 6.7264e+07 | N*s/m |
| Optimization Success | True | - |
| Function Value | 1.0161e+00 | - |
| Iterations | 32 | - |
| Evaluations | 63 | - |
| Execution Time | 6.2727e+00 | s |
+---------------------------+---------------------------+---------------------------+
=====================================================================================
Optionally, to display the results only for a given speed (e.g., 900 rpm), use:
plainjournal._print_single_frequency_results(plainjournal.frequency[0])
=====================================================================================
PLAIN JOURNAL RESULTS - 900.0 RPM
=====================================================================================
+---------------------------+---------------------------+---------------------------+
| Parameter | Value | Unit |
+---------------------------+---------------------------+---------------------------+
| Operating Speed | 900.0 | RPM |
| Eccentricity Ratio | 6.5259e-01 | - |
| Attitude Angle | -4.3884e+01 | deg |
| Load Fx | 0.0000e+00 | N |
| Load Fy | -1.1281e+05 | N |
| kxx (Stiffness) | 2.3587e+09 | N/m |
| kxy (Stiffness) | 7.2577e+08 | N/m |
| kyx (Stiffness) | -3.0873e+09 | N/m |
| kyy (Stiffness) | 2.5047e+09 | N/m |
| cxx (Damping) | 3.3502e+07 | N*s/m |
| cxy (Damping) | -3.1875e+07 | N*s/m |
| cyx (Damping) | -3.8124e+07 | N*s/m |
| cyy (Damping) | 9.5236e+07 | N*s/m |
| Optimization Success | True | - |
| Function Value | 1.0492e+00 | - |
| Iterations | 58 | - |
| Evaluations | 109 | - |
| Execution Time | 1.4276e+01 | s |
+---------------------------+---------------------------+---------------------------+
=====================================================================================
1.9: Final remarks#
To perform quick analyses of short bearings, ROSS provides the CylindricalBearing class, which models the fluid film using a closed-form solution to the Reynolds equation. An example of how this class is used is shown below:
from ross import CylindricalBearing
from ross.units import Q_
cylindrical = CylindricalBearing(
n=0,
speed=Q_([1500, 2000], "RPM"),
weight=525,
bearing_length=Q_(30, "mm"),
journal_diameter=Q_(100, "mm"),
radial_clearance=Q_(0.1, "mm"),
oil_viscosity=0.1,
)
To display the value of a specific bearing dynamic coefficient (e.g., kxx) for all the analyzed speeds, use:
cylindrical.kxx
array([12807959.57740873, 13005276.62769462])
To summarize the inputs (geometric parameters and operational conditions) and outputs (equilibrium position and dynamic coefficients) of the CylindricalBraring class, use:
cylindrical._hover_info()
([0],
'Cylindrical Bearing at Node: 0<br>Length: 0.0300 m<br>Journal Diameter: 0.1000 m<br>Radial Clearance: 1.0000e-04 m<br>Oil Viscosity: 0.1000 Pa·s<br>Load: 525.00 N<br>Speed: 1500.00 ... 2000.00 RPM<br>Eccentricity: 0.2663 ... 0.2126<br>Attitude Angle: 0.1989 ... 0.1617 rad<br>Sommerfeld: 3.571e+00 ... 4.762e+00<br>Kxx: 1.281e+07 ... 1.301e+07 N/m<br>Kyy: 8.815e+06 ... 8.021e+06 N/m<br>Cxx: 2.329e+05 ... 2.249e+05 N·s/m<br>Cyy: 2.949e+05 ... 2.620e+05 N·s/m<br>')
Section 2: TiltingPad#
This tutorial provides a guide on how to simulate the thermo-hydrodynamic behavior of tilting pad journal bearings and compute their dynamic coefficients using the ROSS library.
Tilting pad bearings differ from fixed-geometry bearings because each pad can rotate independently around a pivot, allowing better stability and reducing the risk of oil whirl and oil whip instabilities. The analysis typically involves solving the Reynolds equation coupled with thermal effects.
Tilting pad bearings can be analyzed using two different thermal modeling approaches available in the ROSS library:
Adiabatic model (thermal_type=”adiabatic”): Only the oil film temperature is computed. Heat transfer within the pad is neglected, resulting in a faster solution suitable for preliminary analyses.
Full thermo-hydrodynamic model (thermal_type=”full”): Includes thermal coupling between the oil film and the pad material. This approach accounts for heat conduction inside the pads and requires iterative convergence, providing more accurate results for high-load or high-speed applications.
2.1: Adiabatic example of application#
First, to run a test case using an adiabatic thermal model (where only the oil film temperature is computed), use the following code:
from ross.bearings.tilting_pad import TiltingPad
from ross.units import Q_
import plotly.io as pio
pio.renderers.default = "notebook"
tilting_pad_adiabatic = TiltingPad(
n=1, # Nodal location of the bearing in the FE model of the rotor
frequency=Q_([3000, 5000], "RPM"), # Journal rotational speed
load=[884.05, -2670.4], # External load applied to the shaft [Fx, Fy]
journal_diameter=101.6e-3, # Journal diameter
radial_clearance=74.9e-6, # Radial clearance
pad_thickness=12.7e-3, # Pad thickness
pivot_angle=Q_([18, 90, 162, 234, 306], "deg"), # Pivot angular positions
pad_arc=Q_([60, 60, 60, 60, 60], "deg"), # Angular extent of each pad
pad_axial_length=Q_([50.8e-3] * 5, "m"), # Axial length of each pad
pre_load=[0.5] * 5, # Preload factor for each pad
offset=[0.5] * 5, # Pivot offset factor for each pad
lubricant="ISOVG32", # Lubricant type
oil_supply_temperature=Q_(40, "degC"), # Oil supply temperature
nx=30, # Number of elements in circumferential direction
nz=30, # Number of elements in axial direction
thermal_type="adiabatic", # Only oil film temperature is computed
equilibrium_type="match_eccentricity", # Method to determine equilibrium
eccentricity=0.35, # Prescribed eccentricity
attitude_angle=Q_(287.5, "deg"), # Attitude angle
initial_pads_angles=[
1.0747e-03,
7.2199e-04,
2.9340e-04,
3.4885e-04,
8.1554e-04,
], # Initial guesses to find the equilibrium position of the pads using an optimization process
solver_options={"xtol": 1e-2, "ftol": 1e-2, "maxiter": 1000}, # Solver options
)
2.2: Equilibrium calculation options#
In addition to the thermal model, the TiltingPad class requires defining how the equilibrium position of the rotor is determined. This is controlled by the equilibrium_type parameter.
Two options are available:
“match_eccentricity”: Uses the prescribed eccentricity and attitude angle, optimizing only the pad rotation angles to satisfy moment equilibrium.
“determine_eccentricity”: Computes the full equilibrium condition by simultaneously determining the eccentricity, attitude angle, and pad rotation angles required to balance the applied loads.
Using prescribed eccentricity (match_eccentricity)
In this approach, the eccentricity and attitude angle are prescribed, and an optimization process is performed to determine the pad rotation angles that satisfy the moment equilibrium condition.
tilting_pad_match = TiltingPad(
n=1,
frequency=Q_([3000], "RPM"),
load=[884.05, -2670.4],
journal_diameter=101.6e-3,
radial_clearance=74.9e-6,
pad_thickness=12.7e-3,
pivot_angle=Q_([18, 90, 162, 234, 306], "deg"),
pad_arc=Q_([60] * 5, "deg"),
pad_axial_length=Q_([50.8e-3] * 5, "m"),
pre_load=[0.5] * 5,
offset=[0.5] * 5,
lubricant="ISOVG32",
oil_supply_temperature=Q_(40, "degC"),
nx=30,
nz=30,
thermal_type="adiabatic",
equilibrium_type="match_eccentricity",
eccentricity=0.35,
attitude_angle=Q_(287.5, "deg"),
)
The convergence of the optimization process can be visualized using:
tilting_pad_match.show_optimization_convergence(show_plots=True)
=====================================================================================
OPTIMIZATION CONVERGENCE - 3000.0 RPM
=====================================================================================
+---------------------------+---------------------------+---------------------------+
| Pad | Iterations | Final Residual [N] |
+---------------------------+---------------------------+---------------------------+
| 1 | 20 | 0.000803 |
| 2 | 20 | 0.034964 |
| 3 | 20 | 0.010224 |
| 4 | 20 | 0.012591 |
| 5 | 24 | 0.056801 |
+---------------------------+---------------------------+---------------------------+
=====================================================================================
This plot shows the evolution of the residual during the optimization, indicating how the pad angles are adjusted to achieve equilibrium.
Determining full equilibrium (determine_eccentricity)
In this approach, the solver determines the full equilibrium condition by optimizing the eccentricity, attitude angle, and pad rotation angles simultaneously, based on the applied load.
tilting_pad_determine = TiltingPad(
n=1,
frequency=Q_([3000], "RPM"),
load=[884.05, -2670.4],
journal_diameter=101.6e-3,
radial_clearance=74.9e-6,
pad_thickness=12.7e-3,
pivot_angle=Q_([18, 90, 162, 234, 306], "deg"),
pad_arc=Q_([60] * 5, "deg"),
pad_axial_length=Q_([50.8e-3] * 5, "m"),
pre_load=[0.5] * 5,
offset=[0.5] * 5,
lubricant="ISOVG32",
oil_supply_temperature=Q_(40, "degC"),
nx=30,
nz=30,
thermal_type="adiabatic",
equilibrium_type="determine_eccentricity",
eccentricity=0.3, # Initial guess
attitude_angle=Q_(270, "deg"), # Initial guess
initial_pads_angles=[0.0] * 5,
)
To visualize the convergence of the optimization:
tilting_pad_determine.show_optimization_convergence(show_plots=True)
=========================================================
OPTIMIZATION CONVERGENCE - 3000.0 RPM
=========================================================
+---------------------------+---------------------------+
| Iteration | Residual [N] |
+---------------------------+---------------------------+
| 0 | 2770.810677 |
| 1 | 2780.392526 |
| 2 | 2661.294934 |
| 3 | 2771.258293 |
| 4 | 2872.874155 |
| 5 | 2960.375583 |
| 6 | 2620.192591 |
| 7 | 2594.481099 |
| 8 | 2524.842309 |
| 9 | 2405.028476 |
| 10 | 2484.103864 |
| 11 | 2427.080226 |
| 12 | 2311.514395 |
| 13 | 2053.043086 |
| 14 | 2078.491734 |
| 15 | 2098.667242 |
| 16 | 1913.618108 |
| 17 | 1447.117915 |
| 18 | 1943.786089 |
| 19 | 1680.720692 |
| 20 | 1451.306012 |
| 21 | 1019.271769 |
| 22 | 379.010791 |
| 23 | 798.605874 |
| 24 | 412.499676 |
| 25 | 212.147849 |
| 26 | 2235.065051 |
| 27 | 1213.811075 |
| 28 | 1454.880423 |
| 29 | 331.346950 |
| 30 | 1384.786589 |
| 31 | 2207.768759 |
| 32 | 468.002916 |
| 33 | 1255.726259 |
| 34 | 726.704157 |
| 35 | 1252.781757 |
| 36 | 233.590113 |
| 37 | 751.952600 |
| 38 | 381.345124 |
| 39 | 914.500512 |
| 40 | 364.209966 |
| 41 | 903.080263 |
| 42 | 188.730308 |
| 43 | 608.134921 |
| 44 | 228.683778 |
| 45 | 222.720928 |
| 46 | 372.030791 |
| 47 | 174.325952 |
| 48 | 388.388727 |
| 49 | 195.430043 |
| 50 | 442.784824 |
| 51 | 114.053909 |
| 52 | 445.461044 |
| 53 | 89.743458 |
| 54 | 95.900774 |
| 55 | 217.613250 |
| 56 | 133.316235 |
| 57 | 246.885105 |
| 58 | 38.235126 |
| 59 | 224.523689 |
| 60 | 72.783497 |
| 61 | 401.249632 |
| 62 | 101.702035 |
| 63 | 215.136496 |
| 64 | 108.068922 |
| 65 | 70.122885 |
| 66 | 160.828214 |
| 67 | 39.132544 |
| 68 | 72.593907 |
| 69 | 157.362003 |
| 70 | 46.815766 |
| 71 | 125.144235 |
| 72 | 57.265317 |
| 73 | 119.335986 |
| 74 | 42.002555 |
| 75 | 115.877363 |
| 76 | 41.798721 |
| 77 | 61.943611 |
| 78 | 112.051542 |
| 79 | 21.240040 |
| 80 | 79.148738 |
| 81 | 44.007986 |
| 82 | 39.463615 |
| 83 | 100.185947 |
| 84 | 11.087302 |
| 85 | 57.655429 |
| 86 | 14.002011 |
| 87 | 32.150575 |
| 88 | 56.366509 |
| 89 | 16.484114 |
| 90 | 41.011490 |
| 91 | 22.795981 |
| 92 | 52.279747 |
| 93 | 35.106189 |
| 94 | 90.817538 |
| 95 | 19.187340 |
| 96 | 54.013515 |
| 97 | 26.236288 |
| 98 | 35.748815 |
| 99 | 23.382044 |
| 100 | 22.305119 |
| 101 | 12.912899 |
| 102 | 17.333035 |
| 103 | 13.190541 |
| 104 | 35.148031 |
| 105 | 9.622784 |
| 106 | 44.096384 |
| 107 | 6.345590 |
| 108 | 13.599627 |
| 109 | 29.884186 |
| 110 | 10.354460 |
| 111 | 17.015504 |
| 112 | 8.182687 |
| 113 | 16.372179 |
| 114 | 8.643211 |
| 115 | 19.142136 |
| 116 | 10.893331 |
| 117 | 11.282043 |
| 118 | 7.647039 |
| 119 | 21.361550 |
| 120 | 8.312675 |
| 121 | 12.116053 |
| 122 | 7.050927 |
| 123 | 8.583360 |
| 124 | 19.626086 |
| 125 | 6.238154 |
| 126 | 9.952924 |
| 127 | 6.794607 |
| 128 | 8.148123 |
| 129 | 13.647781 |
| 130 | 6.553083 |
| 131 | 9.352311 |
| 132 | 6.420322 |
| 133 | 9.930604 |
| 134 | 6.893616 |
| 135 | 11.405913 |
| 136 | 6.172205 |
| 137 | 7.048799 |
| 138 | 7.254029 |
| 139 | 6.296680 |
| 140 | 9.838039 |
| 141 | 6.328158 |
| 142 | 6.308447 |
| 143 | 6.420237 |
| 144 | 6.391880 |
| 145 | 6.236287 |
| 146 | 11.014215 |
| 147 | 6.440194 |
| 148 | 8.428903 |
| 149 | 6.526327 |
| 150 | 6.160907 |
| 151 | 6.358681 |
| 152 | 6.383541 |
| 153 | 6.159232 |
| 154 | 6.265824 |
| 155 | 6.273348 |
| 156 | 6.439182 |
| 157 | 6.620094 |
| 158 | 6.201371 |
| 159 | 6.166974 |
| 160 | 6.624769 |
| 161 | 6.183597 |
| 162 | 6.167783 |
| 163 | 6.160166 |
| 164 | 6.184269 |
| 165 | 6.173375 |
| 166 | 6.170500 |
| 167 | 6.191186 |
| 168 | 6.146080 |
| 169 | 6.259152 |
| 170 | 6.167392 |
| 171 | 6.199255 |
| 172 | 6.155295 |
| 173 | 6.211538 |
| 174 | 6.144368 |
| 175 | 6.164724 |
| 176 | 6.171613 |
| 177 | 6.155398 |
| 178 | 6.155135 |
| 179 | 6.146956 |
| 180 | 6.172091 |
| 181 | 6.159716 |
| 182 | 6.156999 |
| 183 | 6.144505 |
| 184 | 6.149163 |
| 185 | 6.166883 |
| 186 | 6.150304 |
| 187 | 6.156033 |
| 188 | 6.147160 |
| 189 | 6.161184 |
| 190 | 6.147614 |
| 191 | 6.148490 |
| 192 | 6.160704 |
| 193 | 6.146608 |
| 194 | 6.150590 |
| 195 | 6.145218 |
| 196 | 6.146638 |
| 197 | 6.149235 |
| 198 | 6.145619 |
| 199 | 6.147338 |
| 200 | 6.153499 |
| 201 | 6.145189 |
| 202 | 6.145298 |
| 203 | 6.145352 |
| 204 | 6.147760 |
| 205 | 6.145724 |
| 206 | 6.146338 |
| 207 | 6.145255 |
| 208 | 6.144550 |
| 209 | 6.145488 |
| 210 | 6.145563 |
| 211 | 6.145356 |
| 212 | 6.144503 |
| 213 | 6.147092 |
+---------------------------+---------------------------+
=========================================================
In this case, the optimization process is more complex, as it involves additional variables (eccentricity and attitude angle), which may lead to a higher number of iterations.
The convergence behavior depends on the initial guesses and operating conditions. Poor initial estimates may increase the number of iterations or affect convergence stability, especially when using determine_eccentricity.
Each entry corresponds to a frequency and contains the residual evolution during the optimization process.
2.3: Execution time#
To display the total execution time of the thermo-hydrodynamic analysis, use:
tilting_pad_adiabatic.show_execution_time()
Execution time: 3.11 seconds
2.4: Pressure and temperature distributions#
The command (plot_results) generates a set of plots showing the pressure and temperature distributions for all pads, including scatter and contour visualizations.
The freq_index parameter specifies which operating speed to visualize. When multiple rotational speeds are defined in the frequency input, each one is stored internally as an index:
freq_index=0 → first speed
freq_index=1 → second speed
and so on
In this example, since multiple rotational speeds are defined, freq_index=0 corresponds to the first speed and freq_index=1 to the second one.
tilting_pad_adiabatic.plot_results(show_plots=True, freq_index=1)
2.5: Pressure distribution (single pad)#
This plot displays the pressure distribution for a single pad. The selected pad is automatically defined by the model and typically corresponds to the most loaded pad (i.e., the pad with the highest pressure levels).
fig = tilting_pad_adiabatic.plot_pressure_distribution()
fig.show()
This plot shows the pressure field in the lubricant film for the selected pad, where:
The circumferential direction (θ) is shown along the X-axis
The axial direction is shown along the Y-axis
The pressure magnitude is represented along the Z-axis
Note that this plot displays the results for a single pad. The selected pad is internally defined by the model.
2.6: Thermal distribution (single pad)#
In addition to the global thermal visualization, it is also possible to analyze the temperature distribution for a specific pad at a given operating condition.
This is particularly useful when multiple rotational speeds are defined, as it allows isolating the thermal behavior for a selected frequency and pad.
tilting_pad_adiabatic.plot_thermal_pad_results(freq_index=0, pad_index=0)
In this example:
freq_index=0 selects the first rotational speed
pad_index=0 selects the first pad (indexing starts at 0)
This plot provides a detailed contour representation of the temperature field within the lubricant film for the selected pad, enabling a more localized thermal analysis.
2.7: Film average temperature#
To analyze the average temperature of the lubricant film along each pad, the axial-averaged temperature can be plotted as a function of the local pad angle:
fig = tilting_pad_adiabatic.plot_film_average_temperature()
fig.show()
This plot shows the variation of the average oil film temperature along the pad arc for each pad.
The temperature is averaged along the axial direction, allowing a clearer visualization of circumferential thermal gradients and the identification of hot spots.
This representation is particularly useful for comparing the thermal behavior between pads and understanding how heat is distributed along the bearing.
2.8: Coefficients comparison#
The stiffness and damping coefficients for each rotational speed can be directly accessed as:
tilting_pad_adiabatic.show_coefficients_comparison()
=============================================================================================================================================================================================================================================================
DYNAMIC COEFFICIENTS COMPARISON TABLE
=============================================================================================================================================================================================================================================================
+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+
| Frequency [RPM] | kxx [N/m] | kxy [N/m] | kyx [N/m] | kyy [N/m] | cxx [N*s/m] | cxy [N*s/m] | cyx [N*s/m] | cyy [N*s/m] |
+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+
| 3000.0 | 8.9482e+07 | -3.3834e+07 | -3.3834e+07 | 1.1652e+08 | 2.8196e+05 | -2.8321e+04 | -2.8321e+04 | 3.3166e+05 |
| 5000.0 | 1.4566e+08 | -5.2125e+07 | -5.2125e+07 | 1.8782e+08 | 2.7119e+05 | -2.5064e+04 | -2.5064e+04 | 3.1678e+05 |
+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+---------------------------+
=============================================================================================================================================================================================================================================================
2.9: Show results#
To display in table format the optimization results and the bearing dynamic coefficients for each speed analyzed, use the following command:
tilting_pad_adiabatic.show_results()
======================================================================
TILTING PAD BEARING RESULTS - 3000.0 RPM
======================================================================
+----------------------+----------------------+----------------------+
| Parameter | Value | Unit |
+----------------------+----------------------+----------------------+
| Operating Speed | 3000.0 | RPM |
| Equilibrium Type | match_eccentricity | - |
| Number of Pads | 5 | - |
| Maximum Pressure | 1.8953e+06 | Pa |
| Maximum Temperature | 56.71 | °C |
| Minimum Film | 5.0040e-05 | m |
| Thickness | | |
| Eccentricity | 0.3500 | - |
| Attitude Angle | 287.50 | ° |
| Total Force X | -8.8446e+02 | N |
| Total Force Y | 2.6730e+03 | N |
| kxx (Stiffness) | 8.9482e+07 | N/m |
| kxy (Stiffness) | -3.3834e+07 | N/m |
| kyx (Stiffness) | -3.3834e+07 | N/m |
| kyy (Stiffness) | 1.1652e+08 | N/m |
| cxx (Damping) | 2.8196e+05 | N*s/m |
| cxy (Damping) | -2.8321e+04 | N*s/m |
| cyx (Damping) | -2.8321e+04 | N*s/m |
| cyy (Damping) | 3.3166e+05 | N*s/m |
+----------------------+----------------------+----------------------+
----------------------------------------------------------------------
PAD ROTATION ANGLES
----------------------------------------------------------------------
+----------------+----------------+----------------+----------------+
| Pad # | Moment [N·m] | Angle [rad] | Angle [°] |
+----------------+----------------+----------------+----------------+
| 1 | -2.7358e-02 | 1.0727e-03 | 6.1461e-02 |
| 2 | 3.4964e-02 | 7.2933e-04 | 4.1788e-02 |
| 3 | 1.0409e-02 | 2.9538e-04 | 1.6924e-02 |
| 4 | -1.7846e-02 | 3.4920e-04 | 2.0007e-02 |
| 5 | 5.6801e-02 | 8.1735e-04 | 4.6831e-02 |
+----------------+----------------+----------------+----------------+
======================================================================
======================================================================
TILTING PAD BEARING RESULTS - 5000.0 RPM
======================================================================
+----------------------+----------------------+----------------------+
| Parameter | Value | Unit |
+----------------------+----------------------+----------------------+
| Operating Speed | 5000.0 | RPM |
| Equilibrium Type | match_eccentricity | - |
| Number of Pads | 5 | - |
| Maximum Pressure | 2.7977e+06 | Pa |
| Maximum Temperature | 64.89 | °C |
| Minimum Film | 5.0040e-05 | m |
| Thickness | | |
| Eccentricity | 0.3500 | - |
| Attitude Angle | 287.50 | ° |
| Total Force X | -1.4315e+03 | N |
| Total Force Y | 4.2994e+03 | N |
| kxx (Stiffness) | 1.4566e+08 | N/m |
| kxy (Stiffness) | -5.2125e+07 | N/m |
| kyx (Stiffness) | -5.2125e+07 | N/m |
| kyy (Stiffness) | 1.8782e+08 | N/m |
| cxx (Damping) | 2.7119e+05 | N*s/m |
| cxy (Damping) | -2.5064e+04 | N*s/m |
| cyx (Damping) | -2.5064e+04 | N*s/m |
| cyy (Damping) | 3.1678e+05 | N*s/m |
+----------------------+----------------------+----------------------+
----------------------------------------------------------------------
PAD ROTATION ANGLES
----------------------------------------------------------------------
+----------------+----------------+----------------+----------------+
| Pad # | Moment [N·m] | Angle [rad] | Angle [°] |
+----------------+----------------+----------------+----------------+
| 1 | -1.3181e-02 | 1.0978e-03 | 6.2898e-02 |
| 2 | -5.4086e-02 | 7.4670e-04 | 4.2783e-02 |
| 3 | 3.7730e-02 | 3.1599e-04 | 1.8105e-02 |
| 4 | 2.7346e-02 | 3.7790e-04 | 2.1652e-02 |
| 5 | -7.9700e-03 | 8.4953e-04 | 4.8674e-02 |
+----------------+----------------+----------------+----------------+
======================================================================
2.10: Thermal example of application#
On the other hand, to include thermal coupling between the oil film, shaft geometry, and the pad material, use:
from ross.bearings.tilting_pad import TiltingPad
from ross.units import Q_
tilting_pad_fullThermo = TiltingPad(
n=1, # Nodal location of the bearing in the FE model of the rotor
frequency=Q_([3000], "RPM"), # Journal rotational speed
load=[884.05, -2670.4], # External load applied to the shaft [Fx, Fy]
journal_diameter=101.6e-3, # Journal diameter
radial_clearance=74.9e-6, # Radial clearance
pad_thickness=12.7e-3, # Pad thickness
pivot_angle=Q_([18, 90, 162, 234, 306], "deg"), # Pivot angular positions
pad_arc=Q_([60, 60, 60, 60, 60], "deg"), # Angular extent of each pad
pad_axial_length=Q_([50.8e-3] * 5, "m"), # Axial length of each pad
pre_load=[0.5] * 5, # Preload factor for each pad
offset=[0.5] * 5, # Pivot offset factor for each pad
lubricant="ISOVG32", # Lubricant type
oil_supply_temperature=Q_(40, "degC"), # Oil supply temperature
nx=30, # Number of elements in circumferential direction
nz=30, # Number of elements in axial direction
thermal_type="full", # Includes pad thermal conduction
nr_pad=20, # Number of elements in radial direction inside the pad
k_pad=45.0, # Pad thermal conductivity [W/(m·K)]
h_edge=2000.0, # Heat transfer coefficient at pad edges [W/(m²·K)]
equilibrium_type="match_eccentricity", # Method to determine equilibrium
eccentricity=0.35, # Prescribed eccentricity
attitude_angle=Q_(287.5, "deg"), # Attitude angle
initial_pads_angles=[
1.0747e-03,
7.2199e-04,
2.9340e-04,
3.4885e-04,
8.1554e-04,
], # Initial guesses to find the equilibrium position of the pads using an optimization process
solver_options={"xtol": 1e-2, "ftol": 1e-2, "maxiter": 1000}, # Solver options
hot_oil_carry_over=0.8, # Hot oil carry-over factor
journal_temperature=Q_(25, "degC"), # Journal temperature
max_jtemp_iter=150, # Maximum number of iterations for journal temperature convergence
jtemp_error=0.5, # Convergence tolerance for journal temperature [°C]
max_inlet_iterations=40, # Maximum number of outer iterations for inlet temperature convergence
inlet_temperature_tolerance=0.3, # Convergence tolerance for pad inlet temperatures [°C]
)
The same commands presented in the adiabatic example are also valid for the full thermal model. This includes:
Running the thermo-hydrodynamic analysis
Running the optimization options and accessing their convergence
Displaying execution time
Visualizing pressure and temperature fields
Obtaining dynamic coefficients and results
However, when using the full thermal model, additional thermal effects are considered, allowing for a more detailed analysis of temperature distribution and heat transfer within the pads.
2.11: Babbitt surface temperature#
To evaluate the temperature at the pad surface (babbitt layer), use:
fig = tilting_pad_fullThermo.plot_babbitt_surface_temperature()
fig.show()
This plot represents the temperature distribution at the pad–film interface, which is critical for assessing thermal loading and potential material degradation.
2.12: Solid pad temperature distribution#
To visualize the temperature distribution inside a specific pad (index = 0) material (solid domain):
fig = tilting_pad_fullThermo.plot_solid_pad_results(pad_index=0)
fig.show()
These contour plots show the temperature variation across the pad thickness, providing insight into heat conduction and thermal gradients within the pad material.
Section 3: ThrustPad#
This tutorial demonstrates how to use the ThrustPad class from ROSS to simulate the thermo-hydrodynamic behavior of tilting-pad thrust bearings and compute their dynamic coefficients for different operating conditions.
Thrust bearings are designed to support axial (thrust) loads and are commonly found in rotating machinery such as turbines, compressors, and pumps. The tilting-pad design allows each pad to tilt around a pivot point, creating a converging oil film that generates the hydrodynamic pressure needed to support the load.
The ThrustPad class implements a Thermo-Hydro-Dynamic (THD) model that solves the coupled Reynolds (pressure) and energy (temperature) equations for each pad. The main outputs are:
Pressure and temperature distributions across each pad
Film thickness at the pivot point
Axial stiffness coefficient
kzz[N/m]Axial damping coefficient
czz[N·s/m]
There are two modes for determining the pad equilibrium position:
calculate: The solver finds the equilibrium film thickness and inclination angles simultaneously by minimizing residual forces and moments.imposed: The film thickness at the pivot is fixed by the user; only the inclination angles are optimized.
3.1: Calculate mode#
In the calculate mode, the equilibrium film thickness is determined automatically by the solver together with the radial and circumferential inclination angles. Initial estimates must be provided as starting points for the optimization.
The example below represents a large industrial thrust bearing (e.g., in a hydraulic turbine) with 12 pads, operating at 90 RPM under an axial load of 13.32 MN:
from ross.bearings.thrust_pad import ThrustPad
from ross.units import Q_
import plotly.io as pio
pio.renderers.default = "notebook"
bearing_calculate = ThrustPad(
n=1, # Node in the FE model of the rotor
pad_inner_radius=Q_(1150, "mm"), # Inner radius of the pad
pad_outer_radius=Q_(1725, "mm"), # Outer radius of the pad
pad_pivot_radius=Q_(1442.5, "mm"), # Radius of the pivot point
pad_arc_length=Q_(26, "deg"), # Arc length of each pad
angular_pivot_position=Q_(15, "deg"), # Angular position of the pivot
oil_supply_temperature=Q_(40, "degC"), # Oil inlet temperature
lubricant="ISOVG68", # Lubricant type
n_pad=12, # Number of pads
n_theta=10, # Mesh elements in circumferential dir.
n_radial=10, # Mesh elements in radial dir.
frequency=Q_([90], "RPM"), # Rotor speed
equilibrium_position_mode="calculate", # Auto-calculate film thickness
axial_load=13.320e6, # Total axial load [N]
radial_inclination_angle=Q_(-2.75e-04, "rad"), # Initial radial inclination angle
circumferential_inclination_angle=Q_(
-1.70e-05, "rad"
), # Initial circumferential inclination
initial_film_thickness=Q_(0.2, "mm"), # Initial film thickness estimate
)
The ThrustPad class automatically runs the full THD analysis during instantiation. The inner loop solves the Reynolds and energy equations iteratively until the viscosity field converges, while the outer loop finds the equilibrium position using Nelder-Mead optimization.
3.2: Optimization convergence#
To display — in both table and graph formats — the evolution of the residual force/moment during the optimization process, use:
bearing_calculate.show_optimization_convergence(show_plots=True)
=========================================================
OPTIMIZATION CONVERGENCE - 90.0 RPM
=========================================================
+---------------------------+---------------------------+
| Iteration | Residual [N] |
+---------------------------+---------------------------+
| 1 | 378947.112688 |
| 2 | 27724.023669 |
| 3 | 4753.901941 |
| 4 | 2267.382187 |
| 5 | 381.146245 |
| 6 | 25.796494 |
| 7 | 6.502821 |
| 8 | 0.143314 |
| 9 | 0.162634 |
| 10 | 0.113231 |
| 11 | 0.084976 |
+---------------------------+---------------------------+
=========================================================
Each row in the table corresponds to one iteration of the outer convergence loop. The residual represents the net force/moment imbalance [N] at the current pad equilibrium estimate. The optimization converges when this residual falls below the tolerance_force_moment threshold (default: 0.1 N).
3.3: Execution time#
To display the total wall-clock time spent during the analysis, use:
bearing_calculate.show_execution_time()
Execution time: 22.34 seconds
3.4: Pressure and temperature distributions#
The plot_results method generates both 3D surface plots and 2D contour plots for the pressure and temperature fields across the bearing pad.
The freq_index parameter selects which operating speed to visualize. Since a single speed was defined (90 RPM), freq_index=0 corresponds to that speed:
figs = bearing_calculate.plot_results(show_plots=True, freq_index=0)