feat!: add minimum maneuver duration constraint (#541)

* feat: add minimum maneuver duration constraint

* fix: print test

* feat: address feedback

Author: vishwa2710

bindings/python/src/OpenSpaceToolkitAstrodynamicsPy/Trajectory/Sequence.cpp

@@ -209,16 +209,18 @@ inline void OpenSpaceToolkitAstrodynamicsPy_Trajectory_Sequence(pybind11::module
                     const NumericalSolver&,
                     const Array<Shared<Dynamics>>&,
                     const Duration&,
+                    const Duration&,
                     const Size&>(),
                 R"doc(
                     Construct a new `Sequence` object.
 
                     Args:
-                    segments (list[Segment], optional): The segments.
-                    numerical_solver (NumericalSolver, optional): The numerical solver.
-                    dynamics (list[Dynamics], optional): The dynamics.
-                    maximum_propagation_duration (Duration, optional): The maximum propagation duration.
-                    verbosity (int, optional): The verbosity level.
+                    segments (list[Segment], optional): The segments. Defaults to an empty list.
+                    numerical_solver (NumericalSolver, optional): The numerical solver. Defaults to the default conditional numerical solver.
+                    dynamics (list[Dynamics], optional): The dynamics. Defaults to an empty list.
+                    maximum_propagation_duration (Duration, optional): The maximum propagation duration. Defaults to 30 days.
+                    minimum_maneuver_duration (Duration, optional): The minimum maneuver duration. Defaults to Undefined. If defined, maneuvers less than this duration will be skipped.
+                    verbosity (int, optional): The verbosity level. Defaults to 1.
 
                     Returns:
                         Sequence: The new `Sequence` object.
@@ -229,7 +231,8 @@ inline void OpenSpaceToolkitAstrodynamicsPy_Trajectory_Sequence(pybind11::module
                     "numerical_solver", NumericalSolver::DefaultConditional(), "NumericalSolver.default_conditional()"
                 ),
                 arg_v("dynamics", Array<Shared<Dynamics>>::Empty(), "[]"),
-                arg_v("maximum_propagation_duration", Duration::Days(30.0), "duration.days(30.0)"),
+                arg_v("maximum_propagation_duration", Duration::Days(30.0), "Duration.days(30.0)"),
+                arg_v("minimum_maneuver_duration", Duration::Undefined(), "Duration.undefined()"),
                 arg("verbosity") = 1
             )
 
@@ -280,6 +283,17 @@ inline void OpenSpaceToolkitAstrodynamicsPy_Trajectory_Sequence(pybind11::module
 
                 )doc"
             )
+            .def(
+                "get_minimum_maneuver_duration",
+                &Sequence::getMinimumManeuverDuration,
+                R"doc(
+                    Get the minimum maneuver duration.
+
+                    Returns:
+                        Duration: The minimum maneuver duration.
+
+                )doc"
+            )
 
             .def(
                 "add_segment",

bindings/python/test/trajectory/test_sequence.py

@@ -156,17 +156,7 @@ def state(
 
     return State(
         instant,
-        [
-            717094.039086306,
-            -6872433.2241124,
-            46175.9696673281,
-            -970.650826004612,
-            -45.4598114773158,
-            7529.82424886455,
-            dry_mass.in_kilograms() + wet_mass.in_kilograms(),
-            cross_sectional_surface_area,
-            drag_coefficient,
-        ],
+        coordinates,
         frame,
         coordinate_broker,
     )
@@ -351,6 +341,28 @@ def sequence_solution(
     )
 
 
+@pytest.fixture
+def minimum_maneuver_duration():
+    return Duration.minutes(1.0)
+
+
+@pytest.fixture
+def sequence_with_minimum_maneuver_duration(
+    segments: list[Segment],
+    numerical_solver: NumericalSolver,
+    dynamics: list,
+    maximum_propagation_duration: Duration,
+    minimum_maneuver_duration: Duration,
+):
+    return Sequence(
+        segments=segments,
+        dynamics=dynamics,
+        numerical_solver=numerical_solver,
+        maximum_propagation_duration=maximum_propagation_duration,
+        minimum_maneuver_duration=minimum_maneuver_duration,
+    )
+
+
 class TestSequenceSolution:
     def test_properties(
         self,
@@ -421,6 +433,16 @@ def test_get_maximum_propagation_duration(
     ):
         assert sequence.get_maximum_propagation_duration() == maximum_propagation_duration
 
+    def test_get_minimum_maneuver_duration(
+        self,
+        sequence_with_minimum_maneuver_duration: Sequence,
+        minimum_maneuver_duration: Duration,
+    ):
+        assert (
+            sequence_with_minimum_maneuver_duration.get_minimum_maneuver_duration()
+            == minimum_maneuver_duration
+        )
+
     def test_add_segment(
         self,
         sequence: Sequence,

include/OpenSpaceToolkit/Astrodynamics/Trajectory/Segment.hpp

@@ -68,14 +68,14 @@ class Segment
         ///
         /// @param aName Name of the segment
         /// @param aDynamicsArray Array of dynamics
-        /// @param aStates Array of states for the segment
+        /// @param aStateArray Array of states for the segment
         /// @param aConditionIsSatisfied True if the event condition is satisfied
         /// @param aSegmentType Type of segment
         /// @return An instance of Solution
         Solution(
             const String& aName,
             const Array<Shared<Dynamics>>& aDynamicsArray,
-            const Array<State>& aStates,
+            const Array<State>& aStateArray,
             const bool& aConditionIsSatisfied,
             const Segment::Type& aSegmentType
         );

include/OpenSpaceToolkit/Astrodynamics/Trajectory/Sequence.hpp

@@ -125,25 +125,29 @@ class Sequence
     ///                  const Array<Shared<Dynamics>> dynamicsArray =
     ///                  {std::make_shared<CentralBodyGravity>(Earth::GravitationalParameter())};
     ///                  const Duration maximumPropagationDuration = Duration::Days(7.0);
+    ///                  const Duration minimumManeuverDuration = Duration::Zero();
     ///                  const Size verbosity = 0;
     ///
     ///                  Sequence sequence = {segmentArray, numericalSolver, dynamicsArray,
-    ///                  maximumPropagationDuration, verbosity};
+    ///                  maximumPropagationDuration, minimumManeuverDuration, verbosity};
     ///
     /// @endcode
     ///
     /// @param aSegmentArray An array of segments. Defaults to empty.
-    /// @param aNumericalSolver A Numerical Solver. Defaults to Undefined.
+    /// @param aNumericalSolver A Numerical Solver. Defaults to NumericalSolver::DefaultConditional().
     /// @param aDynamicsArray An array of shared dynamics. Defaults to empty.
-    /// @param segmentPropagationDurationLimit Maximum duration for propagation. Defaults to 7.0
+    /// @param aSegmentPropagationDurationLimit Maximum duration for propagation. Defaults to 30.0
     /// days.
-    /// @param verbosity Verbosity level for the solver [0 (low) - 5 (high)]. Defaults to 0.
+    /// @param aMinimumManeuverDuration Minimum duration for maneuver, maneuvers less than this duration
+    /// will be skipped. Defaults to Undefined.
+    /// @param aVerbosityLevel Verbosity level for the solver [0 (low) - 5 (high)]. Defaults to 0.
     Sequence(
         const Array<Segment>& aSegmentArray = Array<Segment>::Empty(),
-        const NumericalSolver& aNumericalSolver = NumericalSolver::Undefined(),
+        const NumericalSolver& aNumericalSolver = NumericalSolver::DefaultConditional(),
         const Array<Shared<Dynamics>>& aDynamicsArray = Array<Shared<Dynamics>>::Empty(),
-        const Duration& segmentPropagationDurationLimit = Duration::Days(7.0),
-        const Size& verbosity = 0
+        const Duration& aSegmentPropagationDurationLimit = Duration::Days(30.0),
+        const Duration& aMinimumManeuverDuration = Duration::Undefined(),
+        const Size& aVerbosityLevel = 0
     );
 
     /// @brief Output stream operator.
@@ -173,6 +177,11 @@ class Sequence
     /// @return Maximum propagation duration.
     Duration getMaximumPropagationDuration() const;
 
+    /// @brief Get minimum maneuver duration.
+    ///
+    /// @return Minimum maneuver duration.
+    Duration getMinimumManeuverDuration() const;
+
     /// @brief Add a trajectory segment.
     ///
     /// @param aTrajectorySegment A trajectory segment.
@@ -224,6 +233,7 @@ class Sequence
     NumericalSolver numericalSolver_;
     Array<Shared<Dynamics>> dynamics_;
     Duration segmentPropagationDurationLimit_;
+    Duration minimumManeuverDuration_;
 };
 
 }  // namespace trajectory

src/OpenSpaceToolkit/Astrodynamics/Trajectory/Segment.cpp

@@ -35,13 +35,13 @@ using ostk::astrodynamics::trajectory::StateBuilder;
 Segment::Solution::Solution(
     const String& aName,
     const Array<Shared<Dynamics>>& aDynamicsArray,
-    const Array<State>& aStates,
+    const Array<State>& aStateArray,
     const bool& aConditionIsSatisfied,
     const Segment::Type& aSegmentType
 )
     : name(aName),
       dynamics(aDynamicsArray),
-      states(aStates),
+      states(aStateArray),
       conditionIsSatisfied(aConditionIsSatisfied),
       segmentType(aSegmentType)
 {
@@ -505,10 +505,21 @@ Segment::Solution Segment::solve(const State& aState, const Duration& maximumPro
         aState, aState.accessInstant() + maximumPropagationDuration, *eventCondition_
     );
 
+    // Expand states based on input state
+    const StateBuilder stateBuilder = {aState};
+
+    Array<State> states = Array<State>::Empty();
+    states.reserve(propagator.accessNumericalSolver().accessObservedStates().getSize());
+
+    for (const State& state : propagator.accessNumericalSolver().accessObservedStates())
+    {
+        states.add(stateBuilder.expand(state.inFrame(aState.accessFrame()), aState));
+    }
+
     return {
         name_,
         dynamics_,
-        propagator.accessNumericalSolver().accessObservedStates(),
+        states,
         conditionSolution.conditionIsSatisfied,
         type_,
     };

src/OpenSpaceToolkit/Astrodynamics/Trajectory/Sequence.cpp

@@ -183,35 +183,47 @@ Sequence::Sequence(
     const Array<Segment>& aSegmentArray,
     const NumericalSolver& aNumericalSolver,
     const Array<Shared<Dynamics>>& aDynamicsArray,
-    const Duration& maximumPropagationDuration,
-    const Size& verbosity
+    const Duration& aMaximumPropagationDuration,
+    const Duration& aMinimumManeuverDuration,
+    const Size& aVerbosityLevel
 )
     : segments_(aSegmentArray),
       numericalSolver_(aNumericalSolver),
       dynamics_(aDynamicsArray),
-      segmentPropagationDurationLimit_(maximumPropagationDuration)
+      segmentPropagationDurationLimit_(aMaximumPropagationDuration),
+      minimumManeuverDuration_(aMinimumManeuverDuration)
 {
-    if (verbosity == 5)
+    if (aMaximumPropagationDuration <= Duration::Zero())
+    {
+        throw ostk::core::error::RuntimeError("Maximum propagation duration must be strictly positive.");
+    }
+
+    if (aMinimumManeuverDuration.isDefined() && aMinimumManeuverDuration <= Duration::Zero())
+    {
+        throw ostk::core::error::RuntimeError("Minimum maneuver duration must be strictly positive.");
+    }
+
+    if (aVerbosityLevel == 5)
     {
         boost::log::core::get()->set_filter(boost::log::trivial::severity >= boost::log::trivial::trace);
     }
-    else if (verbosity == 4)
+    else if (aVerbosityLevel == 4)
     {
         boost::log::core::get()->set_filter(boost::log::trivial::severity >= boost::log::trivial::debug);
     }
-    else if (verbosity == 3)
+    else if (aVerbosityLevel == 3)
     {
         boost::log::core::get()->set_filter(boost::log::trivial::severity >= boost::log::trivial::info);
     }
-    else if (verbosity == 2)
+    else if (aVerbosityLevel == 2)
     {
         boost::log::core::get()->set_filter(boost::log::trivial::severity >= boost::log::trivial::warning);
     }
-    else if (verbosity == 1)
+    else if (aVerbosityLevel == 1)
     {
         boost::log::core::get()->set_filter(boost::log::trivial::severity >= boost::log::trivial::error);
     }
-    else if (verbosity == 0)
+    else if (aVerbosityLevel == 0)
     {
         boost::log::core::get()->set_filter(boost::log::trivial::severity >= boost::log::trivial::fatal);
     }
@@ -248,6 +260,11 @@ Duration Sequence::getMaximumPropagationDuration() const
     return segmentPropagationDurationLimit_;
 }
 
+Duration Sequence::getMinimumManeuverDuration() const
+{
+    return minimumManeuverDuration_;
+}
+
 void Sequence::addSegment(const Segment& aSegment)
 {
     segments_.add(aSegment);
@@ -295,6 +312,18 @@ Sequence::Solution Sequence::solve(const State& aState, const Size& aRepetitionC
 
             BOOST_LOG_TRIVIAL(debug) << "\n" << segmentSolution << std::endl;
 
+            if (segment.getType() == Segment::Type::Maneuver && minimumManeuverDuration_.isDefined())
+            {
+                if (segmentSolution.getPropagationDuration() < minimumManeuverDuration_)
+                {
+                    BOOST_LOG_TRIVIAL(debug)
+                        << "Maneuver duration is less than the minimum maneuver duration. Skipping this maneuver."
+                        << std::endl;
+
+                    continue;
+                }
+            }
+
             segmentSolutions.add(segmentSolution);
 
             // Terminate Sequence unsuccessfully if the segment condition was not satisfied
@@ -346,6 +375,19 @@ Sequence::Solution Sequence::solveToCondition(
 
             BOOST_LOG_TRIVIAL(debug) << "\n" << segmentSolution << std::endl;
 
+            // Skip maneuver if it is less than the minimum maneuver duration
+            if (segment.getType() == Segment::Type::Maneuver && minimumManeuverDuration_.isDefined())
+            {
+                if (segmentSolution.getPropagationDuration() < minimumManeuverDuration_)
+                {
+                    BOOST_LOG_TRIVIAL(debug)
+                        << "Maneuver duration is less than the minimum maneuver duration. Skipping this maneuver."
+                        << std::endl;
+
+                    continue;
+                }
+            }
+
             segmentSolutions.add(segmentSolution);
 
             // Terminate Sequence unsuccessfully if the segment condition was not satisfied
@@ -405,6 +447,10 @@ void Sequence::print(std::ostream& anOutputStream, bool displayDecorator) const
     ostk::core::utils::Print::Line(anOutputStream)
         << "Maximum Propagation Duration:" << segmentPropagationDurationLimit_.toString();
 
+    ostk::core::utils::Print::Line(anOutputStream)
+        << "Minimum Maneuver Duration:"
+        << (minimumManeuverDuration_.isDefined() ? minimumManeuverDuration_.toString() : "Undefined");
+
     if (displayDecorator)
     {
         ostk::core::utils::Print::Footer(anOutputStream);