Skip to content

Commit

Permalink
Fix Hohmann bugs and implement more features
Browse files Browse the repository at this point in the history 
Previous PR was just broke, fixes a few things, implements
planning the capture burn, doing rendezvous/transfer, the
"period offset" is now just a "lag time" (chasing time?)

Signed-off-by: Lamont Granquist <lamont@scriptkiddie.org>
  • Loading branch information
@lamont-granquist
lamont-granquist committed Sep 28, 2023
1 parent a926641 commit 78f7d3c
Show file tree
Hide file tree
Showing 7 changed files with 82 additions and 100 deletions.
10 changes: 5 additions & 5 deletions Localization/en-us.cfg
View file
Expand Up @@ -115,12 +115,12 @@ Localization

//Hohmann transfer
#MechJeb_Hohm_title = bi-impulsive (Hohmann) transfer to target
#MechJeb_Hohm_intercept_only = intercept only, no capture burn (impact/flyby)
#MechJeb_Hohm_simpleTransfer = simple coplanar Hohmann transfer
#MechJeb_Hohm_Label1 = fractional target period offset
#MechJeb_Hohm_intercept_only = no insertion burn (impact/flyby)
#MechJeb_Hohm_simpleTransfer = coplanar maneuver
#MechJeb_Hohm_Label1 = lag time offset

#MechJeb_Hohm_Exception1 = must select a target for the bi-impulsive transfer.
#MechJeb_Hohm_Exception2 = target for bi-impulsive transfer must be in the same sphere of influence.
#MechJeb_Hohm_Exception1 = must select a target for the maneuver.
#MechJeb_Hohm_Exception2 = target for maneuver must be in the same sphere of influence.
#MechJeb_Hohm_Exception3 = ascending node with target doesn't exist.
#MechJeb_Hohm_Exception4 = descending node with target doesn't exist.
#MechJeb_Hohm_Exception5 = neither ascending nor descending node with target exists.
Expand Down
29 changes: 11 additions & 18 deletions MechJeb2/Maneuver/OperationTransfer.cs
View file
Expand Up @@ -25,7 +25,7 @@ public class OperationGeneric : Operation

[UsedImplicitly]
[Persistent(pass = (int)Pass.GLOBAL)]
public EditableDouble PeriodOffset = 0;
public EditableDouble LagTime = 0;

[Persistent(pass = (int)Pass.GLOBAL)]
public EditableTime MinDepartureUT = 0;
Expand All @@ -51,15 +51,17 @@ public override void DoParametersGUI(Orbit o, double universalTime, MechJebModul
Capture =
!GUILayout.Toggle(!Capture, Localizer.Format("#MechJeb_Hohm_intercept_only")); //no capture burn (impact/flyby)
if (Capture)
PlanCapture = GUILayout.Toggle(PlanCapture, "Plan capture burn");
PlanCapture = GUILayout.Toggle(PlanCapture, "Plan insertion burn");
Coplanar = GUILayout.Toggle(Coplanar, Localizer.Format("#MechJeb_Hohm_simpleTransfer")); //coplanar maneuver
GuiUtils.SimpleTextBox(Localizer.Format("#MechJeb_Hohm_Label1"), PeriodOffset); //fractional target period offset
GUILayout.BeginHorizontal();
if (GUILayout.Toggle(Rendezvous, "Rendezvous"))
Rendezvous = true;
if (GUILayout.Toggle(!Rendezvous, "Transfer"))
Rendezvous = false;
GUILayout.EndHorizontal();
if (Rendezvous)
GuiUtils.SimpleTextBox(Localizer.Format("#MechJeb_Hohm_Label1"), LagTime, "sec"); //fractional target period offset

/*
if (!Coplanar)
{
Expand All @@ -70,31 +72,20 @@ public override void DoParametersGUI(Orbit o, double universalTime, MechJebModul

protected override List<ManeuverParameters> MakeNodesImpl(Orbit o, double universalTime, MechJebModuleTargetController target)
{
double ut;

if (!target.NormalTargetExists)
{
throw new OperationException(Localizer.Format("#MechJeb_Hohm_Exception1")); //must select a target for the bi-impulsive transfer.
}

if (o.referenceBody != target.TargetOrbit.referenceBody)
{
throw
new OperationException(
Localizer.Format("#MechJeb_Hohm_Exception2")); //target for bi-impulsive transfer must be in the same sphere of influence.
}

Vector3d dV;

Orbit targetOrbit = target.TargetOrbit;

if (PeriodOffset != 0)
{
targetOrbit = target.TargetOrbit.Clone();
targetOrbit.MutatedOrbit(PeriodOffset);
}
double lagTime = Rendezvous ? LagTime.val : 0;

(dV, ut) = OrbitalManeuverCalculator.DeltaVAndTimeForHohmannTransfer(o, targetOrbit, universalTime, coplanar: Coplanar, rendezvous: Rendezvous);
(Vector3d dV1, double ut1, Vector3d dV2, double ut2) =
OrbitalManeuverCalculator.DeltaVAndTimeForHohmannTransfer(o, targetOrbit, universalTime, lagTime, Coplanar, Rendezvous, Capture);

/*
else
Expand Down Expand Up @@ -132,7 +123,9 @@ protected override List<ManeuverParameters> MakeNodesImpl(Orbit o, double univer
}
*/

return new List<ManeuverParameters> { new ManeuverParameters(dV, ut) };
if (Capture && PlanCapture)
return new List<ManeuverParameters> { new ManeuverParameters(dV1, ut1), new ManeuverParameters(dV2, ut2) };
return new List<ManeuverParameters> { new ManeuverParameters(dV1, ut1) };
}
}
}
4 changes: 0 additions & 4 deletions MechJeb2/MechJeb2.csproj
View file
Expand Up @@ -56,9 +56,6 @@
<Compile Include="FlyingSim\SimulatedParachute.cs"/>
<Compile Include="FlyingSim\SimulatedPart.cs"/>
<Compile Include="FlyingSim\SimulatedVessel.cs"/>
<Compile Include="FuelFlowSimulation.cs"/>
<Compile Include="FuelNode.cs"/>
<Compile Include="FuelStats.cs"/>
<Compile Include="GlobalSuppressions.cs"/>
<Compile Include="GLUtils.cs"/>
<Compile Include="GuiUtils.cs"/>
Expand Down Expand Up @@ -228,7 +225,6 @@
<Compile Include="MechJebModuleSpaceplaneAutopilot.cs"/>
<Compile Include="MechJebModuleSpaceplaneGuidance.cs"/>
<Compile Include="MechJebModuleSpinupController.cs"/>
<Compile Include="MechJebModuleStageStatsOld.cs"/>
<Compile Include="MechJebModuleStagingController.cs"/>
<Compile Include="MechJebModuleSuicideTimer.cs"/>
<Compile Include="MechJebModuleTargetController.cs"/>
Expand Down
83 changes: 53 additions & 30 deletions MechJeb2/MechJebLib/Maneuvers/CoplanarTransfer.cs
View file
Expand Up @@ -11,11 +11,11 @@ public class CoplanarTransfer
{
private struct Args
{
public V3 R1;
public V3 V1;
public V3 R2;
public V3 V2;
public bool Rendezvous;
public V3 R1;
public V3 V1;
public V3 R2;
public V3 V2;
public bool Capture;
}

private static (V3 dv1, V3 dv2) LambertFunction(double dt, double tt, double dt2, Args args)
Expand All @@ -26,10 +26,7 @@ private static (V3 dv1, V3 dv2) LambertFunction(double dt, double tt, double dt2
V3 v2 = args.V2;

(V3 rburn, V3 vburn) = Shepperd.Solve(1.0, dt, r1, v1);
double targdt = dt + tt;
if (!args.Rendezvous)
targdt += dt2;
(V3 rf2, V3 vf2) = Shepperd.Solve(1.0, targdt, r2, v2);
(V3 rf2, V3 vf2) = Shepperd.Solve(1.0, dt + tt + dt2, r2, v2);
(V3 vi, V3 vf) = Gooding.Solve(1.0, rburn, vburn, rf2, tt, 0);

return (vi - vburn, vf2 - vf);
Expand All @@ -51,16 +48,16 @@ private static void NLPFunction(double[] x, ref double func, object obj)

(V3 dv1, V3 dv2) = LambertFunction(dt, tt, dt2, args);

func = dv1.magnitude + dv2.magnitude;

if (args.Rendezvous) // drive unused dt2 to zero
func += dt2 * dt2;
if (args.Capture)
func = dv1.magnitude + dv2.magnitude;
else
func = dv1.magnitude;
}

public static (V3 dv1, double dt, V3 dv2, double tt) Maneuver(double mu, V3 r1, V3 v1, V3 r2, V3 v2, double dtguess, double dt2guess,
bool coplanar = true, bool rendezvous = true, double dtmin = double.NegativeInfinity, double dtmax = double.PositiveInfinity,
double ttmin = double.NegativeInfinity,
double ttmax = double.PositiveInfinity, bool optguard = false)
public static (V3 dv1, double dt, V3 dv2, double tt) Maneuver(double mu, V3 r1, V3 v1, V3 r2, V3 v2, double dtguess, double dt2Guess,
bool coplanar = true, bool rendezvous = true, bool capture = true, double dtmin = double.NegativeInfinity,
double dtmax = double.PositiveInfinity,
double ttmin = double.NegativeInfinity, double ttmax = double.PositiveInfinity, double dt2Min = double.NegativeInfinity, double dt2Max = double.PositiveInfinity, bool optguard = false)
{
if (!mu.IsFinite())
throw new ArgumentException("bad mu in ChangeOrbitalElement");
Expand Down Expand Up @@ -99,29 +96,31 @@ private static void NLPFunction(double[] x, ref double func, object obj)
r2 /= scale.LengthScale;
v2 /= scale.VelocityScale;
dtguess /= scale.TimeScale;
dt2guess /= scale.TimeScale;
dt2Guess /= scale.TimeScale;
dtmin /= scale.TimeScale;
dtmax /= scale.TimeScale;
ttmin /= scale.TimeScale;
ttmax /= scale.TimeScale;
dt2Min /= scale.TimeScale;
dt2Max /= scale.TimeScale;

(_, _, double ttguess, _) = Maths.HohmannTransferParameters(1.0, r1, r2);

x[0] = dtguess;
x[1] = ttguess;
x[2] = dt2guess;
x[2] = dt2Guess;

var args = new Args
{
R1 = r1,
V1 = v1,
R2 = r2,
V2 = v2,
Rendezvous = rendezvous
Capture = capture
};

double[] bndl = { dtmin, ttmin, double.NegativeInfinity };
double[] bndu = { dtmax, ttmax, double.PositiveInfinity };
double[] bndl = { dtmin, ttmin, dt2Min };
double[] bndu = { dtmax, ttmax, dt2Max};

alglib.minbleiccreatef(x, DIFFSTEP, out alglib.minbleicstate state);
alglib.minbleicsetbc(state, bndl, bndu);
Expand Down Expand Up @@ -158,24 +157,48 @@ private static void NLPFunction(double[] x, ref double func, object obj)
}

public static (V3 dv1, double dt, V3 dv2, double tt) NextManeuver(double mu, V3 r1, V3 v1, V3 r2, V3 v2, int maxiter = 50,
bool coplanar = true, bool rendezvous = true, bool optguard = false)
double lagTime = double.NaN, bool coplanar = true, bool rendezvous = true, bool capture = true, bool optguard = false)
{
double synodicPeriod = Maths.SynodicPeriod(mu, r1, v1, r2, v2);
double targetPeriod = Maths.PeriodFromStateVectors(mu, r2, v2);
double dtguess = 0;

for (int iter = 0; iter < maxiter; iter++)
{
(V3 dv1, double dt, V3 dv2, double tt) =
Maneuver(mu, r1, v1, r2, v2, dtguess, 0, coplanar, optguard: optguard, rendezvous: rendezvous);

// we have to try the other side of the target orbit since we might get eg. the DN instead of the AN when the AN is closer
// (this may be insufficient and may need more of a search box but then we're O(N^2) and i think basinhopping or porkchop
// plots will be the better solution)
if (!rendezvous)
V3 dv1, dv2;
double dt, tt;
double dt2Guess = 0;

if (rendezvous)
{
double dt2Min = 0;
double dt2Max = 0;
if (lagTime.IsFinite())
{
dt2Min = -lagTime;
dt2Max = -lagTime;
dt2Guess = -lagTime;
}

(dv1, dt, dv2, tt) =
Maneuver(mu, r1, v1, r2, v2, dtguess, dt2Guess, dt2Min: dt2Min, dt2Max: dt2Max, coplanar: coplanar, rendezvous: true, capture: capture,
optguard: optguard);
}

else
{
(dv1, dt, dv2, tt) =
Maneuver(mu, r1, v1, r2, v2, dtguess, dt2Guess, coplanar, rendezvous: false, capture,
optguard: optguard);

// we have to try the other side of the target orbit since we might get eg. the DN instead of the AN when the AN is closer
// (this may be insufficient and may need more of a search box but then we're O(N^2) and i think basinhopping or porkchop
// plots will be the better solution)

(V3 a, double b, V3 c, double d) =
Maneuver(mu, r1, v1, r2, v2, dtguess, targetPeriod * 0.5, coplanar, optguard: optguard, rendezvous: rendezvous);
Maneuver(mu, r1, v1, r2, v2, dtguess, targetPeriod * 0.5, coplanar: coplanar, optguard: optguard,
rendezvous: false, capture: capture);
if (b > 0 && (b < dt || dt < 0))
{
dv1 = a;
Expand Down
4 changes: 2 additions & 2 deletions MechJeb2/MechJebModuleRendezvousAutopilot.cs
View file
Expand Up @@ -135,8 +135,8 @@ public override void Drive(FlightCtrlState s)
{
//We're not on an intercept course, but we have a circular orbit in the right plane.

(Vector3d hohmannDV, double hohmannUT) =
OrbitalManeuverCalculator.DeltaVAndTimeForHohmannTransfer(Orbit, Core.Target.TargetOrbit, VesselState.time, false);
(Vector3d hohmannDV, double hohmannUT, _, _) =
OrbitalManeuverCalculator.DeltaVAndTimeForHohmannTransfer(Orbit, Core.Target.TargetOrbit, VesselState.time, coplanar: false);

double numPhasingOrbits = (hohmannUT - VesselState.time) / Orbit.period;

Expand Down
4 changes: 2 additions & 2 deletions MechJeb2/MechJebModuleRendezvousGuidance.cs
View file
Expand Up @@ -109,8 +109,8 @@ protected override void WindowGUI(int windowID)

if (GUILayout.Button(Localizer.Format("#MechJeb_RZplan_button3"))) //"Intercept with Hohmann transfer"
{
(Vector3d dV, double UT) =
OrbitalManeuverCalculator.DeltaVAndTimeForHohmannTransfer(Orbit, Core.Target.TargetOrbit, VesselState.time, false);
(Vector3d dV, double UT, _, _) =
OrbitalManeuverCalculator.DeltaVAndTimeForHohmannTransfer(Orbit, Core.Target.TargetOrbit, VesselState.time, coplanar: false);
Vessel.RemoveAllManeuverNodes();
Vessel.PlaceManeuverNode(Orbit, dV, UT);
}
Expand Down
48 changes: 9 additions & 39 deletions MechJeb2/OrbitalManeuverCalculator.cs
View file
@@ -1,12 +1,10 @@
using System;
using MechJebLib.Core;
using MechJebLib.Core.TwoBody;
using MechJebLib.Maneuvers;
using MechJebLib.Primitives;
using Smooth.Pools;
using UnityEngine;
using static MechJebLib.Statics;
using Debug = UnityEngine.Debug;

namespace MuMech
{
Expand Down Expand Up @@ -155,15 +153,19 @@ public static Vector3d DeltaVAndTimeToMatchPlanesDescending(Orbit o, Orbit targe
//The output burnUT will be the first transfer window found after the given UT.
//Assumes o and target are in approximately the same plane, and orbiting in the same direction.
//Also assumes that o is a perfectly circular orbit (though result should be OK for small eccentricity).
public static ( Vector3d dV, double burnUT) DeltaVAndTimeForHohmannTransfer(Orbit o, Orbit target, double ut, bool coplanar = true, bool rendezvous = true)
public static ( Vector3d dV1, double UT1, Vector3d dV2, double UT2) DeltaVAndTimeForHohmannTransfer(Orbit o, Orbit target, double ut,
double lagTime = double.NaN,
bool coplanar = true, bool rendezvous = true, bool capture = true)
{
(V3 r1, V3 v1) = o.RightHandedStateVectorsAtUT(ut);
(V3 r2, V3 v2) = target.RightHandedStateVectorsAtUT(ut);

(V3 dv, double dt, _, _) =
CoplanarTransfer.NextManeuver(o.referenceBody.gravParameter, r1, v1, r2, v2, coplanar: coplanar, rendezvous: rendezvous);
(V3 dv1, double dt1, V3 dv2, double dt2) =
CoplanarTransfer.NextManeuver(o.referenceBody.gravParameter, r1, v1, r2, v2, lagTime: lagTime, coplanar: coplanar,
rendezvous: rendezvous,
capture: capture);

return (dv.V3ToWorld(), dt);
return (dv1.V3ToWorld(), ut + dt1, dv2.V3ToWorld(), ut + dt2);
}

// Computes the delta-V of a burn at a given time that will put an object with a given orbit on a
Expand Down Expand Up @@ -191,38 +193,6 @@ public static ( Vector3d dV, double burnUT) DeltaVAndTimeForHohmannTransfer(Orbi
return ((transferVi - vi).V3ToWorld(), (vf - transferVf).V3ToWorld());
}

// Lambert Solver Driver function.
//
// This uses Shepperd's method instead of using KSP's Orbit class.
//
// The reference time is usually 'now' or the first time the burn can start.
//
// GM - grav parameter of the celestial
// pos - position of the source orbit at a reference time
// vel - velocity of the source orbit at a reference time
// tpos - position of the target orbit at a reference time
// tvel - velocity of the target orbit at a reference time
// DT - time of the burn in seconds after the reference time
// TT - transfer time of the burn in seconds after the burn time
// secondDV - the second burn dV
// returns - the first burn dV
//
private static Vector3d DeltaVToInterceptAtTime(double GM, Vector3d pos, Vector3d vel, Vector3d tpos, Vector3d tvel, double dt, double tt,
out Vector3d secondDV, bool posigrade = true)
{
// advance the source orbit to ref + DT
(V3 pos1, V3 vel1) = Shepperd.Solve(GM, dt, pos.ToV3(), vel.ToV3());

// advance the target orbit to ref + DT + TT
(V3 pos2, V3 vel2) = Shepperd.Solve(GM, dt + tt, tpos.ToV3(), tvel.ToV3());

(V3 transferVi, V3 transferVf) = Gooding.Solve(GM, pos1, vel1, pos2, posigrade ? tt : -tt, 0);

secondDV = (vel2 - transferVf).ToVector3d();

return (transferVi - vel1).ToVector3d();
}

// This does a line-search to find the burnUT for the cheapest course correction that will intercept exactly
public static Vector3d DeltaVAndTimeForCheapestCourseCorrection(Orbit o, double UT, Orbit target, out double burnUT)
{
Expand Down Expand Up @@ -316,7 +286,7 @@ public static Vector3d DeltaVAndTimeForCheapestCourseCorrection(Orbit o, double
if (syncPhaseAngle)
{
//time the ejection burn to intercept the target
(idealDeltaV, idealBurnUT) = DeltaVAndTimeForHohmannTransfer(planetOrbit, target, UT);
(idealDeltaV, idealBurnUT, _, _) = DeltaVAndTimeForHohmannTransfer(planetOrbit, target, UT);
}
else
{
Expand Down

0 comments on commit 78f7d3c

Please sign in to comment.