pub struct IntegratorOptionsBuilder<TypedBuilderFields = ((), (), (), (), (), (), (), ())> { /* private fields */ }
Expand description
Builder for IntegratorOptions
instances.
See IntegratorOptions::builder()
for more info.
Implementations§
Source§impl<__min_step, __max_step, __tolerance, __attempts, __fixed_step, __error_ctrl, __integration_frame> IntegratorOptionsBuilder<((), __min_step, __max_step, __tolerance, __attempts, __fixed_step, __error_ctrl, __integration_frame)>
impl<__min_step, __max_step, __tolerance, __attempts, __fixed_step, __error_ctrl, __integration_frame> IntegratorOptionsBuilder<((), __min_step, __max_step, __tolerance, __attempts, __fixed_step, __error_ctrl, __integration_frame)>
pub fn init_step( self, init_step: Duration, ) -> IntegratorOptionsBuilder<((Duration,), __min_step, __max_step, __tolerance, __attempts, __fixed_step, __error_ctrl, __integration_frame)>
Source§impl<__init_step, __max_step, __tolerance, __attempts, __fixed_step, __error_ctrl, __integration_frame> IntegratorOptionsBuilder<(__init_step, (), __max_step, __tolerance, __attempts, __fixed_step, __error_ctrl, __integration_frame)>
impl<__init_step, __max_step, __tolerance, __attempts, __fixed_step, __error_ctrl, __integration_frame> IntegratorOptionsBuilder<(__init_step, (), __max_step, __tolerance, __attempts, __fixed_step, __error_ctrl, __integration_frame)>
Sourcepub fn min_step(
self,
min_step: Duration,
) -> IntegratorOptionsBuilder<(__init_step, (Duration,), __max_step, __tolerance, __attempts, __fixed_step, __error_ctrl, __integration_frame)>
pub fn min_step( self, min_step: Duration, ) -> IntegratorOptionsBuilder<(__init_step, (Duration,), __max_step, __tolerance, __attempts, __fixed_step, __error_ctrl, __integration_frame)>
Examples found in repository?
examples/03_geo_analysis/stationkeeping.rs (line 106)
28fn main() -> Result<(), Box<dyn Error>> {
29 pel::init();
30 // Set up the dynamics like in the orbit raise.
31 let almanac = Arc::new(MetaAlmanac::latest().map_err(Box::new)?);
32 let epoch = Epoch::from_gregorian_utc_hms(2024, 2, 29, 12, 13, 14);
33
34 // Define the GEO orbit, and we're just going to maintain it very tightly.
35 let earth_j2000 = almanac.frame_from_uid(EARTH_J2000)?;
36 let orbit = Orbit::try_keplerian(42164.0, 1e-5, 0., 163.0, 75.0, 0.0, epoch, earth_j2000)?;
37 println!("{orbit:x}");
38
39 let sc = Spacecraft::builder()
40 .orbit(orbit)
41 .mass(Mass::from_dry_and_prop_masses(1000.0, 1000.0)) // 1000 kg of dry mass and prop, totalling 2.0 tons
42 .srp(SRPData::from_area(3.0 * 6.0)) // Assuming 1 kW/m^2 or 18 kW, giving a margin of 4.35 kW for on-propulsion consumption
43 .thruster(Thruster {
44 // "NEXT-STEP" row in Table 2
45 isp_s: 4435.0,
46 thrust_N: 0.472,
47 })
48 .mode(GuidanceMode::Thrust) // Start thrusting immediately.
49 .build();
50
51 // Set up the spacecraft dynamics like in the orbit raise example.
52
53 let prop_time = 30.0 * Unit::Day;
54
55 // Define the guidance law -- we're just using a Ruggiero controller as demonstrated in AAS-2004-5089.
56 let objectives = &[
57 Objective::within_tolerance(StateParameter::SMA, 42_164.0, 5.0), // 5 km
58 Objective::within_tolerance(StateParameter::Eccentricity, 0.001, 5e-5),
59 Objective::within_tolerance(StateParameter::Inclination, 0.05, 1e-2),
60 ];
61
62 let ruggiero_ctrl = Ruggiero::from_max_eclipse(objectives, sc, 0.2)?;
63 println!("{ruggiero_ctrl}");
64
65 let mut orbital_dyn = OrbitalDynamics::point_masses(vec![MOON, SUN]);
66
67 let mut jgm3_meta = MetaFile {
68 uri: "http://public-data.nyxspace.com/nyx/models/JGM3.cof.gz".to_string(),
69 crc32: Some(0xF446F027), // Specifying the CRC32 avoids redownloading it if it's cached.
70 };
71 jgm3_meta.process(true)?;
72
73 let harmonics = Harmonics::from_stor(
74 almanac.frame_from_uid(IAU_EARTH_FRAME)?,
75 HarmonicsMem::from_cof(&jgm3_meta.uri, 8, 8, true)?,
76 );
77 orbital_dyn.accel_models.push(harmonics);
78
79 let srp_dyn = SolarPressure::default(EARTH_J2000, almanac.clone())?;
80 let sc_dynamics = SpacecraftDynamics::from_model(orbital_dyn, srp_dyn)
81 .with_guidance_law(ruggiero_ctrl.clone());
82
83 println!("{sc_dynamics}");
84
85 // Finally, let's use the Monte Carlo framework built into Nyx to propagate spacecraft.
86
87 // Let's start by defining the dispersion.
88 // The MultivariateNormal structure allows us to define the dispersions in any of the orbital parameters, but these are applied directly in the Cartesian state space.
89 // Note that additional validation on the MVN is in progress -- https://github.com/nyx-space/nyx/issues/339.
90 let mc_rv = MvnSpacecraft::new(
91 sc,
92 vec![StateDispersion::zero_mean(StateParameter::SMA, 3.0)],
93 )?;
94
95 let my_mc = MonteCarlo::new(
96 sc, // Nominal state
97 mc_rv,
98 "03_geo_sk".to_string(), // Scenario name
99 None, // No specific seed specified, so one will be drawn from the computer's entropy.
100 );
101
102 // Build the propagator setup.
103 let setup = Propagator::rk89(
104 sc_dynamics.clone(),
105 IntegratorOptions::builder()
106 .min_step(10.0_f64.seconds())
107 .error_ctrl(ErrorControl::RSSCartesianStep)
108 .build(),
109 );
110
111 let num_runs = 25;
112 let rslts = my_mc.run_until_epoch(setup, almanac.clone(), sc.epoch() + prop_time, num_runs);
113
114 assert_eq!(rslts.runs.len(), num_runs);
115
116 // For all of the resulting trajectories, we'll want to compute the percentage of penumbra and umbra.
117
118 rslts.to_parquet(
119 "03_geo_sk.parquet",
120 Some(vec![
121 &EclipseLocator::cislunar(almanac.clone()).to_penumbra_event()
122 ]),
123 ExportCfg::default(),
124 almanac,
125 )?;
126
127 Ok(())
128}
More examples
examples/03_geo_analysis/raise.rs (line 120)
27fn main() -> Result<(), Box<dyn Error>> {
28 pel::init();
29
30 // Dynamics models require planetary constants and ephemerides to be defined.
31 // Let's start by grabbing those by using ANISE's latest MetaAlmanac.
32 // This will automatically download the DE440s planetary ephemeris,
33 // the daily-updated Earth Orientation Parameters, the high fidelity Moon orientation
34 // parameters (for the Moon Mean Earth and Moon Principal Axes frames), and the PCK11
35 // planetary constants kernels.
36 // For details, refer to https://github.com/nyx-space/anise/blob/master/data/latest.dhall.
37 // Note that we place the Almanac into an Arc so we can clone it cheaply and provide read-only
38 // references to many functions.
39 let almanac = Arc::new(MetaAlmanac::latest().map_err(Box::new)?);
40 // Fetch the EME2000 frame from the Almabac
41 let eme2k = almanac.frame_from_uid(EARTH_J2000).unwrap();
42 // Define the orbit epoch
43 let epoch = Epoch::from_gregorian_utc_hms(2024, 2, 29, 12, 13, 14);
44
45 // Build the spacecraft itself.
46 // Using slide 6 of https://aerospace.org/sites/default/files/2018-11/Davis-Mayberry_HPSEP_11212018.pdf
47 // for the "next gen" SEP characteristics.
48
49 // GTO start
50 let orbit = Orbit::keplerian(24505.9, 0.725, 7.05, 0.0, 0.0, 0.0, epoch, eme2k);
51
52 let sc = Spacecraft::builder()
53 .orbit(orbit)
54 .mass(Mass::from_dry_and_prop_masses(1000.0, 1000.0)) // 1000 kg of dry mass and prop, totalling 2.0 tons
55 .srp(SRPData::from_area(3.0 * 6.0)) // Assuming 1 kW/m^2 or 18 kW, giving a margin of 4.35 kW for on-propulsion consumption
56 .thruster(Thruster {
57 // "NEXT-STEP" row in Table 2
58 isp_s: 4435.0,
59 thrust_N: 0.472,
60 })
61 .mode(GuidanceMode::Thrust) // Start thrusting immediately.
62 .build();
63
64 let prop_time = 180.0 * Unit::Day;
65
66 // Define the guidance law -- we're just using a Ruggiero controller as demonstrated in AAS-2004-5089.
67 let objectives = &[
68 Objective::within_tolerance(StateParameter::SMA, 42_165.0, 20.0),
69 Objective::within_tolerance(StateParameter::Eccentricity, 0.001, 5e-5),
70 Objective::within_tolerance(StateParameter::Inclination, 0.05, 1e-2),
71 ];
72
73 // Ensure that we only thrust if we have more than 20% illumination.
74 let ruggiero_ctrl = Ruggiero::from_max_eclipse(objectives, sc, 0.2).unwrap();
75 println!("{ruggiero_ctrl}");
76
77 // Define the high fidelity dynamics
78
79 // Set up the spacecraft dynamics.
80
81 // Specify that the orbital dynamics must account for the graviational pull of the Moon and the Sun.
82 // The gravity of the Earth will also be accounted for since the spaceraft in an Earth orbit.
83 let mut orbital_dyn = OrbitalDynamics::point_masses(vec![MOON, SUN]);
84
85 // We want to include the spherical harmonics, so let's download the gravitational data from the Nyx Cloud.
86 // We're using the JGM3 model here, which is the default in GMAT.
87 let mut jgm3_meta = MetaFile {
88 uri: "http://public-data.nyxspace.com/nyx/models/JGM3.cof.gz".to_string(),
89 crc32: Some(0xF446F027), // Specifying the CRC32 avoids redownloading it if it's cached.
90 };
91 // And let's download it if we don't have it yet.
92 jgm3_meta.process(true)?;
93
94 // Build the spherical harmonics.
95 // The harmonics must be computed in the body fixed frame.
96 // We're using the long term prediction of the Earth centered Earth fixed frame, IAU Earth.
97 let harmonics = Harmonics::from_stor(
98 almanac.frame_from_uid(IAU_EARTH_FRAME)?,
99 HarmonicsMem::from_cof(&jgm3_meta.uri, 8, 8, true).unwrap(),
100 );
101
102 // Include the spherical harmonics into the orbital dynamics.
103 orbital_dyn.accel_models.push(harmonics);
104
105 // We define the solar radiation pressure, using the default solar flux and accounting only
106 // for the eclipsing caused by the Earth.
107 let srp_dyn = SolarPressure::default(EARTH_J2000, almanac.clone())?;
108
109 // Finalize setting up the dynamics, specifying the force models (orbital_dyn) separately from the
110 // acceleration models (SRP in this case). Use `from_models` to specify multiple accel models.
111 let sc_dynamics = SpacecraftDynamics::from_model(orbital_dyn, srp_dyn)
112 .with_guidance_law(ruggiero_ctrl.clone());
113
114 println!("{orbit:x}");
115
116 // We specify a minimum step in the propagator because the Ruggiero control would otherwise drive this step very low.
117 let (final_state, traj) = Propagator::rk89(
118 sc_dynamics.clone(),
119 IntegratorOptions::builder()
120 .min_step(10.0_f64.seconds())
121 .error_ctrl(ErrorControl::RSSCartesianStep)
122 .build(),
123 )
124 .with(sc, almanac.clone())
125 .for_duration_with_traj(prop_time)?;
126
127 let prop_usage = sc.mass.prop_mass_kg - final_state.mass.prop_mass_kg;
128 println!("{:x}", final_state.orbit);
129 println!("prop usage: {prop_usage:.3} kg");
130
131 // Finally, export the results for analysis, including the penumbra percentage throughout the orbit raise.
132 traj.to_parquet(
133 "./03_geo_raise.parquet",
134 Some(vec![
135 &EclipseLocator::cislunar(almanac.clone()).to_penumbra_event()
136 ]),
137 ExportCfg::default(),
138 almanac,
139 )?;
140
141 for status_line in ruggiero_ctrl.status(&final_state) {
142 println!("{status_line}");
143 }
144
145 ruggiero_ctrl
146 .achieved(&final_state)
147 .expect("objective not achieved");
148
149 Ok(())
150}
Source§impl<__init_step, __min_step, __tolerance, __attempts, __fixed_step, __error_ctrl, __integration_frame> IntegratorOptionsBuilder<(__init_step, __min_step, (), __tolerance, __attempts, __fixed_step, __error_ctrl, __integration_frame)>
impl<__init_step, __min_step, __tolerance, __attempts, __fixed_step, __error_ctrl, __integration_frame> IntegratorOptionsBuilder<(__init_step, __min_step, (), __tolerance, __attempts, __fixed_step, __error_ctrl, __integration_frame)>
Sourcepub fn max_step(
self,
max_step: Duration,
) -> IntegratorOptionsBuilder<(__init_step, __min_step, (Duration,), __tolerance, __attempts, __fixed_step, __error_ctrl, __integration_frame)>
pub fn max_step( self, max_step: Duration, ) -> IntegratorOptionsBuilder<(__init_step, __min_step, (Duration,), __tolerance, __attempts, __fixed_step, __error_ctrl, __integration_frame)>
Examples found in repository?
examples/05_cislunar_spacecraft_link_od/main.rs (line 86)
34fn main() -> Result<(), Box<dyn Error>> {
35 pel::init();
36
37 // ====================== //
38 // === ALMANAC SET UP === //
39 // ====================== //
40
41 let manifest_dir =
42 PathBuf::from(std::env::var("CARGO_MANIFEST_DIR").unwrap_or(".".to_string()));
43
44 let out = manifest_dir.join("data/04_output/");
45
46 let almanac = Arc::new(
47 Almanac::new(
48 &manifest_dir
49 .join("data/01_planetary/pck08.pca")
50 .to_string_lossy(),
51 )
52 .unwrap()
53 .load(
54 &manifest_dir
55 .join("data/01_planetary/de440s.bsp")
56 .to_string_lossy(),
57 )
58 .unwrap(),
59 );
60
61 let eme2k = almanac.frame_from_uid(EARTH_J2000).unwrap();
62 let moon_iau = almanac.frame_from_uid(IAU_MOON_FRAME).unwrap();
63
64 let epoch = Epoch::from_gregorian_tai(2021, 5, 29, 19, 51, 16, 852_000);
65 let nrho = Orbit::cartesian(
66 166_473.631_302_239_7,
67 -274_715.487_253_382_7,
68 -211_233.210_176_686_7,
69 0.933_451_604_520_018_4,
70 0.436_775_046_841_900_9,
71 -0.082_211_021_250_348_95,
72 epoch,
73 eme2k,
74 );
75
76 let tx_nrho_sc = Spacecraft::from(nrho);
77
78 let state_luna = almanac.transform_to(nrho, MOON_J2000, None).unwrap();
79 println!("Start state (dynamics: Earth, Moon, Sun gravity):\n{state_luna}");
80
81 let bodies = vec![EARTH, SUN];
82 let dynamics = SpacecraftDynamics::new(OrbitalDynamics::point_masses(bodies));
83
84 let setup = Propagator::rk89(
85 dynamics,
86 IntegratorOptions::builder().max_step(0.5.minutes()).build(),
87 );
88
89 /* == Propagate the NRHO vehicle == */
90 let prop_time = 1.1 * state_luna.period().unwrap();
91
92 let (nrho_final, mut tx_traj) = setup
93 .with(tx_nrho_sc, almanac.clone())
94 .for_duration_with_traj(prop_time)
95 .unwrap();
96
97 tx_traj.name = Some("NRHO Tx SC".to_string());
98
99 println!("{tx_traj}");
100
101 /* == Propagate an LLO vehicle == */
102 let llo_orbit =
103 Orbit::try_keplerian_altitude(110.0, 1e-4, 90.0, 0.0, 0.0, 0.0, epoch, moon_iau).unwrap();
104
105 let llo_sc = Spacecraft::builder().orbit(llo_orbit).build();
106
107 let (_, llo_traj) = setup
108 .with(llo_sc, almanac.clone())
109 .until_epoch_with_traj(nrho_final.epoch())
110 .unwrap();
111
112 // Export the subset of the first two hours.
113 llo_traj
114 .clone()
115 .filter_by_offset(..2.hours())
116 .to_parquet_simple(out.join("05_caps_llo_truth.pq"), almanac.clone())?;
117
118 /* == Setup the interlink == */
119
120 let mut measurement_types = IndexSet::new();
121 measurement_types.insert(MeasurementType::Range);
122 measurement_types.insert(MeasurementType::Doppler);
123
124 let mut stochastics = IndexMap::new();
125
126 let sa45_csac_allan_dev = 1e-11;
127
128 stochastics.insert(
129 MeasurementType::Range,
130 StochasticNoise::from_hardware_range_km(
131 sa45_csac_allan_dev,
132 10.0.seconds(),
133 link_specific::ChipRate::StandardT4B,
134 link_specific::SN0::Average,
135 ),
136 );
137
138 stochastics.insert(
139 MeasurementType::Doppler,
140 StochasticNoise::from_hardware_doppler_km_s(
141 sa45_csac_allan_dev,
142 10.0.seconds(),
143 link_specific::CarrierFreq::SBand,
144 link_specific::CN0::Average,
145 ),
146 );
147
148 let interlink = InterlinkTxSpacecraft {
149 traj: tx_traj,
150 measurement_types,
151 integration_time: None,
152 timestamp_noise_s: None,
153 ab_corr: Aberration::LT,
154 stochastic_noises: Some(stochastics),
155 };
156
157 // Devices are the transmitter, which is our NRHO vehicle.
158 let mut devices = BTreeMap::new();
159 devices.insert("NRHO Tx SC".to_string(), interlink);
160
161 let mut configs = BTreeMap::new();
162 configs.insert(
163 "NRHO Tx SC".to_string(),
164 TrkConfig::builder()
165 .strands(vec![Strand {
166 start: epoch,
167 end: nrho_final.epoch(),
168 }])
169 .build(),
170 );
171
172 let mut trk_sim =
173 TrackingArcSim::with_seed(devices.clone(), llo_traj.clone(), configs, 0).unwrap();
174 println!("{trk_sim}");
175
176 let trk_data = trk_sim.generate_measurements(almanac.clone()).unwrap();
177 println!("{trk_data}");
178
179 trk_data
180 .to_parquet_simple(out.clone().join("nrho_interlink_msr.pq"))
181 .unwrap();
182
183 // Run a truth OD where we estimate the LLO position
184 let llo_uncertainty = SpacecraftUncertainty::builder()
185 .nominal(llo_sc)
186 .x_km(1.0)
187 .y_km(1.0)
188 .z_km(1.0)
189 .vx_km_s(1e-3)
190 .vy_km_s(1e-3)
191 .vz_km_s(1e-3)
192 .build();
193
194 let mut proc_devices = devices.clone();
195
196 // Define the initial estimate, randomized, seed for reproducibility
197 let mut initial_estimate = llo_uncertainty.to_estimate_randomized(Some(0)).unwrap();
198 // Inflate the covariance -- https://github.com/nyx-space/nyx/issues/339
199 initial_estimate.covar *= 2.5;
200
201 // Increase the noise in the devices to accept more measurements.
202
203 for link in proc_devices.values_mut() {
204 for noise in &mut link.stochastic_noises.as_mut().unwrap().values_mut() {
205 *noise.white_noise.as_mut().unwrap() *= 3.0;
206 }
207 }
208
209 let init_err = initial_estimate
210 .orbital_state()
211 .ric_difference(&llo_orbit)
212 .unwrap();
213
214 println!("initial estimate:\n{initial_estimate}");
215 println!("RIC errors = {init_err}",);
216
217 let odp = InterlinkKalmanOD::new(
218 setup.clone(),
219 KalmanVariant::ReferenceUpdate,
220 Some(ResidRejectCrit::default()),
221 proc_devices,
222 almanac.clone(),
223 );
224
225 // Shrink the data to process.
226 let arc = trk_data.filter_by_offset(..2.hours());
227
228 let od_sol = odp.process_arc(initial_estimate, &arc).unwrap();
229
230 println!("{od_sol}");
231
232 od_sol
233 .to_parquet(
234 out.join(format!("05_caps_interlink_od_sol.pq")),
235 ExportCfg::default(),
236 )
237 .unwrap();
238
239 let od_traj = od_sol.to_traj().unwrap();
240
241 od_traj
242 .ric_diff_to_parquet(
243 &llo_traj,
244 out.join(format!("05_caps_interlink_llo_est_error.pq")),
245 ExportCfg::default(),
246 )
247 .unwrap();
248
249 let final_est = od_sol.estimates.last().unwrap();
250 assert!(final_est.within_3sigma(), "should be within 3 sigma");
251
252 println!("ESTIMATE\n{final_est:x}\n");
253 let truth = llo_traj.at(final_est.epoch()).unwrap();
254 println!("TRUTH\n{truth:x}");
255
256 let final_err = truth
257 .orbit
258 .ric_difference(&final_est.orbital_state())
259 .unwrap();
260 println!("ERROR {final_err}");
261
262 // Build the residuals versus reference plot.
263 let rvr_sol = odp
264 .process_arc(initial_estimate, &arc.resid_vs_ref_check())
265 .unwrap();
266
267 rvr_sol
268 .to_parquet(
269 out.join(format!("05_caps_interlink_resid_v_ref.pq")),
270 ExportCfg::default(),
271 )
272 .unwrap();
273
274 let final_rvr = rvr_sol.estimates.last().unwrap();
275
276 println!("RMAG error {:.3} m", final_err.rmag_km() * 1e3);
277 println!(
278 "Pure prop error {:.3} m",
279 final_rvr
280 .orbital_state()
281 .ric_difference(&final_est.orbital_state())
282 .unwrap()
283 .rmag_km()
284 * 1e3
285 );
286
287 Ok(())
288}
Source§impl<__init_step, __min_step, __max_step, __attempts, __fixed_step, __error_ctrl, __integration_frame> IntegratorOptionsBuilder<(__init_step, __min_step, __max_step, (), __attempts, __fixed_step, __error_ctrl, __integration_frame)>
impl<__init_step, __min_step, __max_step, __attempts, __fixed_step, __error_ctrl, __integration_frame> IntegratorOptionsBuilder<(__init_step, __min_step, __max_step, (), __attempts, __fixed_step, __error_ctrl, __integration_frame)>
pub fn tolerance( self, tolerance: f64, ) -> IntegratorOptionsBuilder<(__init_step, __min_step, __max_step, (f64,), __attempts, __fixed_step, __error_ctrl, __integration_frame)>
Source§impl<__init_step, __min_step, __max_step, __tolerance, __fixed_step, __error_ctrl, __integration_frame> IntegratorOptionsBuilder<(__init_step, __min_step, __max_step, __tolerance, (), __fixed_step, __error_ctrl, __integration_frame)>
impl<__init_step, __min_step, __max_step, __tolerance, __fixed_step, __error_ctrl, __integration_frame> IntegratorOptionsBuilder<(__init_step, __min_step, __max_step, __tolerance, (), __fixed_step, __error_ctrl, __integration_frame)>
pub fn attempts( self, attempts: u8, ) -> IntegratorOptionsBuilder<(__init_step, __min_step, __max_step, __tolerance, (u8,), __fixed_step, __error_ctrl, __integration_frame)>
Source§impl<__init_step, __min_step, __max_step, __tolerance, __attempts, __error_ctrl, __integration_frame> IntegratorOptionsBuilder<(__init_step, __min_step, __max_step, __tolerance, __attempts, (), __error_ctrl, __integration_frame)>
impl<__init_step, __min_step, __max_step, __tolerance, __attempts, __error_ctrl, __integration_frame> IntegratorOptionsBuilder<(__init_step, __min_step, __max_step, __tolerance, __attempts, (), __error_ctrl, __integration_frame)>
pub fn fixed_step( self, fixed_step: bool, ) -> IntegratorOptionsBuilder<(__init_step, __min_step, __max_step, __tolerance, __attempts, (bool,), __error_ctrl, __integration_frame)>
Source§impl<__init_step, __min_step, __max_step, __tolerance, __attempts, __fixed_step, __integration_frame> IntegratorOptionsBuilder<(__init_step, __min_step, __max_step, __tolerance, __attempts, __fixed_step, (), __integration_frame)>
impl<__init_step, __min_step, __max_step, __tolerance, __attempts, __fixed_step, __integration_frame> IntegratorOptionsBuilder<(__init_step, __min_step, __max_step, __tolerance, __attempts, __fixed_step, (), __integration_frame)>
Sourcepub fn error_ctrl(
self,
error_ctrl: ErrorControl,
) -> IntegratorOptionsBuilder<(__init_step, __min_step, __max_step, __tolerance, __attempts, __fixed_step, (ErrorControl,), __integration_frame)>
pub fn error_ctrl( self, error_ctrl: ErrorControl, ) -> IntegratorOptionsBuilder<(__init_step, __min_step, __max_step, __tolerance, __attempts, __fixed_step, (ErrorControl,), __integration_frame)>
Examples found in repository?
examples/03_geo_analysis/stationkeeping.rs (line 107)
28fn main() -> Result<(), Box<dyn Error>> {
29 pel::init();
30 // Set up the dynamics like in the orbit raise.
31 let almanac = Arc::new(MetaAlmanac::latest().map_err(Box::new)?);
32 let epoch = Epoch::from_gregorian_utc_hms(2024, 2, 29, 12, 13, 14);
33
34 // Define the GEO orbit, and we're just going to maintain it very tightly.
35 let earth_j2000 = almanac.frame_from_uid(EARTH_J2000)?;
36 let orbit = Orbit::try_keplerian(42164.0, 1e-5, 0., 163.0, 75.0, 0.0, epoch, earth_j2000)?;
37 println!("{orbit:x}");
38
39 let sc = Spacecraft::builder()
40 .orbit(orbit)
41 .mass(Mass::from_dry_and_prop_masses(1000.0, 1000.0)) // 1000 kg of dry mass and prop, totalling 2.0 tons
42 .srp(SRPData::from_area(3.0 * 6.0)) // Assuming 1 kW/m^2 or 18 kW, giving a margin of 4.35 kW for on-propulsion consumption
43 .thruster(Thruster {
44 // "NEXT-STEP" row in Table 2
45 isp_s: 4435.0,
46 thrust_N: 0.472,
47 })
48 .mode(GuidanceMode::Thrust) // Start thrusting immediately.
49 .build();
50
51 // Set up the spacecraft dynamics like in the orbit raise example.
52
53 let prop_time = 30.0 * Unit::Day;
54
55 // Define the guidance law -- we're just using a Ruggiero controller as demonstrated in AAS-2004-5089.
56 let objectives = &[
57 Objective::within_tolerance(StateParameter::SMA, 42_164.0, 5.0), // 5 km
58 Objective::within_tolerance(StateParameter::Eccentricity, 0.001, 5e-5),
59 Objective::within_tolerance(StateParameter::Inclination, 0.05, 1e-2),
60 ];
61
62 let ruggiero_ctrl = Ruggiero::from_max_eclipse(objectives, sc, 0.2)?;
63 println!("{ruggiero_ctrl}");
64
65 let mut orbital_dyn = OrbitalDynamics::point_masses(vec![MOON, SUN]);
66
67 let mut jgm3_meta = MetaFile {
68 uri: "http://public-data.nyxspace.com/nyx/models/JGM3.cof.gz".to_string(),
69 crc32: Some(0xF446F027), // Specifying the CRC32 avoids redownloading it if it's cached.
70 };
71 jgm3_meta.process(true)?;
72
73 let harmonics = Harmonics::from_stor(
74 almanac.frame_from_uid(IAU_EARTH_FRAME)?,
75 HarmonicsMem::from_cof(&jgm3_meta.uri, 8, 8, true)?,
76 );
77 orbital_dyn.accel_models.push(harmonics);
78
79 let srp_dyn = SolarPressure::default(EARTH_J2000, almanac.clone())?;
80 let sc_dynamics = SpacecraftDynamics::from_model(orbital_dyn, srp_dyn)
81 .with_guidance_law(ruggiero_ctrl.clone());
82
83 println!("{sc_dynamics}");
84
85 // Finally, let's use the Monte Carlo framework built into Nyx to propagate spacecraft.
86
87 // Let's start by defining the dispersion.
88 // The MultivariateNormal structure allows us to define the dispersions in any of the orbital parameters, but these are applied directly in the Cartesian state space.
89 // Note that additional validation on the MVN is in progress -- https://github.com/nyx-space/nyx/issues/339.
90 let mc_rv = MvnSpacecraft::new(
91 sc,
92 vec![StateDispersion::zero_mean(StateParameter::SMA, 3.0)],
93 )?;
94
95 let my_mc = MonteCarlo::new(
96 sc, // Nominal state
97 mc_rv,
98 "03_geo_sk".to_string(), // Scenario name
99 None, // No specific seed specified, so one will be drawn from the computer's entropy.
100 );
101
102 // Build the propagator setup.
103 let setup = Propagator::rk89(
104 sc_dynamics.clone(),
105 IntegratorOptions::builder()
106 .min_step(10.0_f64.seconds())
107 .error_ctrl(ErrorControl::RSSCartesianStep)
108 .build(),
109 );
110
111 let num_runs = 25;
112 let rslts = my_mc.run_until_epoch(setup, almanac.clone(), sc.epoch() + prop_time, num_runs);
113
114 assert_eq!(rslts.runs.len(), num_runs);
115
116 // For all of the resulting trajectories, we'll want to compute the percentage of penumbra and umbra.
117
118 rslts.to_parquet(
119 "03_geo_sk.parquet",
120 Some(vec![
121 &EclipseLocator::cislunar(almanac.clone()).to_penumbra_event()
122 ]),
123 ExportCfg::default(),
124 almanac,
125 )?;
126
127 Ok(())
128}
More examples
examples/03_geo_analysis/raise.rs (line 121)
27fn main() -> Result<(), Box<dyn Error>> {
28 pel::init();
29
30 // Dynamics models require planetary constants and ephemerides to be defined.
31 // Let's start by grabbing those by using ANISE's latest MetaAlmanac.
32 // This will automatically download the DE440s planetary ephemeris,
33 // the daily-updated Earth Orientation Parameters, the high fidelity Moon orientation
34 // parameters (for the Moon Mean Earth and Moon Principal Axes frames), and the PCK11
35 // planetary constants kernels.
36 // For details, refer to https://github.com/nyx-space/anise/blob/master/data/latest.dhall.
37 // Note that we place the Almanac into an Arc so we can clone it cheaply and provide read-only
38 // references to many functions.
39 let almanac = Arc::new(MetaAlmanac::latest().map_err(Box::new)?);
40 // Fetch the EME2000 frame from the Almabac
41 let eme2k = almanac.frame_from_uid(EARTH_J2000).unwrap();
42 // Define the orbit epoch
43 let epoch = Epoch::from_gregorian_utc_hms(2024, 2, 29, 12, 13, 14);
44
45 // Build the spacecraft itself.
46 // Using slide 6 of https://aerospace.org/sites/default/files/2018-11/Davis-Mayberry_HPSEP_11212018.pdf
47 // for the "next gen" SEP characteristics.
48
49 // GTO start
50 let orbit = Orbit::keplerian(24505.9, 0.725, 7.05, 0.0, 0.0, 0.0, epoch, eme2k);
51
52 let sc = Spacecraft::builder()
53 .orbit(orbit)
54 .mass(Mass::from_dry_and_prop_masses(1000.0, 1000.0)) // 1000 kg of dry mass and prop, totalling 2.0 tons
55 .srp(SRPData::from_area(3.0 * 6.0)) // Assuming 1 kW/m^2 or 18 kW, giving a margin of 4.35 kW for on-propulsion consumption
56 .thruster(Thruster {
57 // "NEXT-STEP" row in Table 2
58 isp_s: 4435.0,
59 thrust_N: 0.472,
60 })
61 .mode(GuidanceMode::Thrust) // Start thrusting immediately.
62 .build();
63
64 let prop_time = 180.0 * Unit::Day;
65
66 // Define the guidance law -- we're just using a Ruggiero controller as demonstrated in AAS-2004-5089.
67 let objectives = &[
68 Objective::within_tolerance(StateParameter::SMA, 42_165.0, 20.0),
69 Objective::within_tolerance(StateParameter::Eccentricity, 0.001, 5e-5),
70 Objective::within_tolerance(StateParameter::Inclination, 0.05, 1e-2),
71 ];
72
73 // Ensure that we only thrust if we have more than 20% illumination.
74 let ruggiero_ctrl = Ruggiero::from_max_eclipse(objectives, sc, 0.2).unwrap();
75 println!("{ruggiero_ctrl}");
76
77 // Define the high fidelity dynamics
78
79 // Set up the spacecraft dynamics.
80
81 // Specify that the orbital dynamics must account for the graviational pull of the Moon and the Sun.
82 // The gravity of the Earth will also be accounted for since the spaceraft in an Earth orbit.
83 let mut orbital_dyn = OrbitalDynamics::point_masses(vec![MOON, SUN]);
84
85 // We want to include the spherical harmonics, so let's download the gravitational data from the Nyx Cloud.
86 // We're using the JGM3 model here, which is the default in GMAT.
87 let mut jgm3_meta = MetaFile {
88 uri: "http://public-data.nyxspace.com/nyx/models/JGM3.cof.gz".to_string(),
89 crc32: Some(0xF446F027), // Specifying the CRC32 avoids redownloading it if it's cached.
90 };
91 // And let's download it if we don't have it yet.
92 jgm3_meta.process(true)?;
93
94 // Build the spherical harmonics.
95 // The harmonics must be computed in the body fixed frame.
96 // We're using the long term prediction of the Earth centered Earth fixed frame, IAU Earth.
97 let harmonics = Harmonics::from_stor(
98 almanac.frame_from_uid(IAU_EARTH_FRAME)?,
99 HarmonicsMem::from_cof(&jgm3_meta.uri, 8, 8, true).unwrap(),
100 );
101
102 // Include the spherical harmonics into the orbital dynamics.
103 orbital_dyn.accel_models.push(harmonics);
104
105 // We define the solar radiation pressure, using the default solar flux and accounting only
106 // for the eclipsing caused by the Earth.
107 let srp_dyn = SolarPressure::default(EARTH_J2000, almanac.clone())?;
108
109 // Finalize setting up the dynamics, specifying the force models (orbital_dyn) separately from the
110 // acceleration models (SRP in this case). Use `from_models` to specify multiple accel models.
111 let sc_dynamics = SpacecraftDynamics::from_model(orbital_dyn, srp_dyn)
112 .with_guidance_law(ruggiero_ctrl.clone());
113
114 println!("{orbit:x}");
115
116 // We specify a minimum step in the propagator because the Ruggiero control would otherwise drive this step very low.
117 let (final_state, traj) = Propagator::rk89(
118 sc_dynamics.clone(),
119 IntegratorOptions::builder()
120 .min_step(10.0_f64.seconds())
121 .error_ctrl(ErrorControl::RSSCartesianStep)
122 .build(),
123 )
124 .with(sc, almanac.clone())
125 .for_duration_with_traj(prop_time)?;
126
127 let prop_usage = sc.mass.prop_mass_kg - final_state.mass.prop_mass_kg;
128 println!("{:x}", final_state.orbit);
129 println!("prop usage: {prop_usage:.3} kg");
130
131 // Finally, export the results for analysis, including the penumbra percentage throughout the orbit raise.
132 traj.to_parquet(
133 "./03_geo_raise.parquet",
134 Some(vec![
135 &EclipseLocator::cislunar(almanac.clone()).to_penumbra_event()
136 ]),
137 ExportCfg::default(),
138 almanac,
139 )?;
140
141 for status_line in ruggiero_ctrl.status(&final_state) {
142 println!("{status_line}");
143 }
144
145 ruggiero_ctrl
146 .achieved(&final_state)
147 .expect("objective not achieved");
148
149 Ok(())
150}
Source§impl<__init_step, __min_step, __max_step, __tolerance, __attempts, __fixed_step, __error_ctrl> IntegratorOptionsBuilder<(__init_step, __min_step, __max_step, __tolerance, __attempts, __fixed_step, __error_ctrl, ())>
impl<__init_step, __min_step, __max_step, __tolerance, __attempts, __fixed_step, __error_ctrl> IntegratorOptionsBuilder<(__init_step, __min_step, __max_step, __tolerance, __attempts, __fixed_step, __error_ctrl, ())>
Sourcepub fn integration_frame(
self,
integration_frame: Frame,
) -> IntegratorOptionsBuilder<(__init_step, __min_step, __max_step, __tolerance, __attempts, __fixed_step, __error_ctrl, (Option<Frame>,))>
pub fn integration_frame( self, integration_frame: Frame, ) -> IntegratorOptionsBuilder<(__init_step, __min_step, __max_step, __tolerance, __attempts, __fixed_step, __error_ctrl, (Option<Frame>,))>
If a frame is specified and the propagator state is in a different frame, it it changed to this frame prior to integration.
Note, when setting this, it’s recommended to call strip
on the Frame.
Source§impl<__init_step: Optional<Duration>, __min_step: Optional<Duration>, __max_step: Optional<Duration>, __tolerance: Optional<f64>, __attempts: Optional<u8>, __fixed_step: Optional<bool>, __error_ctrl: Optional<ErrorControl>, __integration_frame: Optional<Option<Frame>>> IntegratorOptionsBuilder<(__init_step, __min_step, __max_step, __tolerance, __attempts, __fixed_step, __error_ctrl, __integration_frame)>
impl<__init_step: Optional<Duration>, __min_step: Optional<Duration>, __max_step: Optional<Duration>, __tolerance: Optional<f64>, __attempts: Optional<u8>, __fixed_step: Optional<bool>, __error_ctrl: Optional<ErrorControl>, __integration_frame: Optional<Option<Frame>>> IntegratorOptionsBuilder<(__init_step, __min_step, __max_step, __tolerance, __attempts, __fixed_step, __error_ctrl, __integration_frame)>
Sourcepub fn build(self) -> IntegratorOptions
pub fn build(self) -> IntegratorOptions
Finalise the builder and create its IntegratorOptions
instance
Examples found in repository?
examples/03_geo_analysis/stationkeeping.rs (line 108)
28fn main() -> Result<(), Box<dyn Error>> {
29 pel::init();
30 // Set up the dynamics like in the orbit raise.
31 let almanac = Arc::new(MetaAlmanac::latest().map_err(Box::new)?);
32 let epoch = Epoch::from_gregorian_utc_hms(2024, 2, 29, 12, 13, 14);
33
34 // Define the GEO orbit, and we're just going to maintain it very tightly.
35 let earth_j2000 = almanac.frame_from_uid(EARTH_J2000)?;
36 let orbit = Orbit::try_keplerian(42164.0, 1e-5, 0., 163.0, 75.0, 0.0, epoch, earth_j2000)?;
37 println!("{orbit:x}");
38
39 let sc = Spacecraft::builder()
40 .orbit(orbit)
41 .mass(Mass::from_dry_and_prop_masses(1000.0, 1000.0)) // 1000 kg of dry mass and prop, totalling 2.0 tons
42 .srp(SRPData::from_area(3.0 * 6.0)) // Assuming 1 kW/m^2 or 18 kW, giving a margin of 4.35 kW for on-propulsion consumption
43 .thruster(Thruster {
44 // "NEXT-STEP" row in Table 2
45 isp_s: 4435.0,
46 thrust_N: 0.472,
47 })
48 .mode(GuidanceMode::Thrust) // Start thrusting immediately.
49 .build();
50
51 // Set up the spacecraft dynamics like in the orbit raise example.
52
53 let prop_time = 30.0 * Unit::Day;
54
55 // Define the guidance law -- we're just using a Ruggiero controller as demonstrated in AAS-2004-5089.
56 let objectives = &[
57 Objective::within_tolerance(StateParameter::SMA, 42_164.0, 5.0), // 5 km
58 Objective::within_tolerance(StateParameter::Eccentricity, 0.001, 5e-5),
59 Objective::within_tolerance(StateParameter::Inclination, 0.05, 1e-2),
60 ];
61
62 let ruggiero_ctrl = Ruggiero::from_max_eclipse(objectives, sc, 0.2)?;
63 println!("{ruggiero_ctrl}");
64
65 let mut orbital_dyn = OrbitalDynamics::point_masses(vec![MOON, SUN]);
66
67 let mut jgm3_meta = MetaFile {
68 uri: "http://public-data.nyxspace.com/nyx/models/JGM3.cof.gz".to_string(),
69 crc32: Some(0xF446F027), // Specifying the CRC32 avoids redownloading it if it's cached.
70 };
71 jgm3_meta.process(true)?;
72
73 let harmonics = Harmonics::from_stor(
74 almanac.frame_from_uid(IAU_EARTH_FRAME)?,
75 HarmonicsMem::from_cof(&jgm3_meta.uri, 8, 8, true)?,
76 );
77 orbital_dyn.accel_models.push(harmonics);
78
79 let srp_dyn = SolarPressure::default(EARTH_J2000, almanac.clone())?;
80 let sc_dynamics = SpacecraftDynamics::from_model(orbital_dyn, srp_dyn)
81 .with_guidance_law(ruggiero_ctrl.clone());
82
83 println!("{sc_dynamics}");
84
85 // Finally, let's use the Monte Carlo framework built into Nyx to propagate spacecraft.
86
87 // Let's start by defining the dispersion.
88 // The MultivariateNormal structure allows us to define the dispersions in any of the orbital parameters, but these are applied directly in the Cartesian state space.
89 // Note that additional validation on the MVN is in progress -- https://github.com/nyx-space/nyx/issues/339.
90 let mc_rv = MvnSpacecraft::new(
91 sc,
92 vec![StateDispersion::zero_mean(StateParameter::SMA, 3.0)],
93 )?;
94
95 let my_mc = MonteCarlo::new(
96 sc, // Nominal state
97 mc_rv,
98 "03_geo_sk".to_string(), // Scenario name
99 None, // No specific seed specified, so one will be drawn from the computer's entropy.
100 );
101
102 // Build the propagator setup.
103 let setup = Propagator::rk89(
104 sc_dynamics.clone(),
105 IntegratorOptions::builder()
106 .min_step(10.0_f64.seconds())
107 .error_ctrl(ErrorControl::RSSCartesianStep)
108 .build(),
109 );
110
111 let num_runs = 25;
112 let rslts = my_mc.run_until_epoch(setup, almanac.clone(), sc.epoch() + prop_time, num_runs);
113
114 assert_eq!(rslts.runs.len(), num_runs);
115
116 // For all of the resulting trajectories, we'll want to compute the percentage of penumbra and umbra.
117
118 rslts.to_parquet(
119 "03_geo_sk.parquet",
120 Some(vec![
121 &EclipseLocator::cislunar(almanac.clone()).to_penumbra_event()
122 ]),
123 ExportCfg::default(),
124 almanac,
125 )?;
126
127 Ok(())
128}
More examples
examples/03_geo_analysis/raise.rs (line 122)
27fn main() -> Result<(), Box<dyn Error>> {
28 pel::init();
29
30 // Dynamics models require planetary constants and ephemerides to be defined.
31 // Let's start by grabbing those by using ANISE's latest MetaAlmanac.
32 // This will automatically download the DE440s planetary ephemeris,
33 // the daily-updated Earth Orientation Parameters, the high fidelity Moon orientation
34 // parameters (for the Moon Mean Earth and Moon Principal Axes frames), and the PCK11
35 // planetary constants kernels.
36 // For details, refer to https://github.com/nyx-space/anise/blob/master/data/latest.dhall.
37 // Note that we place the Almanac into an Arc so we can clone it cheaply and provide read-only
38 // references to many functions.
39 let almanac = Arc::new(MetaAlmanac::latest().map_err(Box::new)?);
40 // Fetch the EME2000 frame from the Almabac
41 let eme2k = almanac.frame_from_uid(EARTH_J2000).unwrap();
42 // Define the orbit epoch
43 let epoch = Epoch::from_gregorian_utc_hms(2024, 2, 29, 12, 13, 14);
44
45 // Build the spacecraft itself.
46 // Using slide 6 of https://aerospace.org/sites/default/files/2018-11/Davis-Mayberry_HPSEP_11212018.pdf
47 // for the "next gen" SEP characteristics.
48
49 // GTO start
50 let orbit = Orbit::keplerian(24505.9, 0.725, 7.05, 0.0, 0.0, 0.0, epoch, eme2k);
51
52 let sc = Spacecraft::builder()
53 .orbit(orbit)
54 .mass(Mass::from_dry_and_prop_masses(1000.0, 1000.0)) // 1000 kg of dry mass and prop, totalling 2.0 tons
55 .srp(SRPData::from_area(3.0 * 6.0)) // Assuming 1 kW/m^2 or 18 kW, giving a margin of 4.35 kW for on-propulsion consumption
56 .thruster(Thruster {
57 // "NEXT-STEP" row in Table 2
58 isp_s: 4435.0,
59 thrust_N: 0.472,
60 })
61 .mode(GuidanceMode::Thrust) // Start thrusting immediately.
62 .build();
63
64 let prop_time = 180.0 * Unit::Day;
65
66 // Define the guidance law -- we're just using a Ruggiero controller as demonstrated in AAS-2004-5089.
67 let objectives = &[
68 Objective::within_tolerance(StateParameter::SMA, 42_165.0, 20.0),
69 Objective::within_tolerance(StateParameter::Eccentricity, 0.001, 5e-5),
70 Objective::within_tolerance(StateParameter::Inclination, 0.05, 1e-2),
71 ];
72
73 // Ensure that we only thrust if we have more than 20% illumination.
74 let ruggiero_ctrl = Ruggiero::from_max_eclipse(objectives, sc, 0.2).unwrap();
75 println!("{ruggiero_ctrl}");
76
77 // Define the high fidelity dynamics
78
79 // Set up the spacecraft dynamics.
80
81 // Specify that the orbital dynamics must account for the graviational pull of the Moon and the Sun.
82 // The gravity of the Earth will also be accounted for since the spaceraft in an Earth orbit.
83 let mut orbital_dyn = OrbitalDynamics::point_masses(vec![MOON, SUN]);
84
85 // We want to include the spherical harmonics, so let's download the gravitational data from the Nyx Cloud.
86 // We're using the JGM3 model here, which is the default in GMAT.
87 let mut jgm3_meta = MetaFile {
88 uri: "http://public-data.nyxspace.com/nyx/models/JGM3.cof.gz".to_string(),
89 crc32: Some(0xF446F027), // Specifying the CRC32 avoids redownloading it if it's cached.
90 };
91 // And let's download it if we don't have it yet.
92 jgm3_meta.process(true)?;
93
94 // Build the spherical harmonics.
95 // The harmonics must be computed in the body fixed frame.
96 // We're using the long term prediction of the Earth centered Earth fixed frame, IAU Earth.
97 let harmonics = Harmonics::from_stor(
98 almanac.frame_from_uid(IAU_EARTH_FRAME)?,
99 HarmonicsMem::from_cof(&jgm3_meta.uri, 8, 8, true).unwrap(),
100 );
101
102 // Include the spherical harmonics into the orbital dynamics.
103 orbital_dyn.accel_models.push(harmonics);
104
105 // We define the solar radiation pressure, using the default solar flux and accounting only
106 // for the eclipsing caused by the Earth.
107 let srp_dyn = SolarPressure::default(EARTH_J2000, almanac.clone())?;
108
109 // Finalize setting up the dynamics, specifying the force models (orbital_dyn) separately from the
110 // acceleration models (SRP in this case). Use `from_models` to specify multiple accel models.
111 let sc_dynamics = SpacecraftDynamics::from_model(orbital_dyn, srp_dyn)
112 .with_guidance_law(ruggiero_ctrl.clone());
113
114 println!("{orbit:x}");
115
116 // We specify a minimum step in the propagator because the Ruggiero control would otherwise drive this step very low.
117 let (final_state, traj) = Propagator::rk89(
118 sc_dynamics.clone(),
119 IntegratorOptions::builder()
120 .min_step(10.0_f64.seconds())
121 .error_ctrl(ErrorControl::RSSCartesianStep)
122 .build(),
123 )
124 .with(sc, almanac.clone())
125 .for_duration_with_traj(prop_time)?;
126
127 let prop_usage = sc.mass.prop_mass_kg - final_state.mass.prop_mass_kg;
128 println!("{:x}", final_state.orbit);
129 println!("prop usage: {prop_usage:.3} kg");
130
131 // Finally, export the results for analysis, including the penumbra percentage throughout the orbit raise.
132 traj.to_parquet(
133 "./03_geo_raise.parquet",
134 Some(vec![
135 &EclipseLocator::cislunar(almanac.clone()).to_penumbra_event()
136 ]),
137 ExportCfg::default(),
138 almanac,
139 )?;
140
141 for status_line in ruggiero_ctrl.status(&final_state) {
142 println!("{status_line}");
143 }
144
145 ruggiero_ctrl
146 .achieved(&final_state)
147 .expect("objective not achieved");
148
149 Ok(())
150}
examples/05_cislunar_spacecraft_link_od/main.rs (line 86)
34fn main() -> Result<(), Box<dyn Error>> {
35 pel::init();
36
37 // ====================== //
38 // === ALMANAC SET UP === //
39 // ====================== //
40
41 let manifest_dir =
42 PathBuf::from(std::env::var("CARGO_MANIFEST_DIR").unwrap_or(".".to_string()));
43
44 let out = manifest_dir.join("data/04_output/");
45
46 let almanac = Arc::new(
47 Almanac::new(
48 &manifest_dir
49 .join("data/01_planetary/pck08.pca")
50 .to_string_lossy(),
51 )
52 .unwrap()
53 .load(
54 &manifest_dir
55 .join("data/01_planetary/de440s.bsp")
56 .to_string_lossy(),
57 )
58 .unwrap(),
59 );
60
61 let eme2k = almanac.frame_from_uid(EARTH_J2000).unwrap();
62 let moon_iau = almanac.frame_from_uid(IAU_MOON_FRAME).unwrap();
63
64 let epoch = Epoch::from_gregorian_tai(2021, 5, 29, 19, 51, 16, 852_000);
65 let nrho = Orbit::cartesian(
66 166_473.631_302_239_7,
67 -274_715.487_253_382_7,
68 -211_233.210_176_686_7,
69 0.933_451_604_520_018_4,
70 0.436_775_046_841_900_9,
71 -0.082_211_021_250_348_95,
72 epoch,
73 eme2k,
74 );
75
76 let tx_nrho_sc = Spacecraft::from(nrho);
77
78 let state_luna = almanac.transform_to(nrho, MOON_J2000, None).unwrap();
79 println!("Start state (dynamics: Earth, Moon, Sun gravity):\n{state_luna}");
80
81 let bodies = vec![EARTH, SUN];
82 let dynamics = SpacecraftDynamics::new(OrbitalDynamics::point_masses(bodies));
83
84 let setup = Propagator::rk89(
85 dynamics,
86 IntegratorOptions::builder().max_step(0.5.minutes()).build(),
87 );
88
89 /* == Propagate the NRHO vehicle == */
90 let prop_time = 1.1 * state_luna.period().unwrap();
91
92 let (nrho_final, mut tx_traj) = setup
93 .with(tx_nrho_sc, almanac.clone())
94 .for_duration_with_traj(prop_time)
95 .unwrap();
96
97 tx_traj.name = Some("NRHO Tx SC".to_string());
98
99 println!("{tx_traj}");
100
101 /* == Propagate an LLO vehicle == */
102 let llo_orbit =
103 Orbit::try_keplerian_altitude(110.0, 1e-4, 90.0, 0.0, 0.0, 0.0, epoch, moon_iau).unwrap();
104
105 let llo_sc = Spacecraft::builder().orbit(llo_orbit).build();
106
107 let (_, llo_traj) = setup
108 .with(llo_sc, almanac.clone())
109 .until_epoch_with_traj(nrho_final.epoch())
110 .unwrap();
111
112 // Export the subset of the first two hours.
113 llo_traj
114 .clone()
115 .filter_by_offset(..2.hours())
116 .to_parquet_simple(out.join("05_caps_llo_truth.pq"), almanac.clone())?;
117
118 /* == Setup the interlink == */
119
120 let mut measurement_types = IndexSet::new();
121 measurement_types.insert(MeasurementType::Range);
122 measurement_types.insert(MeasurementType::Doppler);
123
124 let mut stochastics = IndexMap::new();
125
126 let sa45_csac_allan_dev = 1e-11;
127
128 stochastics.insert(
129 MeasurementType::Range,
130 StochasticNoise::from_hardware_range_km(
131 sa45_csac_allan_dev,
132 10.0.seconds(),
133 link_specific::ChipRate::StandardT4B,
134 link_specific::SN0::Average,
135 ),
136 );
137
138 stochastics.insert(
139 MeasurementType::Doppler,
140 StochasticNoise::from_hardware_doppler_km_s(
141 sa45_csac_allan_dev,
142 10.0.seconds(),
143 link_specific::CarrierFreq::SBand,
144 link_specific::CN0::Average,
145 ),
146 );
147
148 let interlink = InterlinkTxSpacecraft {
149 traj: tx_traj,
150 measurement_types,
151 integration_time: None,
152 timestamp_noise_s: None,
153 ab_corr: Aberration::LT,
154 stochastic_noises: Some(stochastics),
155 };
156
157 // Devices are the transmitter, which is our NRHO vehicle.
158 let mut devices = BTreeMap::new();
159 devices.insert("NRHO Tx SC".to_string(), interlink);
160
161 let mut configs = BTreeMap::new();
162 configs.insert(
163 "NRHO Tx SC".to_string(),
164 TrkConfig::builder()
165 .strands(vec![Strand {
166 start: epoch,
167 end: nrho_final.epoch(),
168 }])
169 .build(),
170 );
171
172 let mut trk_sim =
173 TrackingArcSim::with_seed(devices.clone(), llo_traj.clone(), configs, 0).unwrap();
174 println!("{trk_sim}");
175
176 let trk_data = trk_sim.generate_measurements(almanac.clone()).unwrap();
177 println!("{trk_data}");
178
179 trk_data
180 .to_parquet_simple(out.clone().join("nrho_interlink_msr.pq"))
181 .unwrap();
182
183 // Run a truth OD where we estimate the LLO position
184 let llo_uncertainty = SpacecraftUncertainty::builder()
185 .nominal(llo_sc)
186 .x_km(1.0)
187 .y_km(1.0)
188 .z_km(1.0)
189 .vx_km_s(1e-3)
190 .vy_km_s(1e-3)
191 .vz_km_s(1e-3)
192 .build();
193
194 let mut proc_devices = devices.clone();
195
196 // Define the initial estimate, randomized, seed for reproducibility
197 let mut initial_estimate = llo_uncertainty.to_estimate_randomized(Some(0)).unwrap();
198 // Inflate the covariance -- https://github.com/nyx-space/nyx/issues/339
199 initial_estimate.covar *= 2.5;
200
201 // Increase the noise in the devices to accept more measurements.
202
203 for link in proc_devices.values_mut() {
204 for noise in &mut link.stochastic_noises.as_mut().unwrap().values_mut() {
205 *noise.white_noise.as_mut().unwrap() *= 3.0;
206 }
207 }
208
209 let init_err = initial_estimate
210 .orbital_state()
211 .ric_difference(&llo_orbit)
212 .unwrap();
213
214 println!("initial estimate:\n{initial_estimate}");
215 println!("RIC errors = {init_err}",);
216
217 let odp = InterlinkKalmanOD::new(
218 setup.clone(),
219 KalmanVariant::ReferenceUpdate,
220 Some(ResidRejectCrit::default()),
221 proc_devices,
222 almanac.clone(),
223 );
224
225 // Shrink the data to process.
226 let arc = trk_data.filter_by_offset(..2.hours());
227
228 let od_sol = odp.process_arc(initial_estimate, &arc).unwrap();
229
230 println!("{od_sol}");
231
232 od_sol
233 .to_parquet(
234 out.join(format!("05_caps_interlink_od_sol.pq")),
235 ExportCfg::default(),
236 )
237 .unwrap();
238
239 let od_traj = od_sol.to_traj().unwrap();
240
241 od_traj
242 .ric_diff_to_parquet(
243 &llo_traj,
244 out.join(format!("05_caps_interlink_llo_est_error.pq")),
245 ExportCfg::default(),
246 )
247 .unwrap();
248
249 let final_est = od_sol.estimates.last().unwrap();
250 assert!(final_est.within_3sigma(), "should be within 3 sigma");
251
252 println!("ESTIMATE\n{final_est:x}\n");
253 let truth = llo_traj.at(final_est.epoch()).unwrap();
254 println!("TRUTH\n{truth:x}");
255
256 let final_err = truth
257 .orbit
258 .ric_difference(&final_est.orbital_state())
259 .unwrap();
260 println!("ERROR {final_err}");
261
262 // Build the residuals versus reference plot.
263 let rvr_sol = odp
264 .process_arc(initial_estimate, &arc.resid_vs_ref_check())
265 .unwrap();
266
267 rvr_sol
268 .to_parquet(
269 out.join(format!("05_caps_interlink_resid_v_ref.pq")),
270 ExportCfg::default(),
271 )
272 .unwrap();
273
274 let final_rvr = rvr_sol.estimates.last().unwrap();
275
276 println!("RMAG error {:.3} m", final_err.rmag_km() * 1e3);
277 println!(
278 "Pure prop error {:.3} m",
279 final_rvr
280 .orbital_state()
281 .ric_difference(&final_est.orbital_state())
282 .unwrap()
283 .rmag_km()
284 * 1e3
285 );
286
287 Ok(())
288}
Trait Implementations§
Auto Trait Implementations§
impl<TypedBuilderFields> Freeze for IntegratorOptionsBuilder<TypedBuilderFields>where
TypedBuilderFields: Freeze,
impl<TypedBuilderFields> RefUnwindSafe for IntegratorOptionsBuilder<TypedBuilderFields>where
TypedBuilderFields: RefUnwindSafe,
impl<TypedBuilderFields> Send for IntegratorOptionsBuilder<TypedBuilderFields>where
TypedBuilderFields: Send,
impl<TypedBuilderFields> Sync for IntegratorOptionsBuilder<TypedBuilderFields>where
TypedBuilderFields: Sync,
impl<TypedBuilderFields> Unpin for IntegratorOptionsBuilder<TypedBuilderFields>where
TypedBuilderFields: Unpin,
impl<TypedBuilderFields> UnwindSafe for IntegratorOptionsBuilder<TypedBuilderFields>where
TypedBuilderFields: UnwindSafe,
Blanket Implementations§
Source§impl<T> BorrowMut<T> for Twhere
T: ?Sized,
impl<T> BorrowMut<T> for Twhere
T: ?Sized,
Source§fn borrow_mut(&mut self) -> &mut T
fn borrow_mut(&mut self) -> &mut T
Mutably borrows from an owned value. Read more
Source§impl<T> CloneToUninit for Twhere
T: Clone,
impl<T> CloneToUninit for Twhere
T: Clone,
§impl<T> Instrument for T
impl<T> Instrument for T
§fn instrument(self, span: Span) -> Instrumented<Self>
fn instrument(self, span: Span) -> Instrumented<Self>
§fn in_current_span(self) -> Instrumented<Self>
fn in_current_span(self) -> Instrumented<Self>
Source§impl<T> IntoEither for T
impl<T> IntoEither for T
Source§fn into_either(self, into_left: bool) -> Either<Self, Self>
fn into_either(self, into_left: bool) -> Either<Self, Self>
Converts
self
into a Left
variant of Either<Self, Self>
if into_left
is true
.
Converts self
into a Right
variant of Either<Self, Self>
otherwise. Read moreSource§fn into_either_with<F>(self, into_left: F) -> Either<Self, Self>
fn into_either_with<F>(self, into_left: F) -> Either<Self, Self>
Converts
self
into a Left
variant of Either<Self, Self>
if into_left(&self)
returns true
.
Converts self
into a Right
variant of Either<Self, Self>
otherwise. Read more§impl<T> Pointable for T
impl<T> Pointable for T
§impl<SS, SP> SupersetOf<SS> for SPwhere
SS: SubsetOf<SP>,
impl<SS, SP> SupersetOf<SS> for SPwhere
SS: SubsetOf<SP>,
§fn to_subset(&self) -> Option<SS>
fn to_subset(&self) -> Option<SS>
The inverse inclusion map: attempts to construct
self
from the equivalent element of its
superset. Read more§fn is_in_subset(&self) -> bool
fn is_in_subset(&self) -> bool
Checks if
self
is actually part of its subset T
(and can be converted to it).§fn to_subset_unchecked(&self) -> SS
fn to_subset_unchecked(&self) -> SS
Use with care! Same as
self.to_subset
but without any property checks. Always succeeds.§fn from_subset(element: &SS) -> SP
fn from_subset(element: &SS) -> SP
The inclusion map: converts
self
to the equivalent element of its superset.