pub struct IntegratorOptions {
pub init_step: Duration,
pub min_step: Duration,
pub max_step: Duration,
pub tolerance: f64,
pub attempts: u8,
pub fixed_step: bool,
pub error_ctrl: ErrorControl,
pub integration_frame: Option<Frame>,
}Expand description
Stores the integrator options, including the minimum and maximum step sizes, and the central body to perform the integration.
Note that different step sizes and max errors are only used for adaptive
methods. To use a fixed step integrator, initialize the options using with_fixed_step, and
use whichever adaptive step integrator is desired. For example, initializing an RK45 with
fixed step options will lead to an RK4 being used instead of an RK45.
Fields§
§init_step: Duration§min_step: Duration§max_step: Duration§tolerance: f64§attempts: u8§fixed_step: bool§error_ctrl: ErrorControl§integration_frame: 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.
Implementations§
Source§impl IntegratorOptions
impl IntegratorOptions
Sourcepub fn builder() -> IntegratorOptionsBuilder<((), (), (), (), (), (), (), ())>
pub fn builder() -> IntegratorOptionsBuilder<((), (), (), (), (), (), (), ())>
Create a builder for building IntegratorOptions.
On the builder, call .init_step(...)(optional), .min_step(...)(optional), .max_step(...)(optional), .tolerance(...)(optional), .attempts(...)(optional), .fixed_step(...)(optional), .error_ctrl(...)(optional), .integration_frame(...)(optional) to set the values of the fields.
Finally, call .build() to create the instance of IntegratorOptions.
Examples found in repository?
129fn evaluate_weights(
130 weights: &[f32],
131 prop_time_days: f64,
132 state: Arc<SharedState>,
133) -> Result<(f64, f64), Box<dyn Error>> {
134 let ηthresholds: Vec<f64> = weights.iter().map(|w| *w as f64).collect();
135
136 let eme2k = state.almanac.frame_info(EARTH_J2000).unwrap();
137 let epoch = Epoch::from_gregorian_utc_hms(2024, 2, 29, 12, 13, 14);
138
139 let orbit = Orbit::keplerian(24505.9, 0.725, 7.05, 0.0, 0.0, 0.0, epoch, eme2k);
140
141 let sc = Spacecraft::builder()
142 .orbit(orbit)
143 .mass(Mass::from_dry_and_prop_masses(1000.0, 1000.0))
144 .srp(SRPData::from_area(3.0 * 6.0))
145 .thruster(Thruster {
146 isp_s: 4435.0,
147 thrust_N: 0.472,
148 })
149 .mode(GuidanceMode::Thrust)
150 .build();
151
152 let prop_time = prop_time_days * Unit::Day;
153
154 let objectives = &[
155 Objective::within_tolerance(
156 StateParameter::Element(OrbitalElement::SemiMajorAxis),
157 30_000.0,
158 20.0,
159 ),
160 Objective::within_tolerance(
161 StateParameter::Element(OrbitalElement::Eccentricity),
162 0.001,
163 5e-5,
164 ),
165 Objective::within_tolerance(
166 StateParameter::Element(OrbitalElement::Inclination),
167 0.05,
168 1e-2,
169 ),
170 ];
171
172 // let kluever_ctrl = Kluever::from_max_eclipse(objectives, &weights_f64, 0.2);
173 let ctrl = Ruggiero::from_ηthresholds(objectives, &ηthresholds, sc)?;
174
175 let mut orbital_dyn = OrbitalDynamics::point_masses(vec![MOON, SUN]);
176 orbital_dyn.accel_models.push(state.harmonics.clone());
177
178 let sc_dynamics = SpacecraftDynamics::from_model(orbital_dyn, state.srp_dyn.clone())
179 .with_guidance_law(ctrl.clone());
180
181 let (final_state, _traj) = Propagator::rk89(
182 sc_dynamics.clone(),
183 IntegratorOptions::builder()
184 .min_step(10.0_f64.seconds())
185 .tolerance(1e-8)
186 .error_ctrl(ErrorControl::RSSCartesianStep)
187 .build(),
188 )
189 .with(sc, state.almanac.clone())
190 .for_duration_with_traj(prop_time)?;
191
192 let prop_usage = sc.mass.prop_mass_kg - final_state.mass.prop_mass_kg;
193
194 let mut penalty = 0.0;
195 for obj in objectives {
196 let (achieved, error) = obj.assess(&final_state)?;
197 if !achieved {
198 penalty += error.abs();
199 }
200 info!("{obj} error: {error:.3}, achieved? {achieved}");
201 }
202
203 info!("{ηthresholds:?} -> {prop_usage:.3} kg\tpenalty = {penalty:.3}");
204
205 Ok((prop_usage, penalty * 1000.0))
206}More examples
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_info(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(
58 StateParameter::Element(OrbitalElement::SemiMajorAxis),
59 42_165.0,
60 20.0,
61 ),
62 Objective::within_tolerance(
63 StateParameter::Element(OrbitalElement::Eccentricity),
64 0.001,
65 5e-5,
66 ),
67 Objective::within_tolerance(
68 StateParameter::Element(OrbitalElement::Inclination),
69 0.05,
70 1e-2,
71 ),
72 ];
73
74 let ruggiero_ctrl = Ruggiero::from_max_eclipse(objectives, sc, 0.2)?;
75 println!("{ruggiero_ctrl}");
76
77 let mut orbital_dyn = OrbitalDynamics::point_masses(vec![MOON, SUN]);
78
79 let mut jgm3_meta = MetaFile {
80 uri: "http://public-data.nyxspace.com/nyx/models/JGM3.cof.gz".to_string(),
81 crc32: Some(0xF446F027), // Specifying the CRC32 avoids redownloading it if it's cached.
82 };
83 jgm3_meta.process(true)?;
84
85 let harmonics = GravityField::from_stor(
86 almanac.frame_info(IAU_EARTH_FRAME)?,
87 GravityFieldData::from_cof(&jgm3_meta.uri, 8, 8, true)?,
88 );
89 orbital_dyn.accel_models.push(harmonics);
90
91 let srp_dyn = SolarPressure::default_flux(EARTH_J2000, almanac.clone())?;
92 let sc_dynamics = SpacecraftDynamics::from_model(orbital_dyn, srp_dyn)
93 .with_guidance_law(ruggiero_ctrl.clone());
94
95 println!("{sc_dynamics}");
96
97 // Finally, let's use the Monte Carlo framework built into Nyx to propagate spacecraft.
98
99 // Let's start by defining the dispersion.
100 // 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.
101 // Note that additional validation on the MVN is in progress -- https://github.com/nyx-space/nyx/issues/339.
102 let mc_rv = MvnSpacecraft::new(
103 sc,
104 vec![StateDispersion::zero_mean(
105 StateParameter::Element(OrbitalElement::SemiMajorAxis),
106 3.0,
107 )],
108 )?;
109
110 let my_mc = MonteCarlo::new(
111 sc, // Nominal state
112 mc_rv,
113 "03_geo_sk".to_string(), // Scenario name
114 None, // No specific seed specified, so one will be drawn from the computer's entropy.
115 );
116
117 // Build the propagator setup.
118 let setup = Propagator::rk89(
119 sc_dynamics.clone(),
120 IntegratorOptions::builder()
121 .min_step(10.0_f64.seconds())
122 .error_ctrl(ErrorControl::RSSCartesianStep)
123 .build(),
124 );
125
126 let num_runs = 25;
127 let rslts = my_mc.run_until_epoch(setup, almanac.clone(), sc.epoch() + prop_time, num_runs);
128
129 assert_eq!(rslts.runs.len(), num_runs);
130
131 rslts.to_parquet("03_geo_sk.parquet", ExportCfg::default())?;
132
133 Ok(())
134}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_info(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(
69 StateParameter::Element(OrbitalElement::SemiMajorAxis),
70 42_165.0,
71 20.0,
72 ),
73 Objective::within_tolerance(
74 StateParameter::Element(OrbitalElement::Eccentricity),
75 0.001,
76 5e-5,
77 ),
78 Objective::within_tolerance(
79 StateParameter::Element(OrbitalElement::Inclination),
80 0.05,
81 1e-2,
82 ),
83 ];
84
85 // Ensure that we only thrust if we have more than 20% illumination.
86 let ruggiero_ctrl = Ruggiero::from_max_eclipse(objectives, sc, 0.2).unwrap();
87 println!("{ruggiero_ctrl}");
88
89 // Define the high fidelity dynamics
90
91 // Set up the spacecraft dynamics.
92
93 // Specify that the orbital dynamics must account for the graviational pull of the Moon and the Sun.
94 // The gravity of the Earth will also be accounted for since the spaceraft in an Earth orbit.
95 let mut orbital_dyn = OrbitalDynamics::point_masses(vec![MOON, SUN]);
96
97 // We want to include the spherical harmonics, so let's download the gravitational data from the Nyx Cloud.
98 // We're using the JGM3 model here, which is the default in GMAT.
99 let mut jgm3_meta = MetaFile {
100 uri: "http://public-data.nyxspace.com/nyx/models/JGM3.cof.gz".to_string(),
101 crc32: Some(0xF446F027), // Specifying the CRC32 avoids redownloading it if it's cached.
102 };
103 // And let's download it if we don't have it yet.
104 jgm3_meta.process(true)?;
105
106 // Build the spherical harmonics.
107 // The harmonics must be computed in the body fixed frame.
108 // We're using the long term prediction of the Earth centered Earth fixed frame, IAU Earth.
109 let harmonics = GravityField::from_stor(
110 almanac.frame_info(IAU_EARTH_FRAME)?,
111 GravityFieldData::from_cof(&jgm3_meta.uri, 8, 8, true).unwrap(),
112 );
113
114 // Include the spherical harmonics into the orbital dynamics.
115 orbital_dyn.accel_models.push(harmonics);
116
117 // We define the solar radiation pressure, using the default solar flux and accounting only
118 // for the eclipsing caused by the Earth.
119 let srp_dyn = SolarPressure::default_flux(EARTH_J2000, almanac.clone())?;
120
121 // Finalize setting up the dynamics, specifying the force models (orbital_dyn) separately from the
122 // acceleration models (SRP in this case). Use `from_models` to specify multiple accel models.
123 let sc_dynamics = SpacecraftDynamics::from_model(orbital_dyn, srp_dyn)
124 .with_guidance_law(ruggiero_ctrl.clone());
125
126 println!("{orbit:x}");
127
128 // We specify a minimum step in the propagator because the Ruggiero control would otherwise drive this step very low.
129 let (final_state, traj) = Propagator::rk89(
130 sc_dynamics.clone(),
131 IntegratorOptions::builder()
132 .min_step(10.0_f64.seconds())
133 .error_ctrl(ErrorControl::RSSCartesianStep)
134 .build(),
135 )
136 .with(sc, almanac.clone())
137 .for_duration_with_traj(prop_time)?;
138
139 let prop_usage = sc.mass.prop_mass_kg - final_state.mass.prop_mass_kg;
140 println!("{:x}", final_state.orbit);
141 println!("prop usage: {prop_usage:.3} kg");
142
143 // Finally, export the results for analysis, including the penumbra percentage throughout the orbit raise.
144 traj.to_parquet("./03_geo_raise.parquet", ExportCfg::default())?;
145
146 for status_line in ruggiero_ctrl.status(&final_state) {
147 println!("{status_line}");
148 }
149
150 ruggiero_ctrl
151 .achieved(&final_state)
152 .expect("objective not achieved");
153
154 Ok(())
155}34fn main() -> Result<(), Box<dyn Error>> {
35 pel::init();
36
37 // ====================== //
38 // === ALMANAC SET UP === //
39 // ====================== //
40
41 let manifest_dir = PathBuf::from(env!("CARGO_MANIFEST_DIR"));
42
43 let out = manifest_dir.join("data/04_output/");
44
45 let almanac = Arc::new(
46 Almanac::new(
47 &manifest_dir
48 .join("data/01_planetary/pck08.pca")
49 .to_string_lossy(),
50 )
51 .unwrap()
52 .load(
53 &manifest_dir
54 .join("data/01_planetary/de440s.bsp")
55 .to_string_lossy(),
56 )
57 .unwrap(),
58 );
59
60 let eme2k = almanac.frame_info(EARTH_J2000).unwrap();
61 let moon_iau = almanac.frame_info(IAU_MOON_FRAME).unwrap();
62
63 let epoch = Epoch::from_gregorian_tai(2021, 5, 29, 19, 51, 16, 852_000);
64 let nrho = Orbit::cartesian(
65 166_473.631_302_239_7,
66 -274_715.487_253_382_7,
67 -211_233.210_176_686_7,
68 0.933_451_604_520_018_4,
69 0.436_775_046_841_900_9,
70 -0.082_211_021_250_348_95,
71 epoch,
72 eme2k,
73 );
74
75 let tx_nrho_sc = Spacecraft::from(nrho);
76
77 let state_luna = almanac.transform_to(nrho, MOON_J2000, None).unwrap();
78 println!("Start state (dynamics: Earth, Moon, Sun gravity):\n{state_luna}");
79
80 let bodies = vec![EARTH, SUN];
81 let dynamics = SpacecraftDynamics::new(OrbitalDynamics::point_masses(bodies));
82
83 let setup = Propagator::rk89(
84 dynamics,
85 IntegratorOptions::builder().max_step(0.5.minutes()).build(),
86 );
87
88 /* == Propagate the NRHO vehicle == */
89 let prop_time = 1.1 * state_luna.period().unwrap();
90
91 let (nrho_final, mut tx_traj) = setup
92 .with(tx_nrho_sc, almanac.clone())
93 .for_duration_with_traj(prop_time)
94 .unwrap();
95
96 tx_traj.name = Some("NRHO Tx SC".to_string());
97
98 println!("{tx_traj}");
99
100 /* == Propagate an LLO vehicle == */
101 let llo_orbit =
102 Orbit::try_keplerian_altitude(110.0, 1e-4, 90.0, 0.0, 0.0, 0.0, epoch, moon_iau).unwrap();
103
104 let llo_sc = Spacecraft::builder().orbit(llo_orbit).build();
105
106 let (_, llo_traj) = setup
107 .with(llo_sc, almanac.clone())
108 .until_epoch_with_traj(nrho_final.epoch())
109 .unwrap();
110
111 // Export the subset of the first two hours.
112 llo_traj
113 .clone()
114 .filter_by_offset(..2.hours())
115 .to_parquet_simple(out.join("05_caps_llo_truth.pq"))?;
116
117 /* == Setup the interlink == */
118
119 let mut measurement_types = IndexSet::new();
120 measurement_types.insert(MeasurementType::Range);
121 measurement_types.insert(MeasurementType::Doppler);
122
123 let mut stochastics = IndexMap::new();
124
125 let sa45_csac_allan_dev = 1e-11;
126
127 stochastics.insert(
128 MeasurementType::Range,
129 StochasticNoise::from_hardware_range_km(
130 sa45_csac_allan_dev,
131 10.0.seconds(),
132 link_specific::ChipRate::StandardT4B,
133 link_specific::SN0::Average,
134 ),
135 );
136
137 stochastics.insert(
138 MeasurementType::Doppler,
139 StochasticNoise::from_hardware_doppler_km_s(
140 sa45_csac_allan_dev,
141 10.0.seconds(),
142 link_specific::CarrierFreq::SBand,
143 link_specific::CN0::Average,
144 ),
145 );
146
147 let interlink = InterlinkTxSpacecraft {
148 traj: tx_traj,
149 measurement_types,
150 integration_time: None,
151 timestamp_noise_s: None,
152 ab_corr: Aberration::LT,
153 stochastic_noises: Some(stochastics),
154 };
155
156 // Devices are the transmitter, which is our NRHO vehicle.
157 let mut devices = BTreeMap::new();
158 devices.insert("NRHO Tx SC".to_string(), interlink);
159
160 let mut configs = BTreeMap::new();
161 configs.insert(
162 "NRHO Tx SC".to_string(),
163 TrkConfig::builder()
164 .strands(vec![Strand {
165 start: epoch,
166 end: nrho_final.epoch(),
167 }])
168 .build(),
169 );
170
171 let mut trk_sim =
172 TrackingArcSim::with_seed(devices.clone(), llo_traj.clone(), configs, 0).unwrap();
173 println!("{trk_sim}");
174
175 let trk_data = trk_sim.generate_measurements(almanac.clone()).unwrap();
176 println!("{trk_data}");
177
178 trk_data
179 .to_parquet_simple(out.clone().join("nrho_interlink_msr.pq"))
180 .unwrap();
181
182 // Run a truth OD where we estimate the LLO position
183 let llo_uncertainty = SpacecraftUncertainty::builder()
184 .nominal(llo_sc)
185 .x_km(1.0)
186 .y_km(1.0)
187 .z_km(1.0)
188 .vx_km_s(1e-3)
189 .vy_km_s(1e-3)
190 .vz_km_s(1e-3)
191 .build();
192
193 let mut proc_devices = devices.clone();
194
195 // Define the initial estimate, randomized, seed for reproducibility
196 let mut initial_estimate = llo_uncertainty.to_estimate_randomized(Some(0)).unwrap();
197 // Inflate the covariance -- https://github.com/nyx-space/nyx/issues/339
198 initial_estimate.covar *= 2.5;
199
200 // Increase the noise in the devices to accept more measurements.
201
202 for link in proc_devices.values_mut() {
203 for noise in &mut link.stochastic_noises.as_mut().unwrap().values_mut() {
204 *noise.white_noise.as_mut().unwrap() *= 3.0;
205 }
206 }
207
208 let init_err = initial_estimate
209 .orbital_state()
210 .ric_difference(&llo_orbit)
211 .unwrap();
212
213 println!("initial estimate:\n{initial_estimate}");
214 println!("RIC errors = {init_err}",);
215
216 let odp = InterlinkKalmanOD::new(
217 setup.clone(),
218 KalmanVariant::ReferenceUpdate,
219 Some(ResidRejectCrit::default()),
220 proc_devices,
221 almanac.clone(),
222 );
223
224 // Shrink the data to process.
225 let arc = trk_data.filter_by_offset(..2.hours());
226
227 let od_sol = odp.process_arc(initial_estimate, &arc).unwrap();
228
229 println!("{od_sol}");
230
231 od_sol
232 .to_parquet(
233 out.join("05_caps_interlink_od_sol.pq"),
234 ExportCfg::default(),
235 )
236 .unwrap();
237
238 let od_traj = od_sol.to_traj().unwrap();
239
240 od_traj
241 .ric_diff_to_parquet(
242 &llo_traj,
243 out.join("05_caps_interlink_llo_est_error.pq"),
244 ExportCfg::default(),
245 )
246 .unwrap();
247
248 let final_est = od_sol.estimates.last().unwrap();
249 assert!(final_est.within_3sigma(), "should be within 3 sigma");
250
251 println!("ESTIMATE\n{final_est:x}\n");
252 let truth = llo_traj.at(final_est.epoch()).unwrap();
253 println!("TRUTH\n{truth:x}");
254
255 let final_err = truth
256 .orbit
257 .ric_difference(&final_est.orbital_state())
258 .unwrap();
259 println!("ERROR {final_err}");
260
261 // Build the residuals versus reference plot.
262 let rvr_sol = odp
263 .process_arc(initial_estimate, &arc.resid_vs_ref_check())
264 .unwrap();
265
266 rvr_sol
267 .to_parquet(
268 out.join("05_caps_interlink_resid_v_ref.pq"),
269 ExportCfg::default(),
270 )
271 .unwrap();
272
273 let final_rvr = rvr_sol.estimates.last().unwrap();
274
275 println!("RMAG error {:.3} m", final_err.rmag_km() * 1e3);
276 println!(
277 "Pure prop error {:.3} m",
278 final_rvr
279 .orbital_state()
280 .ric_difference(&final_est.orbital_state())
281 .unwrap()
282 .rmag_km()
283 * 1e3
284 );
285
286 Ok(())
287}Source§impl IntegratorOptions
impl IntegratorOptions
Sourcepub fn with_adaptive_step(
min_step: Duration,
max_step: Duration,
tolerance: f64,
error_ctrl: ErrorControl,
) -> Self
pub fn with_adaptive_step( min_step: Duration, max_step: Duration, tolerance: f64, error_ctrl: ErrorControl, ) -> Self
with_adaptive_step initializes an PropOpts such that the integrator is used with an
adaptive step size. The number of attempts is currently fixed to 50 (as in GMAT).
pub fn with_adaptive_step_s( min_step: f64, max_step: f64, tolerance: f64, error_ctrl: ErrorControl, ) -> Self
Sourcepub fn with_fixed_step(step: Duration) -> Self
pub fn with_fixed_step(step: Duration) -> Self
with_fixed_step initializes an PropOpts such that the integrator is used with a fixed
step size.
pub fn with_fixed_step_s(step: f64) -> Self
Sourcepub fn with_tolerance(tolerance: f64) -> Self
pub fn with_tolerance(tolerance: f64) -> Self
Returns the default options with a specific tolerance.
Sourcepub fn with_max_step(max_step: Duration) -> Self
pub fn with_max_step(max_step: Duration) -> Self
Creates a propagator with the provided max step, and sets the initial step to that value as well.
Sourcepub fn set_max_step(&mut self, max_step: Duration)
pub fn set_max_step(&mut self, max_step: Duration)
Set the maximum step size and sets the initial step to that value if currently greater
Sourcepub fn set_min_step(&mut self, min_step: Duration)
pub fn set_min_step(&mut self, min_step: Duration)
Set the minimum step size and sets the initial step to that value if currently smaller
Trait Implementations§
Source§impl Clone for IntegratorOptions
impl Clone for IntegratorOptions
Source§fn clone(&self) -> IntegratorOptions
fn clone(&self) -> IntegratorOptions
1.0.0 (const: unstable) · Source§fn clone_from(&mut self, source: &Self)
fn clone_from(&mut self, source: &Self)
source. Read moreSource§impl Debug for IntegratorOptions
impl Debug for IntegratorOptions
Source§impl Default for IntegratorOptions
impl Default for IntegratorOptions
Source§fn default() -> IntegratorOptions
fn default() -> IntegratorOptions
default returns the same default options as GMAT.
Source§impl<'de> Deserialize<'de> for IntegratorOptions
impl<'de> Deserialize<'de> for IntegratorOptions
Source§fn deserialize<__D>(__deserializer: __D) -> Result<Self, __D::Error>where
__D: Deserializer<'de>,
fn deserialize<__D>(__deserializer: __D) -> Result<Self, __D::Error>where
__D: Deserializer<'de>,
Source§impl Display for IntegratorOptions
impl Display for IntegratorOptions
Source§impl PartialEq for IntegratorOptions
impl PartialEq for IntegratorOptions
Source§fn eq(&self, other: &IntegratorOptions) -> bool
fn eq(&self, other: &IntegratorOptions) -> bool
self and other values to be equal, and is used by ==.Source§impl Serialize for IntegratorOptions
impl Serialize for IntegratorOptions
Source§impl StaticType for IntegratorOptions
impl StaticType for IntegratorOptions
Source§fn static_type() -> SimpleType
fn static_type() -> SimpleType
impl Copy for IntegratorOptions
impl StructuralPartialEq for IntegratorOptions
Auto Trait Implementations§
impl Freeze for IntegratorOptions
impl RefUnwindSafe for IntegratorOptions
impl Send for IntegratorOptions
impl Sync for IntegratorOptions
impl Unpin for IntegratorOptions
impl UnsafeUnpin for IntegratorOptions
impl UnwindSafe for IntegratorOptions
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
Source§impl<T> CloneToUninit for Twhere
T: Clone,
impl<T> CloneToUninit for Twhere
T: Clone,
Source§impl<T> FromDhall for Twhere
T: DeserializeOwned,
impl<T> FromDhall for Twhere
T: DeserializeOwned,
fn from_dhall(v: &Value) -> Result<T, Error>
§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>
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>
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>
self from the equivalent element of its
superset. Read more§fn is_in_subset(&self) -> bool
fn is_in_subset(&self) -> bool
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
self.to_subset but without any property checks. Always succeeds.§fn from_subset(element: &SS) -> SP
fn from_subset(element: &SS) -> SP
self to the equivalent element of its superset.