Overhaul of ReturnFromMoon initial guesser

It now uses Brent's minimization algorithm to try to guess the inclination of the return ellipse and if it is type I or type II return orbit. Restored the 1/2 time midpoint constraints since that helps convergence speed in the hyperbolic initial conditions case and seems to improve other ones as well. This also now finds some better solutions as evidenced by the one test which dropped from 1520 m/s to 1092 m/s. Signed-off-by: Lamont Granquist <lamont@scriptkiddie.org>
MuMech · Nov 22, 2023 · a61ca04 · a61ca04
1 parent 423d774
commit a61ca04
Show file tree

Hide file tree

Showing 5 changed files with 320 additions and 176 deletions.
diff --git a/MechJeb2/OrbitalManeuverCalculator.cs b/MechJeb2/OrbitalManeuverCalculator.cs
@@ -354,6 +354,8 @@ public static Vector3d DeltaVAndTimeForCheapestCourseCorrection(Orbit o, double
 
         public static (Vector3d dv, double dt) DeltaVAndTimeForMoonReturnEjection(Orbit o, double ut, double targetPrimaryRadius)
         {
+            var solver = new ReturnFromMoon();
+
             CelestialBody moon = o.referenceBody;
             CelestialBody primary = moon.referenceBody;
             (V3 moonR0, V3 moonV0) = moon.orbit.RightHandedStateVectorsAtUT(ut);
@@ -362,7 +364,7 @@ public static (Vector3d dv, double dt) DeltaVAndTimeForMoonReturnEjection(Orbit
 
             double dtmin = o.eccentricity >= 1 ? 0 : double.NegativeInfinity;
 
-            (V3 dv, double dt, double newPeR) = ReturnFromMoon.NextManeuver(primary.gravParameter, moon.gravParameter, moonR0,
+            (V3 dv, double dt, double newPeR) = solver.NextManeuver(primary.gravParameter, moon.gravParameter, moonR0,
                 moonV0, moonSOI, r0, v0, targetPrimaryRadius, 0, dtmin);
 
             Debug.Log($"Solved PeR from calcluator: {newPeR}");

diff --git a/MechJebLib/Functions/Astro.cs b/MechJebLib/Functions/Astro.cs
@@ -363,6 +363,13 @@ public static V3 VelocityForInclination(V3 r, V3 v, double newInc)
             return vf;
         }
 
+        [MethodImpl(MethodImplOptions.AggressiveInlining)]
+        public static V3 VelocityForInclination(V3 r, double vmag, double newInc)
+        {
+            V3 vf = ENUHeadingForInclination(newInc, r) * vmag;
+            return ENUToECI(r, vf);
+        }
+
         [MethodImpl(MethodImplOptions.AggressiveInlining)]
         public static V3 VelocityForFPA(V3 r, V3 v, double newFPA)
         {
@@ -564,6 +571,14 @@ public static (V3 p, V3 q) PerifocalFromElements(double mu, double l, double ecc
             return (one * l / (1 + ecc * cnu), two * Sqrt(mu / l));
         }
 
+        [MethodImpl(MethodImplOptions.AggressiveInlining)]
+        public static (V3 r, V3 v) StateVectorsAtTrueAnomaly(double mu, V3 r, V3 v, double nu)
+        {
+            // TODO? there may be a more efficient way to do this
+            (double _, double ecc, double inc, double lan, double argp, double _, double l) = KeplerianFromStateVectors(mu, r, v);
+            return StateVectorsFromKeplerian(mu, l, ecc, inc, lan, argp, nu);
+        }
+
         [MethodImpl(MethodImplOptions.AggressiveInlining)]
         public static (V3 r, V3 v) StateVectorsFromKeplerian(double mu, double l, double ecc, double inc, double lan, double argp, double nu)
         {