ENH: CP and Thrust Eccentricity Effects Generate Roll Moment

rocketpy/simulation/flight.py

@@ -1494,6 +1494,9 @@
                 M3 += M3f - M3d
             except AttributeError:
                 pass
+        # Off center moment
+        M3 += self.rocket.cp_eccentricity_x * R2 - self.rocket.cp_eccentricity_y * R1
+
         # Calculate derivatives
         # Angular acceleration
         alpha1 = (
@@ -1651,10 +1654,6 @@
                     R3 = air_brakes_force  # Substitutes rocket drag coefficient
                 else:
                     R3 += air_brakes_force
-        ## Off center moment
-        M1 += self.rocket.cp_eccentricity_y * R3
-        M2 -= self.rocket.cp_eccentricity_x * R3
-
         # Get rocket velocity in body frame
         vB = Kt @ v
         # Calculate lift and moment for each component of the rocket
@@ -1724,11 +1723,24 @@
                 M3 += M3f - M3d
             except AttributeError:
                 pass
+
+        # Off center moment
+        thrust = self.rocket.motor.thrust.get_value_opt(t)
+        M1 += (
+            self.rocket.cp_eccentricity_y * R3
+            + self.rocket.thrust_eccentricity_x * thrust
+        )
+        M2 -= (
+            self.rocket.cp_eccentricity_x * R3
+            - self.rocket.thrust_eccentricity_y * thrust
+        )
+        M3 += self.rocket.cp_eccentricity_x * R2 - self.rocket.cp_eccentricity_y * R1
+
         weightB = Kt @ Vector([0, 0, -total_mass * self.env.gravity.get_value_opt(z)])
         T00 = total_mass * r_CM
         T03 = 2 * total_mass_dot * (r_NOZ - r_CM) - 2 * total_mass * r_CM_dot
         T04 = (
-            Vector([0, 0, self.rocket.motor.thrust.get_value_opt(t)])
+            Vector([0, 0, thrust])
             - total_mass * r_CM_ddot
             - 2 * total_mass_dot * r_CM_dot
             + total_mass_ddot * (r_NOZ - r_CM)

tests/test_flight.py

@@ -295,6 +295,37 @@ def test_rolling_flight(
     assert test_flight.all_info() == None
 
 
+@patch("matplotlib.pyplot.show")
+def test_eccentricity_on_flight(
+    mock_show,
+    example_plain_env,
+    cesaroni_m1670,
+    calisto,
+    calisto_nose_cone,
+    calisto_trapezoidal_fins,
+    calisto_tail,
+):
+    test_rocket = calisto
+
+    test_rocket.set_rail_buttons(0.082, -0.618)
+    test_rocket.add_motor(cesaroni_m1670, position=-1.373)
+    calisto.aerodynamic_surfaces.add(calisto_trapezoidal_fins, position=-1.04956)
+    calisto.aerodynamic_surfaces.add(calisto_nose_cone, 1.160)
+    calisto.aerodynamic_surfaces.add(calisto_tail, -1.313)
+    calisto.add_cm_eccentricity(x=-0.01, y=-0.01)
+
+    test_flight = Flight(
+        rocket=test_rocket,
+        environment=example_plain_env,
+        rail_length=5.2,
+        inclination=85,
+        heading=0,
+        terminate_on_apogee=True,
+    )
+
+    assert test_flight.all_info() == None
+
+
 @patch("matplotlib.pyplot.show")
 def test_air_brakes_flight(mock_show, flight_calisto_air_brakes):
     """Test the flight of a rocket with air brakes. This test only validates

tests/unit/test_rocket.py

@@ -496,3 +496,26 @@ def test_get_inertia_tensor_derivative_at_time(calisto):
     assert pytest.approx(0, atol) == inertia_tensor.y[2]
     assert pytest.approx(0, atol) == inertia_tensor.z[0]
     assert pytest.approx(0, atol) == inertia_tensor.z[1]
+
+
+def test_add_cp_eccentricity(calisto):
+    """Test add_cp_eccentricity method of the Rocket class."""
+    calisto.add_cp_eccentricity(0.1, 0.1)
+    assert calisto.cp_eccentricity_x == 0.1
+    assert calisto.cp_eccentricity_y == 0.1
+
+
+def test_add_thrust_eccentricity(calisto):
+    """Test add_thrust_eccentricity method of the Rocket class."""
+    calisto.add_thrust_eccentricity(0.1, 0.1)
+    assert calisto.thrust_eccentricity_x == 0.1
+    assert calisto.thrust_eccentricity_y == 0.1
+
+
+def test_add_cm_eccentricity(calisto):
+    """Test add_cm_eccentricity method of the Rocket class."""
+    calisto.add_cm_eccentricity(-0.1, -0.1)
+    assert calisto.cp_eccentricity_x == 0.1
+    assert calisto.cp_eccentricity_y == 0.1
+    assert calisto.thrust_eccentricity_x == 0.1
+    assert calisto.thrust_eccentricity_y == 0.1