BUG: Zero mass flow rates being ignored on LiquidMotors.

Author: phmbressan

CHANGELOG.md

@@ -43,6 +43,7 @@ Attention: The newest changes should be on top -->
 
 ### Fixed
 
+- BUG: Zero Mass Flow Rate in Liquid Motors breaks Exhaust Velocity [#677](https://github.com/RocketPy-Team/RocketPy/pull/677)
 - DOC: Fix documentation dependencies [#651](https://github.com/RocketPy-Team/RocketPy/pull/651)
 - DOC: Fix documentation warnings [#645](https://github.com/RocketPy-Team/RocketPy/pull/645)
 - BUG: Rotational EOMs Not Relative To CDM [#674](https://github.com/RocketPy-Team/RocketPy/pull/674)

rocketpy/motors/liquid_motor.py

@@ -266,16 +266,16 @@ def exhaust_velocity(self):
         """
         times, thrusts = self.thrust.source[:, 0], self.thrust.source[:, 1]
         mass_flow_rates = self.mass_flow_rate(times)
+        exhaust_velocity = np.zeros_like(mass_flow_rates)
 
         # Compute exhaust velocity only for non-zero mass flow rates
         valid_indices = mass_flow_rates != 0
-        valid_times = times[valid_indices]
-        valid_thrusts = thrusts[valid_indices]
-        valid_mass_flow_rates = mass_flow_rates[valid_indices]
 
-        ext_vel = -valid_thrusts / valid_mass_flow_rates
+        exhaust_velocity[valid_indices] = (
+            -thrusts[valid_indices] / mass_flow_rates[valid_indices]
+        )
 
-        return np.column_stack([valid_times, ext_vel])
+        return np.column_stack([times, exhaust_velocity])
 
     @funcify_method("Time (s)", "Propellant Mass (kg)")
     def propellant_mass(self):

rocketpy/motors/motor.py

@@ -483,7 +483,8 @@ def total_mass_flow_rate(self):
         rate should not be greater than `total_mass_flow_rate`, otherwise the
         grains mass flow rate will be negative, losing physical meaning.
         """
-        return -1 * self.thrust / self.exhaust_velocity
+        average_exhaust_velocity = self.total_impulse / self.propellant_initial_mass
+        return self.thrust / -average_exhaust_velocity
 
     @property
     @abstractmethod