Repository: ross-rotordynamics/ross (python)

Issue:
Fix bearing type definition in fluid flow
The bearing is currently being defined wrongly. Instead of taking the diameter to calculate the ratio between it and the bearing length, we are taking the radius.

Pull request description:
Fix bearing type definition and add new test
Fix the bearing type definition and also add a new test that checks if a second analytical way to calculate the pressure matrix matches the numerical way.
Close #175 

Patch:
diff --git a/ross/fluid_flow/fluid_flow.py b/ross/fluid_flow/fluid_flow.py
index 40eadf23..516018fd 100644
--- a/ross/fluid_flow/fluid_flow.py
+++ b/ross/fluid_flow/fluid_flow.py
@@ -114,8 +114,8 @@ class PressureMatrix:
     bearing_type: str
         type of structure. 'short_bearing': short; 'long_bearing': long;
         'medium_size': in between short and long.
-        if length/radius_stator <= 1/8 it is short.
-        if length/radius_stator > 4 it is long.
+        if length/radius_stator <= 1/4 it is short.
+        if length/radius_stator > 8 it is long.
     radial_clearance: float
         Difference between both stator and rotor radius, regardless of eccentricity.
     analytical_pressure_matrix_available: bool
@@ -195,9 +195,9 @@ class PressureMatrix:
         self.radial_clearance = self.radius_stator - self.radius_rotor
         self.difference_between_radius = radius_stator - radius_rotor
         self.bearing_type = ""
-        if self.length / self.radius_stator <= 1 / 8:
+        if self.length / self.radius_stator <= 1 / 4:
             self.bearing_type = "short_bearing"
-        elif self.length / self.radius_stator > 4:
+        elif self.length / self.radius_stator > 8:
             self.bearing_type = "long_bearing"
         else:
             self.bearing_type = "medium_size"
@@ -234,7 +234,8 @@ class PressureMatrix:
         self.numerical_pressure_matrix_available = False
 
     def calculate_pressure_matrix_analytical(self):
-        """This function calculates the pressure matrix analytically.
+        """This function calculates the pressure matrix analytically, based on the book Tribology Series vol. 33, by
+        Frene et al., chapter 5.
         """
         if self.bearing_type == "short_bearing":
             for i in range(self.nz):
@@ -248,6 +249,22 @@ class PressureMatrix:
         self.analytical_pressure_matrix_available = True
         return self.p_mat_analytical
 
+    def calculate_pressure_matrix_analytical2(self):
+        """This function calculates the pressure matrix analytically, based on the chapter Linear and Nonlinear
+        Rotordynamics, by Ishida and Yamamoto, from the book Flow-Induced Vibrations.
+        """
+        if self.bearing_type == "short_bearing":
+            for i in range(self.nz):
+                for j in range(self.ntheta):
+                    # fmt: off
+                    self.p_mat_analytical[i][j] = (3 * self.visc / ((self.difference_between_radius ** 2) *
+                                                                    (1. + self.eccentricity_ratio * np.cos(
+                                                                        j * self.dtheta)) ** 3)) * \
+                                                  (-self.eccentricity_ratio * self.omega * np.sin(j * self.dtheta)) * \
+                                                  (((i * self.dz - (self.length / 2)) ** 2) - (self.length ** 2) / 4)
+                    # fmt: on
+        self.analytical_pressure_matrix_available = True
+
     def calculate_coefficients(self):
         """This function calculates the constants that form the Poisson equation
         of the discrete pressure (central differences in the second
@@ -444,7 +461,7 @@ class PressureMatrix:
             warnings.warn(
                 "Function get_rotor_load suitable only for short bearings. "
                 "The ratio between the bearing length and its radius should be less or "
-                "equal to 0.125. Currently we have "
+                "equal to 0.25. Currently we have "
                 + str(self.length / self.radius_stator)
                 + "."
             )
@@ -498,7 +515,7 @@ class PressureMatrix:
             warnings.warn(
                 "Function calculate_eccentricity_ratio suitable only for short bearings. "
                 "The ratio between the bearing length and its radius should be less or "
-                "equal to 0.125. Currently we have "
+                "equal to 0.25. Currently we have "
                 + str(self.length / self.radius_stator)
                 + "."
             )
@@ -528,7 +545,7 @@ class PressureMatrix:
             warnings.warn(
                 "Function get_analytical_stiffness_matrix suitable only for short bearings. "
                 "The ratio between the bearing length and its radius should be less or "
-                "equal to 0.125. Currently we have "
+                "equal to 0.25. Currently we have "
                 + str(self.length / self.radius_stator)
                 + "."
             )
@@ -562,7 +579,7 @@ class PressureMatrix:
             warnings.warn(
                 "Function get_analytical_damping_matrix suitable only for short bearings. "
                 "The ratio between the bearing length and its radius should be less or "
-                "equal to 0.125. Currently we have "
+                "equal to 0.25. Currently we have "
                 + str(self.length / self.radius_stator)
                 + "."
             )


Test patch:
diff --git a/ross/tests/test_fluid_flow.py b/ross/tests/test_fluid_flow.py
index ebc3ea11..f6c62741 100644
--- a/ross/tests/test_fluid_flow.py
+++ b/ross/tests/test_fluid_flow.py
@@ -114,3 +114,13 @@ def test_numerical_abs_error():
     error = np.linalg.norm(bearing.p_mat_analytical[:][int(bearing.nz / 2)] -
                            bearing.p_mat_numerical[:][int(bearing.nz / 2)], ord=np.inf)
     assert math.isclose(error, 0, abs_tol=0.001)
+
+
+def test_numerical_abs_error2():
+    bearing = fluid_flow_numerical()
+    bearing.calculate_pressure_matrix_analytical2()
+    bearing.calculate_pressure_matrix_numerical()
+    error = np.linalg.norm(bearing.p_mat_analytical[:][int(bearing.nz / 2)] -
+                           bearing.p_mat_numerical[:][int(bearing.nz / 2)], ord=np.inf)
+    assert math.isclose(error, 0, abs_tol=0.001)
+