Skip to content

Commit

Permalink
More ReturnFromMoon cleanup
Browse files Browse the repository at this point in the history 
Bunch of reorganization of the files and now the test files match
up better.

Much better construction of a nearly-feasible initial guess that
helps to minimize the errors at the initial point.

Signed-off-by: Lamont Granquist <lamont@scriptkiddie.org>
  • Loading branch information
@lamont-granquist
lamont-granquist committed Aug 25, 2023
1 parent 7624a87 commit d78272a
Show file tree
Hide file tree
Showing 12 changed files with 460 additions and 296 deletions.
1 change: 1 addition & 0 deletions MechJeb2/MechJeb2.csproj
View file
Expand Up @@ -118,6 +118,7 @@
<Compile Include="MechJebLib\Core\TwoBody\Farnocchia.cs"/>
<Compile Include="MechJebLib\Maneuvers\ChangeOrbitalElement.cs"/>
<Compile Include="MechJebLib\Maneuvers\ReturnFromMoon.cs"/>
<Compile Include="MechJebLib\Maneuvers\Simple.cs"/>
<Compile Include="MechJebLib\Primitives\Dual.cs"/>
<Compile Include="MechJebLib\Primitives\DualV3.cs"/>
<Compile Include="MechJebLib\Primitives\Vn.cs"/>
Expand Down
81 changes: 0 additions & 81 deletions MechJeb2/MechJebLib/Core/Functions/Maneuvers.cs
View file
Expand Up @@ -5,90 +5,9 @@

#nullable enable

using MechJebLib.Core.TwoBody;
using MechJebLib.Primitives;
using MechJebLib.Utils;

namespace MechJebLib.Core.Functions
{
public class Maneuvers
{
public static V3 DeltaVRelativeToCircularVelocity(double mu, V3 r, V3 v, double n = 1.0)
{
Check.Finite(mu);
Check.Finite(r);
Check.Finite(v);
Check.Positive(mu);
Check.NonZero(r);

var h = V3.Cross(r, v);
return n * Maths.CircularVelocityFromHvec(mu, r, h) - v;
}

public static V3 DeltaVToCircularize(double mu, V3 r, V3 v)
{
Check.Finite(mu);
Check.Finite(r);
Check.Finite(v);
Check.Positive(mu);
Check.NonZero(r);

var h = V3.Cross(r, v);
return Maths.CircularVelocityFromHvec(mu, r, h) - v;
}

public static V3 DeltaVToCircularizeAfterTime(double mu, V3 r, V3 v, double dt)
{
Check.Finite(mu);
Check.Finite(r);
Check.Finite(v);
Check.Finite(dt);
Check.Positive(mu);
Check.NonZero(r);

(V3 r1, V3 v1) = Farnocchia.Solve(mu, dt, r, v);
var h = V3.Cross(r1, v1);
return Maths.CircularVelocityFromHvec(mu, r, h) - v1;
}

public static V3 DeltaVToCircularizeAtPeriapsis(double mu, V3 r, V3 v)
{
Check.Finite(mu);
Check.Finite(r);
Check.Finite(v);
Check.Positive(mu);
Check.NonZero(r);

double dt = Maths.TimeToNextPeriapsis(mu, r, v);

Check.Finite(dt);

return DeltaVToCircularizeAfterTime(mu, r, v, dt);
}

public static V3 DeltaVToCircularizeAtApoapsis(double mu, V3 r, V3 v)
{
Check.Finite(mu);
Check.Finite(r);
Check.Finite(v);
Check.Positive(mu);
Check.NonZero(r);

double dt = Maths.TimeToNextApoapsis(mu, r, v);

Check.Finite(dt);

return DeltaVToCircularizeAfterTime(mu, r, v, dt);
}

public static V3 DeltaVToChangeInclination(V3 r, V3 v, double newInc)
{
return Maths.VelocityForInclination(r, v, newInc) - v;
}

public static V3 DeltaVToChangeFPA(V3 r, V3 v, double newFPA)
{
return Maths.VelocityForFPA(r, v, newFPA) - v;
}
}
}
14 changes: 10 additions & 4 deletions MechJeb2/MechJebLib/Core/Maths.cs
View file
Expand Up @@ -134,7 +134,8 @@ public static (Dual sma, Dual ecc) SmaEccFromStateVectors(Dual mu, DualV3 r, Dua
{
var h = DualV3.Cross(r, v);
Dual sma = SmaFromStateVectors(mu, r, v);
return (sma, Dual.Sqrt(Dual.Max(1 - h.sqrMagnitude / (sma * mu), 0)));
var ecc = Dual.Sqrt(Dual.Max(1 - h.sqrMagnitude / (sma * mu), 0));
return (sma, ecc);
}

public static double SmaFromStateVectors(double mu, V3 r, V3 v) => mu / (2.0 * mu / r.magnitude - V3.Dot(v, v));
Expand Down Expand Up @@ -401,10 +402,8 @@ public static double TimeToNextPeriapsis(double mu, V3 r, V3 v)
return Clamp2Pi(-manom) / meanMotion;
}

public static double TimeToNextTrueAnomaly(double mu, V3 r, V3 v, double nu2)
public static double TimetoNextTrueAnomaly(double mu, double sma, double ecc, double nu1, double nu2)
{
(double sma, double ecc, _, _, _, double nu1, _) = KeplerianFromStateVectors(mu, r, v);

double meanMotion = MeanMotion(mu, sma);

double manom1 = Angles.MFromNu(nu1, ecc);
Expand All @@ -416,6 +415,13 @@ public static double TimeToNextTrueAnomaly(double mu, V3 r, V3 v, double nu2)
return (manom2 - manom1) / meanMotion;
}

public static double TimeToNextTrueAnomaly(double mu, V3 r, V3 v, double nu2)
{
(double sma, double ecc, _, _, _, double nu1, _) = KeplerianFromStateVectors(mu, r, v);

return TimetoNextTrueAnomaly(mu, sma, ecc, nu1, nu2);
}

public static double TimeToNextRadius(double mu, V3 r, V3 v, double radius)
{
double nu1 = TrueAnomalyFromRadius(mu, r, v, radius);
Expand Down
10 changes: 5 additions & 5 deletions MechJeb2/MechJebLib/Maneuvers/ChangeOrbitalElement.cs
View file
Expand Up @@ -106,7 +106,7 @@ private static void NLPFunction(double[] x, double[] fi, object obj)
}
}

private static V3 DeltaV(double mu, V3 r, V3 v, double value, Type type, bool optguard=false)
private static V3 DeltaV(double mu, V3 r, V3 v, double value, Type type, bool optguard = false)
{
if (!mu.IsFinite())
throw new ArgumentException("bad mu in ChangeOrbitalElement");
Expand All @@ -116,7 +116,7 @@ private static V3 DeltaV(double mu, V3 r, V3 v, double value, Type type, bool op
throw new ArgumentException("bad v in ChangeOrbitalElement");

const double DIFFSTEP = 1e-7;
const double EPSX = 1e-15;
const double EPSX = 1e-7;
const int MAXITS = 1000;

const int NVARIABLES = 2;
Expand All @@ -137,7 +137,7 @@ private static V3 DeltaV(double mu, V3 r, V3 v, double value, Type type, bool op
// changing the ECC is actually a global optimization problem due to basins around parabolic ecc == 1.0
double boost = Sign(value - 1.0) * 0.1 + Sqrt(2);

V3 dv = Core.Functions.Maneuvers.DeltaVRelativeToCircularVelocity(1.0, p, q, boost);
V3 dv = Simple.DeltaVRelativeToCircularVelocity(1.0, p, q, boost);
x[0] = dv.x;
x[1] = dv.y;
}
Expand All @@ -161,7 +161,7 @@ private static V3 DeltaV(double mu, V3 r, V3 v, double value, Type type, bool op
}

alglib.minnlcoptimize(state, NLPFunction2, null, args);
alglib.minnlcresults(state, out x, out alglib.minnlcreport rep);
alglib.minnlcresults(state, out double[] x2, out alglib.minnlcreport rep);

if (rep.terminationtype < 0)
throw new Exception(
Expand All @@ -175,7 +175,7 @@ private static V3 DeltaV(double mu, V3 r, V3 v, double value, Type type, bool op
throw new Exception("alglib optguard caught an error, i should report better on errors now");
}

return rot * new V3(x[0], x[1], 0) * scale.VelocityScale;
return rot * new V3(x2[0], x2[1], 0) * scale.VelocityScale;
}

public static V3 ChangePeriapsis(double mu, V3 r, V3 v, double peR, bool optguard = false)
Expand Down

0 comments on commit d78272a

Please sign in to comment.