From 554946d4b02be5380793558987416d4aaa76b215 Mon Sep 17 00:00:00 2001
From: Ayoub Abgo <ayoubdussap@gmail.com>
Date: Wed, 28 Jan 2026 16:25:34 +0100
Subject: [PATCH 01/10] feat: implement cluster_motor class with dynamic
 inertias

---
 rocketpy/motors/cluster_motor.py | 228 +++++++++++++++++++++++++++++++
 1 file changed, 228 insertions(+)
 create mode 100644 rocketpy/motors/cluster_motor.py

diff --git a/rocketpy/motors/cluster_motor.py b/rocketpy/motors/cluster_motor.py
new file mode 100644
index 000000000..0b8964d6f
--- /dev/null
+++ b/rocketpy/motors/cluster_motor.py
@@ -0,0 +1,228 @@
+import matplotlib.pyplot as plt
+import numpy as np
+from rocketpy import Function
+from rocketpy.motors import Motor
+
+
+class ClusterMotor(Motor):
+    """
+    A class representing a cluster of N identical motors arranged symmetrically.
+
+    This class aggregates the physical properties (thrust, mass, inertia) of
+    multiple motors using the Parallel Axis Theorem (Huygens-Steiner theorem).
+
+    Attributes
+    ----------
+    motor : SolidMotor
+        The single motor instance used in the cluster.
+    number : int
+        The number of motors in the cluster.
+    radius : float
+        The radial distance from the rocket's central axis to the center of each motor.
+    """
+
+    def __init__(self, motor, number, radius):
+        """
+        Initialize the ClusterMotor.
+
+        Parameters
+        ----------
+        motor : SolidMotor
+            The base motor to be clustered.
+        number : int
+            Number of motors. Must be >= 2.
+        radius : float
+            Distance from center of rocket to center of motor (m).
+        """
+        self.motor = motor
+        self.number = number
+        self.radius = radius
+        dry_inertia_cluster = self._calculate_dry_inertia()
+
+        super().__init__(
+            thrust_source=motor.thrust_source,
+            nozzle_radius=motor.nozzle_radius,
+            burn_time=motor.burn_time,
+            dry_mass=motor.dry_mass * number,
+            dry_inertia=dry_inertia_cluster,
+            center_of_dry_mass_position=motor.center_of_dry_mass_position,
+            coordinate_system_orientation=motor.coordinate_system_orientation,
+            interpolation_method="linear",
+        )
+
+        self.throat_radius = motor.throat_radius
+        self.grain_number = motor.grain_number
+        self.grain_density = motor.grain_density
+        self.grain_outer_radius = motor.grain_outer_radius
+        self.grain_initial_inner_radius = motor.grain_initial_inner_radius
+        self.grain_initial_height = motor.grain_initial_height
+        self.grains_center_of_mass_position = motor.grains_center_of_mass_position
+        self._thrust = self.motor.thrust * self.number
+        self._propellant_mass = self.motor.propellant_mass * self.number
+        self._propellant_initial_mass = self.number * self.motor.propellant_initial_mass
+        self._center_of_propellant_mass = self.motor.center_of_propellant_mass
+        Ixx_term1 = self.motor.propellant_I_11 * self.number
+        Ixx_term2 = self.motor.propellant_mass * (0.5 * self.number * self.radius**2)
+        self._propellant_I_11 = Ixx_term1 + Ixx_term2
+        self._propellant_I_22 = self._propellant_I_11
+
+        Izz_term1 = self.motor.propellant_I_33 * self.number
+        Izz_term2 = self.motor.propellant_mass * (self.number * self.radius**2)
+        self._propellant_I_33 = Izz_term1 + Izz_term2
+
+        zero_func = Function(0)
+        self._propellant_I_12 = zero_func
+        self._propellant_I_13 = zero_func
+        self._propellant_I_23 = zero_func
+
+    @property
+    def thrust(self):
+        return self._thrust
+
+    @thrust.setter
+    def thrust(self, value):
+        self._thrust = value
+
+    @property
+    def propellant_mass(self):
+        return self._propellant_mass
+
+    @propellant_mass.setter
+    def propellant_mass(self, value):
+        self._propellant_mass = value
+
+    @property
+    def propellant_initial_mass(self):
+        return self._propellant_initial_mass
+
+    @propellant_initial_mass.setter
+    def propellant_initial_mass(self, value):
+        self._propellant_initial_mass = value
+
+    @property
+    def center_of_propellant_mass(self):
+        return self._center_of_propellant_mass
+
+    @center_of_propellant_mass.setter
+    def center_of_propellant_mass(self, value):
+        self._center_of_propellant_mass = value
+
+    @property
+    def propellant_I_11(self):
+        return self._propellant_I_11
+
+    @propellant_I_11.setter
+    def propellant_I_11(self, value):
+        self._propellant_I_11 = value
+
+    @property
+    def propellant_I_22(self):
+        return self._propellant_I_22
+
+    @propellant_I_22.setter
+    def propellant_I_22(self, value):
+        self._propellant_I_22 = value
+
+    @property
+    def propellant_I_33(self):
+        return self._propellant_I_33
+
+    @propellant_I_33.setter
+    def propellant_I_33(self, value):
+        self._propellant_I_33 = value
+
+    @property
+    def propellant_I_12(self):
+        return self._propellant_I_12
+
+    @propellant_I_12.setter
+    def propellant_I_12(self, value):
+        self._propellant_I_12 = value
+
+    @property
+    def propellant_I_13(self):
+        return self._propellant_I_13
+
+    @propellant_I_13.setter
+    def propellant_I_13(self, value):
+        self._propellant_I_13 = value
+
+    @property
+    def propellant_I_23(self):
+        return self._propellant_I_23
+
+    @propellant_I_23.setter
+    def propellant_I_23(self, value):
+        self._propellant_I_23 = value
+
+    def exhaust_velocity(self, t):
+        return self.motor.exhaust_velocity(t)
+
+    def _calculate_dry_inertia(self):
+        Ixx_loc = self.motor.dry_I_11
+        Iyy_loc = self.motor.dry_I_22
+        Izz_loc = self.motor.dry_I_33
+        m_dry = self.motor.dry_mass
+
+        Izz_cluster = self.number * Izz_loc + self.number * m_dry * (self.radius**2)
+        I_transverse = self.number * Ixx_loc + (self.number / 2) * m_dry * (
+            self.radius**2
+        )
+
+        return (I_transverse, I_transverse, Izz_cluster)
+
+    def info(self):
+        print(f"Cluster Configuration:")
+        print(f" - Motors: {self.number} x {type(self.motor).__name__}")
+        print(f" - Radial Distance: {self.radius} m")
+        return self.motor.info()
+
+    def draw_cluster_layout(self, rocket_radius=None):
+        fig, ax = plt.subplots(figsize=(6, 6))
+        ax.plot(0, 0, "k+", markersize=10, label="Central axis")
+        if rocket_radius:
+            rocket_tube = plt.Circle(
+                (0, 0),
+                rocket_radius,
+                color="black",
+                fill=False,
+                linestyle="--",
+                linewidth=2,
+                label="Rocket",
+            )
+            ax.add_patch(rocket_tube)
+            limit = rocket_radius * 1.2
+        else:
+            limit = self.radius * 2
+        motor_outer_radius = self.grain_outer_radius
+        angles = np.linspace(0, 2 * np.pi, self.number, endpoint=False)
+
+        for i, angle in enumerate(angles):
+            x = self.radius * np.cos(angle)
+            y = self.radius * np.sin(angle)
+            motor_circle = plt.Circle(
+                (x, y),
+                motor_outer_radius,
+                color="red",
+                alpha=0.5,
+                label="Engine" if i == 0 else "",
+            )
+            ax.add_patch(motor_circle)
+            ax.text(
+                x,
+                y,
+                str(i + 1),
+                color="white",
+                ha="center",
+                va="center",
+                fontweight="bold",
+            )
+        ax.set_aspect("equal", "box")
+        ax.set_xlim(-limit, limit)
+        ax.set_ylim(-limit, limit)
+        ax.set_xlabel("Position X (m)")
+        ax.set_ylabel("Position Y (m)")
+        ax.set_title(f"Cluster Configuration : {self.number} engines")
+        ax.grid(True, linestyle=":", alpha=0.6)
+        ax.legend(loc="upper right")
+        plt.show()

From 53ea442250d87605de0e6febba33ebb6c9e94e1f Mon Sep 17 00:00:00 2001
From: Ayoub Abgo <ayoubdussap@gmail.com>
Date: Wed, 28 Jan 2026 16:31:34 +0100
Subject: [PATCH 02/10] feat: implement cluster_motor class with dynamic
 inertias

---
 .../integration/motors/test_cluster_motor.py  | 115 ++++++++++++++++++
 1 file changed, 115 insertions(+)
 create mode 100644 tests/integration/motors/test_cluster_motor.py

diff --git a/tests/integration/motors/test_cluster_motor.py b/tests/integration/motors/test_cluster_motor.py
new file mode 100644
index 000000000..44bdec083
--- /dev/null
+++ b/tests/integration/motors/test_cluster_motor.py
@@ -0,0 +1,115 @@
+import pytest
+import numpy as np
+from rocketpy import SolidMotor, Function
+from rocketpy.motors.cluster_motor import ClusterMotor
+
+@pytest.fixture
+def base_motor():
+    """
+    Creates a simplified SolidMotor for testing purposes.
+    Properties:
+    - Constant Thrust: 1000 N
+    - Burn time: 5 s
+    - Dry mass: 10 kg
+    - Dry Inertia: (1.0, 1.0, 0.1)
+    """
+    thrust_curve = Function(lambda t: 1000 if t < 5 else 0, "Time (s)", "Thrust (N)")
+    
+    return SolidMotor(
+        thrust_source=thrust_curve,
+        burn_time=5,
+        dry_mass=10.0,
+        dry_inertia=(1.0, 1.0, 0.1), # Ixx, Iyy, Izz
+        grain_number=1,
+        grain_density=1000,
+        grain_outer_radius=0.05,
+        grain_initial_inner_radius=0.02,
+        grain_initial_height=0.5,
+        coordinate_system_orientation="nozzle_to_combustion_chamber",
+        nozzle_radius=0.02,
+        grain_separation=0.001,
+        grains_center_of_mass_position=0.25,
+        center_of_dry_mass_position=0.25
+    )
+
+def test_cluster_initialization(base_motor):
+    """
+    Tests if the ClusterMotor initializes basic attributes correctly.
+    """
+    N = 3
+    R = 0.5
+    cluster = ClusterMotor(motor=base_motor, number=N, radius=R)
+    
+    assert cluster.number == N
+    assert cluster.radius == R
+    assert cluster.grain_outer_radius == base_motor.grain_outer_radius
+
+def test_cluster_mass_and_thrust_scaling(base_motor):
+    """
+    Tests if scalar properties (Thrust, Mass) are correctly multiplied by N.
+    """
+    N = 4
+    R = 0.2
+    cluster = ClusterMotor(motor=base_motor, number=N, radius=R)
+    
+    # 1. Check Thrust Scaling
+    # Thrust at t=1 should be N * single_motor_thrust
+    assert np.isclose(cluster.thrust(1), base_motor.thrust(1) * N)
+    
+    # 2. Check Dry Mass Scaling
+    assert np.isclose(cluster.dry_mass, base_motor.dry_mass * N)
+    
+    # 3. Check Propellant Mass Scaling
+    assert np.isclose(cluster.propellant_mass(0), base_motor.propellant_mass(0) * N)
+
+def test_cluster_dry_inertia_steiner_theorem(base_motor):
+    """
+    Tests the implementation of the Parallel Axis Theorem (Huygens-Steiner)
+    for the static (dry) mass of the cluster.
+    
+    Theoretical Formulas:
+    I_zz_cluster = N * I_zz_local + N * m * R^2
+    I_xx_cluster = N * I_xx_local + (N/2) * m * R^2  (Radial symmetry approximation)
+    """
+    N = 3
+    R = 1.0 # 1 meter radius for simpler checking
+    cluster = ClusterMotor(motor=base_motor, number=N, radius=R)
+    
+    m_dry = base_motor.dry_mass
+    Ixx_loc = base_motor.dry_I_11
+    Izz_loc = base_motor.dry_I_33
+    
+    # Expected Izz (Longitudinal / Roll)
+    expected_Izz = N * Izz_loc + N * m_dry * (R**2)
+    
+    # Expected Ixx/Iyy (Transverse / Pitch / Yaw)
+    expected_I_trans = N * Ixx_loc + (N / 2) * m_dry * (R**2)
+    
+    assert np.isclose(cluster.dry_I_33, expected_Izz)
+    assert np.isclose(cluster.dry_I_11, expected_I_trans)
+    assert np.isclose(cluster.dry_I_22, expected_I_trans)
+
+def test_cluster_propellant_inertia_dynamic(base_motor):
+    """
+    Tests if the Steiner theorem is correctly applied dynamically 
+    to the changing propellant mass over time.
+    """
+    N = 2
+    R = 0.5
+    cluster = ClusterMotor(motor=base_motor, number=N, radius=R)
+    
+    t = 0 # Check at t=0
+    
+    m_prop = base_motor.propellant_mass(t)
+    Ixx_prop_loc = base_motor.propellant_I_11(t)
+    Izz_prop_loc = base_motor.propellant_I_33(t)
+    
+    # Expected Dynamic Ixx
+    # Ixx_term1 (Local rotation) + Ixx_term2 (Parallel axis offset)
+    expected_Ixx = (Ixx_prop_loc * N) + (m_prop * 0.5 * N * R**2)
+    
+    # Expected Dynamic Izz
+    expected_Izz = (Izz_prop_loc * N) + (m_prop * N * R**2)
+    
+    assert np.isclose(cluster.propellant_I_11(t), expected_Ixx)
+    assert np.isclose(cluster.propellant_I_33(t), expected_Izz)
\ No newline at end of file

From 3a57b9c45019c1e12a41b2353791e7bcc1973525 Mon Sep 17 00:00:00 2001
From: ayoubdsp <ayoubdussap@gmail.com>
Date: Mon, 16 Mar 2026 21:28:32 +0100
Subject: [PATCH 03/10] fix: scale thrust source before Motor init

Co-authored-by: Copilot Autofix powered by AI <175728472+Copilot@users.noreply.github.com>
---
 rocketpy/motors/cluster_motor.py | 8 ++++++--
 1 file changed, 6 insertions(+), 2 deletions(-)

diff --git a/rocketpy/motors/cluster_motor.py b/rocketpy/motors/cluster_motor.py
index 0b8964d6f..4839e6a84 100644
--- a/rocketpy/motors/cluster_motor.py
+++ b/rocketpy/motors/cluster_motor.py
@@ -39,8 +39,13 @@ def __init__(self, motor, number, radius):
         self.radius = radius
         dry_inertia_cluster = self._calculate_dry_inertia()
 
+        # Use a thrust source scaled by the number of motors so that
+        # all thrust-derived quantities computed by the base Motor class
+        # correspond to the full cluster rather than a single motor.
+        scaled_thrust_source = motor.thrust * number
+
         super().__init__(
-            thrust_source=motor.thrust_source,
+            thrust_source=scaled_thrust_source,
             nozzle_radius=motor.nozzle_radius,
             burn_time=motor.burn_time,
             dry_mass=motor.dry_mass * number,
@@ -57,7 +62,6 @@ def __init__(self, motor, number, radius):
         self.grain_initial_inner_radius = motor.grain_initial_inner_radius
         self.grain_initial_height = motor.grain_initial_height
         self.grains_center_of_mass_position = motor.grains_center_of_mass_position
-        self._thrust = self.motor.thrust * self.number
         self._propellant_mass = self.motor.propellant_mass * self.number
         self._propellant_initial_mass = self.number * self.motor.propellant_initial_mass
         self._center_of_propellant_mass = self.motor.center_of_propellant_mass

From 74022e7fa84091d5172336653c025fb8540586db Mon Sep 17 00:00:00 2001
From: Ayoub Abgo <ayoubdussap@gmail.com>
Date: Mon, 16 Mar 2026 22:43:45 +0100
Subject: [PATCH 04/10] refactor: apply remaining Copilot review suggestions

---
 rocketpy/motors/cluster_motor.py              | 34 ++++++++++++++-----
 .../integration/motors/test_cluster_motor.py  |  5 ++-
 2 files changed, 29 insertions(+), 10 deletions(-)

diff --git a/rocketpy/motors/cluster_motor.py b/rocketpy/motors/cluster_motor.py
index 4839e6a84..8c241fb3b 100644
--- a/rocketpy/motors/cluster_motor.py
+++ b/rocketpy/motors/cluster_motor.py
@@ -34,6 +34,18 @@ def __init__(self, motor, number, radius):
         radius : float
             Distance from center of rocket to center of motor (m).
         """
+        if not isinstance(number, int):
+            raise TypeError(f"number must be an int, got {type(number).__name__}")
+        if number < 2:
+            raise ValueError("number must be >= 2 for a ClusterMotor")
+        if not isinstance(radius, (int, float)):
+            raise TypeError(f"radius must be a real number, got {type(radius).__name__}")
+        if radius < 0:
+            raise ValueError("radius must be non-negative")
+            
+        self.motor = motor
+        self.number = number
+        self.radius = float(radius)
         self.motor = motor
         self.number = number
         self.radius = radius
@@ -158,9 +170,9 @@ def propellant_I_23(self):
     @propellant_I_23.setter
     def propellant_I_23(self, value):
         self._propellant_I_23 = value
-
-    def exhaust_velocity(self, t):
-        return self.motor.exhaust_velocity(t)
+    @property
+    def exhaust_velocity(self):
+        return self.motor.exhaust_velocity
 
     def _calculate_dry_inertia(self):
         Ixx_loc = self.motor.dry_I_11
@@ -169,11 +181,10 @@ def _calculate_dry_inertia(self):
         m_dry = self.motor.dry_mass
 
         Izz_cluster = self.number * Izz_loc + self.number * m_dry * (self.radius**2)
-        I_transverse = self.number * Ixx_loc + (self.number / 2) * m_dry * (
-            self.radius**2
-        )
+        Ixx_cluster = self.number * Ixx_loc + (self.number / 2) * m_dry * (self.radius**2)
+        Iyy_cluster = self.number * Iyy_loc + (self.number / 2) * m_dry * (self.radius**2)
 
-        return (I_transverse, I_transverse, Izz_cluster)
+        return (Ixx_cluster, Iyy_cluster, Izz_cluster)
 
     def info(self):
         print(f"Cluster Configuration:")
@@ -181,7 +192,8 @@ def info(self):
         print(f" - Radial Distance: {self.radius} m")
         return self.motor.info()
 
-    def draw_cluster_layout(self, rocket_radius=None):
+    def draw_cluster_layout(self, rocket_radius=None,show=True):
+        """Draw the geometric layout of the clustered motors."""
         fig, ax = plt.subplots(figsize=(6, 6))
         ax.plot(0, 0, "k+", markersize=10, label="Central axis")
         if rocket_radius:
@@ -229,4 +241,8 @@ def draw_cluster_layout(self, rocket_radius=None):
         ax.set_title(f"Cluster Configuration : {self.number} engines")
         ax.grid(True, linestyle=":", alpha=0.6)
         ax.legend(loc="upper right")
-        plt.show()
+        if show:
+            plt.show()
+        return fig, ax
+    
+        
diff --git a/tests/integration/motors/test_cluster_motor.py b/tests/integration/motors/test_cluster_motor.py
index 44bdec083..7e02b2c0b 100644
--- a/tests/integration/motors/test_cluster_motor.py
+++ b/tests/integration/motors/test_cluster_motor.py
@@ -46,7 +46,7 @@ def test_cluster_initialization(base_motor):
 
 def test_cluster_mass_and_thrust_scaling(base_motor):
     """
-    Tests if scalar properties (Thrust, Mass) are correctly multiplied by N.
+    Tests if scalar and derived properties are correctly multiplied by N and that functional properties preserve their Function behavior
     """
     N = 4
     R = 0.2
@@ -61,6 +61,9 @@ def test_cluster_mass_and_thrust_scaling(base_motor):
     
     # 3. Check Propellant Mass Scaling
     assert np.isclose(cluster.propellant_mass(0), base_motor.propellant_mass(0) * N)
+    assert np.isclose(cluster.total_impulse, base_motor.total_impulse * N)
+    assert np.isclose(cluster.average_thrust, base_motor.average_thrust * N)
+
 
 def test_cluster_dry_inertia_steiner_theorem(base_motor):
     """

From f9727e25872dc98e5e5ba09f557796e5bee325a8 Mon Sep 17 00:00:00 2001
From: Ayoub Abgo <ayoubdussap@gmail.com>
Date: Thu, 19 Mar 2026 22:05:55 +0100
Subject: [PATCH 05/10] style: run black formatter and fix pylint warnings

---
 rocketpy/motors/cluster_motor.py              | 25 ++++++----
 .../integration/motors/test_cluster_motor.py  | 46 ++++++++++---------
 2 files changed, 40 insertions(+), 31 deletions(-)

diff --git a/rocketpy/motors/cluster_motor.py b/rocketpy/motors/cluster_motor.py
index 8c241fb3b..9b2c90a06 100644
--- a/rocketpy/motors/cluster_motor.py
+++ b/rocketpy/motors/cluster_motor.py
@@ -39,10 +39,12 @@ def __init__(self, motor, number, radius):
         if number < 2:
             raise ValueError("number must be >= 2 for a ClusterMotor")
         if not isinstance(radius, (int, float)):
-            raise TypeError(f"radius must be a real number, got {type(radius).__name__}")
+            raise TypeError(
+                f"radius must be a real number, got {type(radius).__name__}"
+            )
         if radius < 0:
             raise ValueError("radius must be non-negative")
-            
+
         self.motor = motor
         self.number = number
         self.radius = float(radius)
@@ -170,6 +172,7 @@ def propellant_I_23(self):
     @propellant_I_23.setter
     def propellant_I_23(self, value):
         self._propellant_I_23 = value
+
     @property
     def exhaust_velocity(self):
         return self.motor.exhaust_velocity
@@ -181,18 +184,22 @@ def _calculate_dry_inertia(self):
         m_dry = self.motor.dry_mass
 
         Izz_cluster = self.number * Izz_loc + self.number * m_dry * (self.radius**2)
-        Ixx_cluster = self.number * Ixx_loc + (self.number / 2) * m_dry * (self.radius**2)
-        Iyy_cluster = self.number * Iyy_loc + (self.number / 2) * m_dry * (self.radius**2)
+        Ixx_cluster = self.number * Ixx_loc + (self.number / 2) * m_dry * (
+            self.radius**2
+        )
+        Iyy_cluster = self.number * Iyy_loc + (self.number / 2) * m_dry * (
+            self.radius**2
+        )
 
         return (Ixx_cluster, Iyy_cluster, Izz_cluster)
 
-    def info(self):
-        print(f"Cluster Configuration:")
+    def info(self, *args, **kwargs):
+        print("Cluster Configuration:")
         print(f" - Motors: {self.number} x {type(self.motor).__name__}")
         print(f" - Radial Distance: {self.radius} m")
-        return self.motor.info()
+        return self.motor.info(*args, **kwargs)
 
-    def draw_cluster_layout(self, rocket_radius=None,show=True):
+    def draw_cluster_layout(self, rocket_radius=None, show=True):
         """Draw the geometric layout of the clustered motors."""
         fig, ax = plt.subplots(figsize=(6, 6))
         ax.plot(0, 0, "k+", markersize=10, label="Central axis")
@@ -244,5 +251,3 @@ def draw_cluster_layout(self, rocket_radius=None,show=True):
         if show:
             plt.show()
         return fig, ax
-    
-        
diff --git a/tests/integration/motors/test_cluster_motor.py b/tests/integration/motors/test_cluster_motor.py
index 7e02b2c0b..537092d40 100644
--- a/tests/integration/motors/test_cluster_motor.py
+++ b/tests/integration/motors/test_cluster_motor.py
@@ -3,6 +3,7 @@
 from rocketpy import SolidMotor, Function
 from rocketpy.motors.cluster_motor import ClusterMotor
 
+
 @pytest.fixture
 def base_motor():
     """
@@ -14,12 +15,12 @@ def base_motor():
     - Dry Inertia: (1.0, 1.0, 0.1)
     """
     thrust_curve = Function(lambda t: 1000 if t < 5 else 0, "Time (s)", "Thrust (N)")
-    
+
     return SolidMotor(
         thrust_source=thrust_curve,
         burn_time=5,
         dry_mass=10.0,
-        dry_inertia=(1.0, 1.0, 0.1), # Ixx, Iyy, Izz
+        dry_inertia=(1.0, 1.0, 0.1),  # Ixx, Iyy, Izz
         grain_number=1,
         grain_density=1000,
         grain_outer_radius=0.05,
@@ -29,9 +30,10 @@ def base_motor():
         nozzle_radius=0.02,
         grain_separation=0.001,
         grains_center_of_mass_position=0.25,
-        center_of_dry_mass_position=0.25
+        center_of_dry_mass_position=0.25,
     )
 
+
 def test_cluster_initialization(base_motor):
     """
     Tests if the ClusterMotor initializes basic attributes correctly.
@@ -39,11 +41,12 @@ def test_cluster_initialization(base_motor):
     N = 3
     R = 0.5
     cluster = ClusterMotor(motor=base_motor, number=N, radius=R)
-    
+
     assert cluster.number == N
     assert cluster.radius == R
     assert cluster.grain_outer_radius == base_motor.grain_outer_radius
 
+
 def test_cluster_mass_and_thrust_scaling(base_motor):
     """
     Tests if scalar and derived properties are correctly multiplied by N and that functional properties preserve their Function behavior
@@ -51,14 +54,14 @@ def test_cluster_mass_and_thrust_scaling(base_motor):
     N = 4
     R = 0.2
     cluster = ClusterMotor(motor=base_motor, number=N, radius=R)
-    
+
     # 1. Check Thrust Scaling
     # Thrust at t=1 should be N * single_motor_thrust
     assert np.isclose(cluster.thrust(1), base_motor.thrust(1) * N)
-    
+
     # 2. Check Dry Mass Scaling
     assert np.isclose(cluster.dry_mass, base_motor.dry_mass * N)
-    
+
     # 3. Check Propellant Mass Scaling
     assert np.isclose(cluster.propellant_mass(0), base_motor.propellant_mass(0) * N)
     assert np.isclose(cluster.total_impulse, base_motor.total_impulse * N)
@@ -69,50 +72,51 @@ def test_cluster_dry_inertia_steiner_theorem(base_motor):
     """
     Tests the implementation of the Parallel Axis Theorem (Huygens-Steiner)
     for the static (dry) mass of the cluster.
-    
+
     Theoretical Formulas:
     I_zz_cluster = N * I_zz_local + N * m * R^2
     I_xx_cluster = N * I_xx_local + (N/2) * m * R^2  (Radial symmetry approximation)
     """
     N = 3
-    R = 1.0 # 1 meter radius for simpler checking
+    R = 1.0  # 1 meter radius for simpler checking
     cluster = ClusterMotor(motor=base_motor, number=N, radius=R)
-    
+
     m_dry = base_motor.dry_mass
     Ixx_loc = base_motor.dry_I_11
     Izz_loc = base_motor.dry_I_33
-    
+
     # Expected Izz (Longitudinal / Roll)
     expected_Izz = N * Izz_loc + N * m_dry * (R**2)
-    
+
     # Expected Ixx/Iyy (Transverse / Pitch / Yaw)
     expected_I_trans = N * Ixx_loc + (N / 2) * m_dry * (R**2)
-    
+
     assert np.isclose(cluster.dry_I_33, expected_Izz)
     assert np.isclose(cluster.dry_I_11, expected_I_trans)
     assert np.isclose(cluster.dry_I_22, expected_I_trans)
 
+
 def test_cluster_propellant_inertia_dynamic(base_motor):
     """
-    Tests if the Steiner theorem is correctly applied dynamically 
+    Tests if the Steiner theorem is correctly applied dynamically
     to the changing propellant mass over time.
     """
     N = 2
     R = 0.5
     cluster = ClusterMotor(motor=base_motor, number=N, radius=R)
-    
-    t = 0 # Check at t=0
-    
+
+    t = 0  # Check at t=0
+
     m_prop = base_motor.propellant_mass(t)
     Ixx_prop_loc = base_motor.propellant_I_11(t)
     Izz_prop_loc = base_motor.propellant_I_33(t)
-    
+
     # Expected Dynamic Ixx
     # Ixx_term1 (Local rotation) + Ixx_term2 (Parallel axis offset)
     expected_Ixx = (Ixx_prop_loc * N) + (m_prop * 0.5 * N * R**2)
-    
+
     # Expected Dynamic Izz
     expected_Izz = (Izz_prop_loc * N) + (m_prop * N * R**2)
-    
+
     assert np.isclose(cluster.propellant_I_11(t), expected_Ixx)
-    assert np.isclose(cluster.propellant_I_33(t), expected_Izz)
\ No newline at end of file
+    assert np.isclose(cluster.propellant_I_33(t), expected_Izz)

From 48d33d695e38395bda73f20bf89c5eb97931da36 Mon Sep 17 00:00:00 2001
From: Ayoub Abgo <ayoubdussap@gmail.com>
Date: Thu, 19 Mar 2026 22:13:07 +0100
Subject: [PATCH 06/10] test: add coverage for validation, setters, and display
 methods

---
 .../integration/motors/test_cluster_motor.py  | 27 +++++++++++++++++++
 1 file changed, 27 insertions(+)

diff --git a/tests/integration/motors/test_cluster_motor.py b/tests/integration/motors/test_cluster_motor.py
index 537092d40..4c7ce20d0 100644
--- a/tests/integration/motors/test_cluster_motor.py
+++ b/tests/integration/motors/test_cluster_motor.py
@@ -120,3 +120,30 @@ def test_cluster_propellant_inertia_dynamic(base_motor):
 
     assert np.isclose(cluster.propellant_I_11(t), expected_Ixx)
     assert np.isclose(cluster.propellant_I_33(t), expected_Izz)
+
+def test_cluster_invalid_inputs(base_motor):
+    """Tests if the validation raises errors for bad inputs."""
+    with pytest.raises(ValueError):
+        ClusterMotor(motor=base_motor, number=1, radius=0.5) # N < 2
+    with pytest.raises(ValueError):
+        ClusterMotor(motor=base_motor, number=2, radius=-1.0) # Radius < 0
+    with pytest.raises(TypeError):
+        ClusterMotor(motor=base_motor, number="two", radius=0.5) # N is string
+
+def test_cluster_methods_and_setters(base_motor):
+    """Touches the display methods and setters to ensure coverage."""
+    cluster = ClusterMotor(motor=base_motor, number=2, radius=0.5)
+    
+    # 1. Touch the info method
+    cluster.info()
+    
+    # 2. Touch the draw method (without showing the plot to avoid blocking tests)
+    cluster.draw_cluster_layout(show=False)
+    cluster.draw_cluster_layout(rocket_radius=0.1, show=False)
+    
+    # 3. Touch a few setters
+    cluster.propellant_mass = 50.0
+    assert cluster.propellant_mass == 50.0
+    
+    cluster.propellant_I_11 = 2.0
+    assert cluster.propellant_I_11 == 2.0
\ No newline at end of file

From c8739011884e79959114166937e42f6384c9f039 Mon Sep 17 00:00:00 2001
From: Ayoub Abgo <ayoubdussap@gmail.com>
Date: Sat, 21 Mar 2026 17:58:38 +0100
Subject: [PATCH 07/10] refactor: resolve pylint too-many-statements and number
 of parameters

---
 rocketpy/motors/cluster_motor.py | 52 +++++++++++++++++++-------------
 1 file changed, 31 insertions(+), 21 deletions(-)

diff --git a/rocketpy/motors/cluster_motor.py b/rocketpy/motors/cluster_motor.py
index 9b2c90a06..6a71dd40c 100644
--- a/rocketpy/motors/cluster_motor.py
+++ b/rocketpy/motors/cluster_motor.py
@@ -1,3 +1,4 @@
+# pylint: disable=invalid-name
 import matplotlib.pyplot as plt
 import numpy as np
 from rocketpy import Function
@@ -48,9 +49,6 @@ def __init__(self, motor, number, radius):
         self.motor = motor
         self.number = number
         self.radius = float(radius)
-        self.motor = motor
-        self.number = number
-        self.radius = radius
         dry_inertia_cluster = self._calculate_dry_inertia()
 
         # Use a thrust source scaled by the number of motors so that
@@ -69,16 +67,14 @@ def __init__(self, motor, number, radius):
             interpolation_method="linear",
         )
 
-        self.throat_radius = motor.throat_radius
-        self.grain_number = motor.grain_number
-        self.grain_density = motor.grain_density
-        self.grain_outer_radius = motor.grain_outer_radius
-        self.grain_initial_inner_radius = motor.grain_initial_inner_radius
-        self.grain_initial_height = motor.grain_initial_height
-        self.grains_center_of_mass_position = motor.grains_center_of_mass_position
+        self._setup_grain_properties()
         self._propellant_mass = self.motor.propellant_mass * self.number
         self._propellant_initial_mass = self.number * self.motor.propellant_initial_mass
         self._center_of_propellant_mass = self.motor.center_of_propellant_mass
+        self._evaluate_propellant_inertia()
+
+    def _evaluate_propellant_inertia(self):
+        """Calculates the dynamic inertia of the propellant using Steiner's theorem."""
         Ixx_term1 = self.motor.propellant_I_11 * self.number
         Ixx_term2 = self.motor.propellant_mass * (0.5 * self.number * self.radius**2)
         self._propellant_I_11 = Ixx_term1 + Ixx_term2
@@ -93,6 +89,16 @@ def __init__(self, motor, number, radius):
         self._propellant_I_13 = zero_func
         self._propellant_I_23 = zero_func
 
+    def _setup_grain_properties(self):
+        """Copies the grain properties from the base motor."""
+        self.throat_radius = self.motor.throat_radius
+        self.grain_number = self.motor.grain_number
+        self.grain_density = self.motor.grain_density
+        self.grain_outer_radius = self.motor.grain_outer_radius
+        self.grain_initial_inner_radius = self.motor.grain_initial_inner_radius
+        self.grain_initial_height = self.motor.grain_initial_height
+        self.grains_center_of_mass_position = self.motor.grains_center_of_mass_position
+
     @property
     def thrust(self):
         return self._thrust
@@ -217,6 +223,21 @@ def draw_cluster_layout(self, rocket_radius=None, show=True):
             limit = rocket_radius * 1.2
         else:
             limit = self.radius * 2
+        self._draw_engines(ax)
+        ax.set_aspect("equal", "box")
+        ax.set_xlim(-limit, limit)
+        ax.set_ylim(-limit, limit)
+        ax.set_xlabel("Position X (m)")
+        ax.set_ylabel("Position Y (m)")
+        ax.set_title(f"Cluster Configuration : {self.number} engines")
+        ax.grid(True, linestyle=":", alpha=0.6)
+        ax.legend(loc="upper right")
+        if show:
+            plt.show()
+        return fig, ax
+
+    def _draw_engines(self, ax):
+        """Draws the individual engines of the cluster."""
         motor_outer_radius = self.grain_outer_radius
         angles = np.linspace(0, 2 * np.pi, self.number, endpoint=False)
 
@@ -240,14 +261,3 @@ def draw_cluster_layout(self, rocket_radius=None, show=True):
                 va="center",
                 fontweight="bold",
             )
-        ax.set_aspect("equal", "box")
-        ax.set_xlim(-limit, limit)
-        ax.set_ylim(-limit, limit)
-        ax.set_xlabel("Position X (m)")
-        ax.set_ylabel("Position Y (m)")
-        ax.set_title(f"Cluster Configuration : {self.number} engines")
-        ax.grid(True, linestyle=":", alpha=0.6)
-        ax.legend(loc="upper right")
-        if show:
-            plt.show()
-        return fig, ax

From eb34298c9cd477ebd9a0c513e11839a7c154d171 Mon Sep 17 00:00:00 2001
From: Ayoub Abgo <ayoubdussap@gmail.com>
Date: Sat, 21 Mar 2026 23:00:23 +0100
Subject: [PATCH 08/10] feat(plots): add visual rendering and layout support
 for ClusterMotor

---
 rocketpy/plots/rocket_plots.py | 91 ++++++++++++++++++++++------------
 1 file changed, 59 insertions(+), 32 deletions(-)

diff --git a/rocketpy/plots/rocket_plots.py b/rocketpy/plots/rocket_plots.py
index 0f087a5d6..ba58ebc42 100644
--- a/rocketpy/plots/rocket_plots.py
+++ b/rocketpy/plots/rocket_plots.py
@@ -1,7 +1,7 @@
 import matplotlib.pyplot as plt
 import numpy as np
 
-from rocketpy.motors import EmptyMotor, HybridMotor, LiquidMotor, SolidMotor
+from rocketpy.motors import HybridMotor, LiquidMotor, SolidMotor
 from rocketpy.rocket.aero_surface import Fins, NoseCone, Tail
 from rocketpy.rocket.aero_surface.generic_surface import GenericSurface
 
@@ -414,57 +414,84 @@ def _draw_tubes(self, ax, drawn_surfaces, vis_args):
     def _draw_motor(self, last_radius, last_x, ax, vis_args):
         """Draws the motor from motor patches"""
         total_csys = self.rocket._csys * self.rocket.motor._csys
+        is_cluster = hasattr(self.rocket.motor, "number")
+        base_motor = self.rocket.motor.motor if is_cluster else self.rocket.motor
+
+        if is_cluster:
+            angles = np.linspace(0, 2 * np.pi, self.rocket.motor.number, endpoint=False)
+            y_offsets = self.rocket.motor.radius * np.cos(angles)
+        else:
+            y_offsets = [0]
         nozzle_position = (
-            self.rocket.motor_position + self.rocket.motor.nozzle_position * total_csys
+            self.rocket.motor_position + base_motor.nozzle_position * total_csys
         )
-
         # Get motor patches translated to the correct position
         motor_patches = self._generate_motor_patches(total_csys, ax)
-
         # Draw patches
-        if not isinstance(self.rocket.motor, EmptyMotor):
-            # Add nozzle last so it is in front of the other patches
-            nozzle = self.rocket.motor.plots._generate_nozzle(
-                translate=(nozzle_position, 0), csys=self.rocket._csys
-            )
-            motor_patches += [nozzle]
+        if type(self.rocket.motor).__name__ != "EmptyMotor":
+            for y_off in y_offsets:
+                nozzle = base_motor.plots._generate_nozzle(
+                    translate=(nozzle_position, y_off), csys=self.rocket._csys
+                )
+                if y_off != y_offsets[0]:
+                    nozzle.set_label("_nolegend_")
+                motor_patches.append(nozzle)
 
-            outline = self.rocket.motor.plots._generate_motor_region(
+            outline = base_motor.plots._generate_motor_region(
                 list_of_patches=motor_patches
             )
-            # add outline first so it is behind the other patches
-            ax.add_patch(outline)
+            if not is_cluster:
+                ax.add_patch(outline)
+
             for patch in motor_patches:
+                if is_cluster:
+                    patch.set_alpha(0.6)
                 ax.add_patch(patch)
-
         self._draw_nozzle_tube(last_radius, last_x, nozzle_position, ax, vis_args)
 
-    def _generate_motor_patches(self, total_csys, ax):  # pylint: disable=unused-argument
+    def _generate_motor_patches(
+        self, total_csys, ax
+    ):  # pylint: disable=unused-argument
         """Generates motor patches for drawing"""
         motor_patches = []
 
-        if isinstance(self.rocket.motor, SolidMotor):
+        is_cluster = hasattr(self.rocket.motor, "number")
+        base_motor = self.rocket.motor.motor if is_cluster else self.rocket.motor
+
+        if isinstance(base_motor, SolidMotor):
+            y_offsets = (
+                self.rocket.motor.radius
+                * np.cos(
+                    np.linspace(0, 2 * np.pi, self.rocket.motor.number, endpoint=False)
+                )
+                if is_cluster
+                else [0]
+            )
             grains_cm_position = (
                 self.rocket.motor_position
-                + self.rocket.motor.grains_center_of_mass_position * total_csys
-            )
-            ax.scatter(
-                grains_cm_position,
-                0,
-                color="brown",
-                label="Grains Center of Mass",
-                s=8,
-                zorder=10,
+                + base_motor.grains_center_of_mass_position * total_csys
             )
+            for y_off in y_offsets:
+                ax.scatter(
+                    grains_cm_position,
+                    y_off,
+                    color="brown",
+                    label="Grains Center of Mass" if y_off == y_offsets[0] else "",
+                    s=8,
+                    zorder=10,
+                )
 
-            chamber = self.rocket.motor.plots._generate_combustion_chamber(
-                translate=(grains_cm_position, 0), label=None
-            )
-            grains = self.rocket.motor.plots._generate_grains(
-                translate=(grains_cm_position, 0)
-            )
+                chamber = base_motor.plots._generate_combustion_chamber(
+                    translate=(grains_cm_position, y_off), label=None
+                )
+                grains = base_motor.plots._generate_grains(
+                    translate=(grains_cm_position, y_off)
+                )
+                if y_off != y_offsets[0]:
+                    for grain in grains:
+                        grain.set_label("_nolegend_")
 
-            motor_patches += [chamber, *grains]
+                motor_patches += [chamber, *grains]
 
         elif isinstance(self.rocket.motor, HybridMotor):
             grains_cm_position = (

From 1cb957d8ff1777025f1e6f8c0285148cbde1a72d Mon Sep 17 00:00:00 2001
From: Ayoub Abgo <ayoubdussap@gmail.com>
Date: Thu, 26 Mar 2026 18:46:03 +0100
Subject: [PATCH 09/10] fix: add missing numpy import for environment
 integration tests (that appeared through th merge

---
 tests/integration/environment/test_environment.py | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/tests/integration/environment/test_environment.py b/tests/integration/environment/test_environment.py
index d919c535d..be8826d0c 100644
--- a/tests/integration/environment/test_environment.py
+++ b/tests/integration/environment/test_environment.py
@@ -1,7 +1,7 @@
 import time
 from datetime import date, datetime, timezone
 from unittest.mock import patch
-
+import numpy as np
 import pytest
 
 

From 91e5a2e0529be8ae3c1010b53f13e964ac30217a Mon Sep 17 00:00:00 2001
From: Ayoub Abgo <ayoubdussap@gmail.com>
Date: Thu, 26 Mar 2026 20:49:22 +0100
Subject: [PATCH 10/10] style: reformat code with ruff

---
 rocketpy/plots/rocket_plots.py                 |  4 +---
 tests/integration/motors/test_cluster_motor.py | 18 ++++++++++--------
 2 files changed, 11 insertions(+), 11 deletions(-)

diff --git a/rocketpy/plots/rocket_plots.py b/rocketpy/plots/rocket_plots.py
index a01a62dfa..2f1f959fa 100644
--- a/rocketpy/plots/rocket_plots.py
+++ b/rocketpy/plots/rocket_plots.py
@@ -455,9 +455,7 @@ def _draw_motor(self, last_radius, last_x, ax, vis_args):
                 ax.add_patch(patch)
         self._draw_nozzle_tube(last_radius, last_x, nozzle_position, ax, vis_args)
 
-    def _generate_motor_patches(
-        self, total_csys, ax
-    ):  # pylint: disable=unused-argument
+    def _generate_motor_patches(self, total_csys, ax):  # pylint: disable=unused-argument
         """Generates motor patches for drawing"""
         motor_patches = []
 
diff --git a/tests/integration/motors/test_cluster_motor.py b/tests/integration/motors/test_cluster_motor.py
index 4c7ce20d0..e5515b5c5 100644
--- a/tests/integration/motors/test_cluster_motor.py
+++ b/tests/integration/motors/test_cluster_motor.py
@@ -121,29 +121,31 @@ def test_cluster_propellant_inertia_dynamic(base_motor):
     assert np.isclose(cluster.propellant_I_11(t), expected_Ixx)
     assert np.isclose(cluster.propellant_I_33(t), expected_Izz)
 
+
 def test_cluster_invalid_inputs(base_motor):
     """Tests if the validation raises errors for bad inputs."""
     with pytest.raises(ValueError):
-        ClusterMotor(motor=base_motor, number=1, radius=0.5) # N < 2
+        ClusterMotor(motor=base_motor, number=1, radius=0.5)  # N < 2
     with pytest.raises(ValueError):
-        ClusterMotor(motor=base_motor, number=2, radius=-1.0) # Radius < 0
+        ClusterMotor(motor=base_motor, number=2, radius=-1.0)  # Radius < 0
     with pytest.raises(TypeError):
-        ClusterMotor(motor=base_motor, number="two", radius=0.5) # N is string
+        ClusterMotor(motor=base_motor, number="two", radius=0.5)  # N is string
+
 
 def test_cluster_methods_and_setters(base_motor):
     """Touches the display methods and setters to ensure coverage."""
     cluster = ClusterMotor(motor=base_motor, number=2, radius=0.5)
-    
+
     # 1. Touch the info method
     cluster.info()
-    
+
     # 2. Touch the draw method (without showing the plot to avoid blocking tests)
     cluster.draw_cluster_layout(show=False)
     cluster.draw_cluster_layout(rocket_radius=0.1, show=False)
-    
+
     # 3. Touch a few setters
     cluster.propellant_mass = 50.0
     assert cluster.propellant_mass == 50.0
-    
+
     cluster.propellant_I_11 = 2.0
-    assert cluster.propellant_I_11 == 2.0
\ No newline at end of file
+    assert cluster.propellant_I_11 == 2.0