1
$\begingroup$

The inclination is not that much as It appears in the image, don't forget to check Z scale.I propagated an orbit with semi-major axis of 7000km and the rest of the keplerian orbital elements 0 on GMAT. Since the orbit is circular and in the Equatorial plane, it should not change over time. Initially I have Z = 0 but overtime the value of Z changes significantly ( approx 2 km on day 1). When I plotted this orbit on MATLAB its inclination angle and semi major axis changed over the period of time. What could be the possible explanation for this?

The following is the GMAT Script:

%General Mission Analysis Tool(GMAT) Script
%Created: 2020-08-08 13:43:20


%----------------------------------------
%---------- Spacecraft
%----------------------------------------

Create Spacecraft DefaultSC;
GMAT DefaultSC.DateFormat = UTCGregorian;
GMAT DefaultSC.Epoch = '21 Mar 2019 00:00:00.000';
GMAT DefaultSC.CoordinateSystem = EarthMJ2000Eq;
GMAT DefaultSC.DisplayStateType = Keplerian;
GMAT DefaultSC.SMA = 7000;
GMAT DefaultSC.ECC = 0;
GMAT DefaultSC.INC = 0;
GMAT DefaultSC.RAAN = 0;
GMAT DefaultSC.AOP = 0;
GMAT DefaultSC.TA = 0;
GMAT DefaultSC.DryMass = 850;
GMAT DefaultSC.Cd = 2.2;
GMAT DefaultSC.Cr = 1.8;
GMAT DefaultSC.DragArea = 15;
GMAT DefaultSC.SRPArea = 1;
GMAT DefaultSC.SPADDragScaleFactor = 1;
GMAT DefaultSC.SPADSRPScaleFactor = 1;
GMAT DefaultSC.NAIFId = -10008001;
GMAT DefaultSC.NAIFIdReferenceFrame = -9008001;
GMAT DefaultSC.OrbitColor = Red;
GMAT DefaultSC.TargetColor = Teal;
GMAT DefaultSC.OrbitErrorCovariance = [ 1e+70 0 0 0 0 0 ; 0 1e+70 0 0 0 0 ; 0 0 1e+70 0 0 0 ; 0 0 0 1e+70 0 0 ; 0 0 0 0 1e+70 0 ; 0 0 0 0 0 1e+70 ];
GMAT DefaultSC.CdSigma = 1e+70;
GMAT DefaultSC.CrSigma = 1e+70;
GMAT DefaultSC.Id = 'SatId';
GMAT DefaultSC.Attitude = CoordinateSystemFixed;
GMAT DefaultSC.SPADSRPInterpolationMethod = Bilinear;
GMAT DefaultSC.SPADSRPScaleFactorSigma = 1e+70;
GMAT DefaultSC.SPADDragInterpolationMethod = Bilinear;
GMAT DefaultSC.SPADDragScaleFactorSigma = 1e+70;
GMAT DefaultSC.ModelFile = 'aura.3ds';
GMAT DefaultSC.ModelOffsetX = 0;
GMAT DefaultSC.ModelOffsetY = 0;
GMAT DefaultSC.ModelOffsetZ = 0;
GMAT DefaultSC.ModelRotationX = 0;
GMAT DefaultSC.ModelRotationY = 0;
GMAT DefaultSC.ModelRotationZ = 0;
GMAT DefaultSC.ModelScale = 1;
GMAT DefaultSC.AttitudeDisplayStateType = 'Quaternion';
GMAT DefaultSC.AttitudeRateDisplayStateType = 'AngularVelocity';
GMAT DefaultSC.AttitudeCoordinateSystem = EarthMJ2000Eq;
GMAT DefaultSC.EulerAngleSequence = '321';





%----------------------------------------
%---------- ForceModels
%----------------------------------------

Create ForceModel DefaultProp_ForceModel;
GMAT DefaultProp_ForceModel.CentralBody = Earth;
GMAT DefaultProp_ForceModel.PrimaryBodies = {Earth};
GMAT DefaultProp_ForceModel.Drag = None;
GMAT DefaultProp_ForceModel.SRP = Off;
GMAT DefaultProp_ForceModel.RelativisticCorrection = Off;
GMAT DefaultProp_ForceModel.ErrorControl = RSSStep;
GMAT DefaultProp_ForceModel.GravityField.Earth.Degree = 2;
GMAT DefaultProp_ForceModel.GravityField.Earth.Order = 0;
GMAT DefaultProp_ForceModel.GravityField.Earth.StmLimit = 100;
GMAT DefaultProp_ForceModel.GravityField.Earth.PotentialFile = 'JGM3.cof';
GMAT DefaultProp_ForceModel.GravityField.Earth.TideModel = 'None';

Create ForceModel InternalODEModel;
GMAT InternalODEModel.CentralBody = Earth;
GMAT InternalODEModel.PrimaryBodies = {Earth};
GMAT InternalODEModel.Drag = None;
GMAT InternalODEModel.SRP = Off;
GMAT InternalODEModel.RelativisticCorrection = Off;
GMAT InternalODEModel.ErrorControl = None;
GMAT InternalODEModel.GravityField.Earth.Degree = 2;
GMAT InternalODEModel.GravityField.Earth.Order = 0;
GMAT InternalODEModel.GravityField.Earth.StmLimit = 100;
GMAT InternalODEModel.GravityField.Earth.PotentialFile = 'JGM2.cof';
GMAT InternalODEModel.GravityField.Earth.TideModel = 'None';

%----------------------------------------
%---------- Propagators
%----------------------------------------

Create Propagator DefaultProp;
GMAT DefaultProp.FM = InternalODEModel;
GMAT DefaultProp.Type = AdamsBashforthMoulton;
GMAT DefaultProp.InitialStepSize = 1;
GMAT DefaultProp.Accuracy = 1e-10;
GMAT DefaultProp.MinStep = 1;
GMAT DefaultProp.MaxStep = 1;
GMAT DefaultProp.MaxStepAttempts = 50;
GMAT DefaultProp.StopIfAccuracyIsViolated = false;
GMAT DefaultProp.LowerError = 1e-13;
GMAT DefaultProp.TargetError = 9.999999999999999e-12;

%----------------------------------------
%---------- Subscribers
%----------------------------------------

Create ReportFile ReportFile1;
GMAT ReportFile1.SolverIterations = Current;
GMAT ReportFile1.UpperLeft = [ 0 0 ];
GMAT ReportFile1.Size = [ 0 0 ];
GMAT ReportFile1.RelativeZOrder = 0;
GMAT ReportFile1.Maximized = false;
GMAT ReportFile1.Filename = <Enter file name and location>;
GMAT ReportFile1.Precision = 20;
GMAT ReportFile1.Add = {DefaultSC.UTCGregorian, DefaultSC.EarthMJ2000Eq.X, DefaultSC.EarthMJ2000Eq.Y, DefaultSC.EarthMJ2000Eq.Z, DefaultSC.EarthMJ2000Eq.VX, DefaultSC.EarthMJ2000Eq.VY, DefaultSC.EarthMJ2000Eq.VZ, DefaultSC.DefaultProp_ForceModel.AccelerationX, DefaultSC.DefaultProp_ForceModel.AccelerationY, DefaultSC.DefaultProp_ForceModel.AccelerationZ};
GMAT ReportFile1.WriteHeaders = true;
GMAT ReportFile1.LeftJustify = On;
GMAT ReportFile1.ZeroFill = Off;
GMAT ReportFile1.FixedWidth = false;
GMAT ReportFile1.Delimiter = ',';
GMAT ReportFile1.ColumnWidth = 30;
GMAT ReportFile1.WriteReport = true;


%----------------------------------------
%---------- Mission Sequence
%----------------------------------------

BeginMissionSequence;
Propagate DefaultProp(DefaultSC) {DefaultSC.ElapsedSecs = 86400};
$\endgroup$
  • $\begingroup$ GMAT DefaultProp.Type = AdamsBashforthMoulton There's your problem. GMAT has lousy propagators. $\endgroup$ – David Hammen Aug 13 at 23:16
  • 2
    $\begingroup$ @DavidHammen An Adams-Bashforth algorithm shouldn't be that bad! I suspect there's a bug in that code somewhere. Which, I guess, makes it a lousy propagator! ;-) $\endgroup$ – Tom Spilker Aug 13 at 23:48
  • $\begingroup$ If no perturbations act along Z there should be no accumulation of numerical error in that direction. I would look at some unexpected perturbation, and in your script it looks like J2 is enabled (GMAT InternalODEModel.GravityField.Earth.Degree = 2, although I'm not confident with GMAT). Have you tried looking at that? $\endgroup$ – LeWavite Aug 14 at 17:23
  • $\begingroup$ It's because you use 'JGM3.cof' or 'JGM2.cof' (I don't know the exact meaning of the GMAT parameters). If you use a spherical Earth, you should see no more than 40 m. $\endgroup$ – Cristiano Aug 15 at 8:28
2
$\begingroup$

In the InternalODEModel of your script, you have J2 enabled. Set the degree to zero and the problem should be fixed.

Create ForceModel InternalODEModel;
GMAT InternalODEModel.CentralBody = Earth;
GMAT InternalODEModel.PrimaryBodies = {Earth};
GMAT InternalODEModel.Drag = None;
GMAT InternalODEModel.SRP = Off;
GMAT InternalODEModel.RelativisticCorrection = Off;
GMAT InternalODEModel.ErrorControl = None;
GMAT InternalODEModel.GravityField.Earth.Degree = 2;  <--------- PROBLEM
GMAT InternalODEModel.GravityField.Earth.Order = 0;
GMAT InternalODEModel.GravityField.Earth.StmLimit = 100;
GMAT InternalODEModel.GravityField.Earth.PotentialFile = 'JGM2.cof';
GMAT InternalODEModel.GravityField.Earth.TideModel = 'None';
| improve this answer | |
$\endgroup$
  • $\begingroup$ I don't understand! The OP says " Since the orbit is circular and in the Equatorial plane..." how could $J_2$ have any effect at all? If it were $J_3$ or tesserals it would be a different story. I can never remember the difference between "Degree" and "Order" but does Degree = 2 really mean $J_2$? If so, I don't see how this is the correct explanation. $\endgroup$ – uhoh Aug 15 at 3:46
  • $\begingroup$ To double check, do GMAT DefaultSC.CoordinateSystem = EarthMJ2000Eq; and GMAT DefaultSC.INC = 0; confirm that this is equatorial? Also, could effects of Sun or Moon be present, or does the OP's script ensure these are turned off? Thanks! $\endgroup$ – uhoh Aug 15 at 3:50
  • 2
    $\begingroup$ @uhoh: yes, J2 is the second-order and zeroth-degree harmonic. Judging from the GMAT script, the OP's orbit is equatorial in the mean equator and equinox frame of J2000. Since the Earth's equator moves due to precession and nutation, there can be an out-of-(equatorial-)plane component of the J2 perturbation in the EMEJ2000 frame. I don't know if this can justify the 2km/day error mentioned by the OP, but an easy way to check would be to set the initial epoch at J2000 and see whether the final Z component decreases. Or just turn off J2 :) $\endgroup$ – LeWavite Aug 15 at 23:20
  • $\begingroup$ @LeWavite Excellent detective work! However if I understand the plot in the question correctly, the ~2 km in ~1 day change demonstrates a drift in inclination and not precession, while the J2000.0 frame issue you mention would result in an inclination of 0.3 degrees in 2020 that would result a precession, and as you suggest that may not be even visible in one day on this plot. But 3D is tricky, maybe what's shown in the image is precession around a 0.3 degree tilted axis and I'm just not seeing it? $\endgroup$ – uhoh Aug 16 at 1:58

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.