Struct IntegratorOptionsBuilder 

pub struct IntegratorOptionsBuilder<TypedBuilderFields = ((), (), (), (), (), (), (), ())> { /* private fields */ }

Builder for IntegratorOptions instances.

See IntegratorOptions::builder() for more info.

Implementations

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)>

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)>

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

nyx-core/examples/03_geo_analysis/raise_optim.rs (line 184)

fn evaluate_weights(
    weights: &[f32],
    prop_time_days: f64,
    state: Arc<SharedState>,
) -> Result<(f64, f64), Box<dyn Error>> {
    let ηthresholds: Vec<f64> = weights.iter().map(|w| *w as f64).collect();

    let eme2k = state.almanac.frame_info(EARTH_J2000).unwrap();
    let epoch = Epoch::from_gregorian_utc_hms(2024, 2, 29, 12, 13, 14);

    let orbit = Orbit::keplerian(24505.9, 0.725, 7.05, 0.0, 0.0, 0.0, epoch, eme2k);

    let sc = Spacecraft::builder()
        .orbit(orbit)
        .mass(Mass::from_dry_and_prop_masses(1000.0, 1000.0))
        .srp(SRPData::from_area(3.0 * 6.0))
        .thruster(Thruster {
            isp_s: 4435.0,
            thrust_N: 0.472,
        })
        .mode(GuidanceMode::Thrust)
        .build();

    let prop_time = prop_time_days * Unit::Day;

    let objectives = &[
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::SemiMajorAxis),
            30_000.0,
            20.0,
        ),
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::Eccentricity),
            0.001,
            5e-5,
        ),
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::Inclination),
            0.05,
            1e-2,
        ),
    ];

    // let kluever_ctrl = Kluever::from_max_eclipse(objectives, &weights_f64, 0.2);
    let ctrl = Ruggiero::from_ηthresholds(objectives, &ηthresholds, sc)?;

    let mut orbital_dyn = OrbitalDynamics::point_masses(vec![MOON, SUN]);
    orbital_dyn.accel_models.push(state.harmonics.clone());

    let sc_dynamics = SpacecraftDynamics::from_model(orbital_dyn, state.srp_dyn.clone())
        .with_guidance_law(ctrl.clone());

    let (final_state, _traj) = Propagator::rk89(
        sc_dynamics.clone(),
        IntegratorOptions::builder()
            .min_step(10.0_f64.seconds())
            .tolerance(1e-8)
            .error_ctrl(ErrorControl::RSSCartesianStep)
            .build(),
    )
    .with(sc, state.almanac.clone())
    .for_duration_with_traj(prop_time)?;

    let prop_usage = sc.mass.prop_mass_kg - final_state.mass.prop_mass_kg;

    let mut penalty = 0.0;
    for obj in objectives {
        let (achieved, error) = obj.assess(&final_state)?;
        if !achieved {
            penalty += error.abs();
        }
        info!("{obj} error: {error:.3}, achieved? {achieved}");
    }

    info!("{ηthresholds:?} -> {prop_usage:.3} kg\tpenalty = {penalty:.3}");

    Ok((prop_usage, penalty * 1000.0))
}

nyx-core/examples/03_geo_analysis/stationkeeping.rs (line 121)

fn main() -> Result<(), Box<dyn Error>> {
    pel::init();
    // Set up the dynamics like in the orbit raise.
    let almanac = Arc::new(MetaAlmanac::latest().map_err(Box::new)?);
    let epoch = Epoch::from_gregorian_utc_hms(2024, 2, 29, 12, 13, 14);

    // Define the GEO orbit, and we're just going to maintain it very tightly.
    let earth_j2000 = almanac.frame_info(EARTH_J2000)?;
    let orbit = Orbit::try_keplerian(42164.0, 1e-5, 0., 163.0, 75.0, 0.0, epoch, earth_j2000)?;
    println!("{orbit:x}");

    let sc = Spacecraft::builder()
        .orbit(orbit)
        .mass(Mass::from_dry_and_prop_masses(1000.0, 1000.0)) // 1000 kg of dry mass and prop, totalling 2.0 tons
        .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
        .thruster(Thruster {
            // "NEXT-STEP" row in Table 2
            isp_s: 4435.0,
            thrust_N: 0.472,
        })
        .mode(GuidanceMode::Thrust) // Start thrusting immediately.
        .build();

    // Set up the spacecraft dynamics like in the orbit raise example.

    let prop_time = 30.0 * Unit::Day;

    // Define the guidance law -- we're just using a Ruggiero controller as demonstrated in AAS-2004-5089.
    let objectives = &[
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::SemiMajorAxis),
            42_165.0,
            20.0,
        ),
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::Eccentricity),
            0.001,
            5e-5,
        ),
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::Inclination),
            0.05,
            1e-2,
        ),
    ];

    let ruggiero_ctrl = Ruggiero::from_max_eclipse(objectives, sc, 0.2)?;
    println!("{ruggiero_ctrl}");

    let mut orbital_dyn = OrbitalDynamics::point_masses(vec![MOON, SUN]);

    let mut jgm3_meta = MetaFile {
        uri: "http://public-data.nyxspace.com/nyx/models/JGM3.cof.gz".to_string(),
        crc32: Some(0xF446F027), // Specifying the CRC32 avoids redownloading it if it's cached.
    };
    jgm3_meta.process(true)?;

    let harmonics = GravityField::from_stor(
        almanac.frame_info(IAU_EARTH_FRAME)?,
        GravityFieldData::from_cof(&jgm3_meta.uri, 8, 8, true)?,
    );
    orbital_dyn.accel_models.push(harmonics);

    let srp_dyn = SolarPressure::default_flux(EARTH_J2000, almanac.clone())?;
    let sc_dynamics = SpacecraftDynamics::from_model(orbital_dyn, srp_dyn)
        .with_guidance_law(ruggiero_ctrl.clone());

    println!("{sc_dynamics}");

    // Finally, let's use the Monte Carlo framework built into Nyx to propagate spacecraft.

    // Let's start by defining the dispersion.
    // 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.
    // Note that additional validation on the MVN is in progress -- https://github.com/nyx-space/nyx/issues/339.
    let mc_rv = MvnSpacecraft::new(
        sc,
        vec![StateDispersion::zero_mean(
            StateParameter::Element(OrbitalElement::SemiMajorAxis),
            3.0,
        )],
    )?;

    let my_mc = MonteCarlo::new(
        sc, // Nominal state
        mc_rv,
        "03_geo_sk".to_string(), // Scenario name
        None, // No specific seed specified, so one will be drawn from the computer's entropy.
    );

    // Build the propagator setup.
    let setup = Propagator::rk89(
        sc_dynamics.clone(),
        IntegratorOptions::builder()
            .min_step(10.0_f64.seconds())
            .error_ctrl(ErrorControl::RSSCartesianStep)
            .build(),
    );

    let num_runs = 25;
    let rslts = my_mc.run_until_epoch(setup, almanac.clone(), sc.epoch() + prop_time, num_runs);

    assert_eq!(rslts.runs.len(), num_runs);

    rslts.to_parquet("03_geo_sk.parquet", ExportCfg::default())?;

    Ok(())
}

nyx-core/examples/03_geo_analysis/raise.rs (line 132)

fn main() -> Result<(), Box<dyn Error>> {
    pel::init();

    // Dynamics models require planetary constants and ephemerides to be defined.
    // Let's start by grabbing those by using ANISE's latest MetaAlmanac.
    // This will automatically download the DE440s planetary ephemeris,
    // the daily-updated Earth Orientation Parameters, the high fidelity Moon orientation
    // parameters (for the Moon Mean Earth and Moon Principal Axes frames), and the PCK11
    // planetary constants kernels.
    // For details, refer to https://github.com/nyx-space/anise/blob/master/data/latest.dhall.
    // Note that we place the Almanac into an Arc so we can clone it cheaply and provide read-only
    // references to many functions.
    let almanac = Arc::new(MetaAlmanac::latest().map_err(Box::new)?);
    // Fetch the EME2000 frame from the Almabac
    let eme2k = almanac.frame_info(EARTH_J2000).unwrap();
    // Define the orbit epoch
    let epoch = Epoch::from_gregorian_utc_hms(2024, 2, 29, 12, 13, 14);

    // Build the spacecraft itself.
    // Using slide 6 of https://aerospace.org/sites/default/files/2018-11/Davis-Mayberry_HPSEP_11212018.pdf
    // for the "next gen" SEP characteristics.

    // GTO start
    let orbit = Orbit::keplerian(24505.9, 0.725, 7.05, 0.0, 0.0, 0.0, epoch, eme2k);

    let sc = Spacecraft::builder()
        .orbit(orbit)
        .mass(Mass::from_dry_and_prop_masses(1000.0, 1000.0)) // 1000 kg of dry mass and prop, totalling 2.0 tons
        .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
        .thruster(Thruster {
            // "NEXT-STEP" row in Table 2
            isp_s: 4435.0,
            thrust_N: 0.472,
        })
        .mode(GuidanceMode::Thrust) // Start thrusting immediately.
        .build();

    let prop_time = 180.0 * Unit::Day;

    // Define the guidance law -- we're just using a Ruggiero controller as demonstrated in AAS-2004-5089.
    let objectives = &[
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::SemiMajorAxis),
            42_165.0,
            20.0,
        ),
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::Eccentricity),
            0.001,
            5e-5,
        ),
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::Inclination),
            0.05,
            1e-2,
        ),
    ];

    // Ensure that we only thrust if we have more than 20% illumination.
    let ruggiero_ctrl = Ruggiero::from_max_eclipse(objectives, sc, 0.2).unwrap();
    println!("{ruggiero_ctrl}");

    // Define the high fidelity dynamics

    // Set up the spacecraft dynamics.

    // Specify that the orbital dynamics must account for the graviational pull of the Moon and the Sun.
    // The gravity of the Earth will also be accounted for since the spaceraft in an Earth orbit.
    let mut orbital_dyn = OrbitalDynamics::point_masses(vec![MOON, SUN]);

    // We want to include the spherical harmonics, so let's download the gravitational data from the Nyx Cloud.
    // We're using the JGM3 model here, which is the default in GMAT.
    let mut jgm3_meta = MetaFile {
        uri: "http://public-data.nyxspace.com/nyx/models/JGM3.cof.gz".to_string(),
        crc32: Some(0xF446F027), // Specifying the CRC32 avoids redownloading it if it's cached.
    };
    // And let's download it if we don't have it yet.
    jgm3_meta.process(true)?;

    // Build the spherical harmonics.
    // The harmonics must be computed in the body fixed frame.
    // We're using the long term prediction of the Earth centered Earth fixed frame, IAU Earth.
    let harmonics = GravityField::from_stor(
        almanac.frame_info(IAU_EARTH_FRAME)?,
        GravityFieldData::from_cof(&jgm3_meta.uri, 8, 8, true).unwrap(),
    );

    // Include the spherical harmonics into the orbital dynamics.
    orbital_dyn.accel_models.push(harmonics);

    // We define the solar radiation pressure, using the default solar flux and accounting only
    // for the eclipsing caused by the Earth.
    let srp_dyn = SolarPressure::default_flux(EARTH_J2000, almanac.clone())?;

    // Finalize setting up the dynamics, specifying the force models (orbital_dyn) separately from the
    // acceleration models (SRP in this case). Use `from_models` to specify multiple accel models.
    let sc_dynamics = SpacecraftDynamics::from_model(orbital_dyn, srp_dyn)
        .with_guidance_law(ruggiero_ctrl.clone());

    println!("{orbit:x}");

    // We specify a minimum step in the propagator because the Ruggiero control would otherwise drive this step very low.
    let (final_state, traj) = Propagator::rk89(
        sc_dynamics.clone(),
        IntegratorOptions::builder()
            .min_step(10.0_f64.seconds())
            .error_ctrl(ErrorControl::RSSCartesianStep)
            .build(),
    )
    .with(sc, almanac.clone())
    .for_duration_with_traj(prop_time)?;

    let prop_usage = sc.mass.prop_mass_kg - final_state.mass.prop_mass_kg;
    println!("{:x}", final_state.orbit);
    println!("prop usage: {prop_usage:.3} kg");

    // Finally, export the results for analysis, including the penumbra percentage throughout the orbit raise.
    traj.to_parquet("./03_geo_raise.parquet", ExportCfg::default())?;

    for status_line in ruggiero_ctrl.status(&final_state) {
        println!("{status_line}");
    }

    ruggiero_ctrl
        .achieved(&final_state)
        .expect("objective not achieved");

    Ok(())
}

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)>

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

nyx-core/examples/05_cislunar_spacecraft_link_od/main.rs (line 85)

fn main() -> Result<(), Box<dyn Error>> {
    pel::init();

    // ====================== //
    // === ALMANAC SET UP === //
    // ====================== //

    let manifest_dir = PathBuf::from(env!("CARGO_MANIFEST_DIR"));

    let out = manifest_dir.join("data/04_output/");

    let almanac = Arc::new(
        Almanac::new(
            &manifest_dir
                .join("data/01_planetary/pck08.pca")
                .to_string_lossy(),
        )
        .unwrap()
        .load(
            &manifest_dir
                .join("data/01_planetary/de440s.bsp")
                .to_string_lossy(),
        )
        .unwrap(),
    );

    let eme2k = almanac.frame_info(EARTH_J2000).unwrap();
    let moon_iau = almanac.frame_info(IAU_MOON_FRAME).unwrap();

    let epoch = Epoch::from_gregorian_tai(2021, 5, 29, 19, 51, 16, 852_000);
    let nrho = Orbit::cartesian(
        166_473.631_302_239_7,
        -274_715.487_253_382_7,
        -211_233.210_176_686_7,
        0.933_451_604_520_018_4,
        0.436_775_046_841_900_9,
        -0.082_211_021_250_348_95,
        epoch,
        eme2k,
    );

    let tx_nrho_sc = Spacecraft::from(nrho);

    let state_luna = almanac.transform_to(nrho, MOON_J2000, None).unwrap();
    println!("Start state (dynamics: Earth, Moon, Sun gravity):\n{state_luna}");

    let bodies = vec![EARTH, SUN];
    let dynamics = SpacecraftDynamics::new(OrbitalDynamics::point_masses(bodies));

    let setup = Propagator::rk89(
        dynamics,
        IntegratorOptions::builder().max_step(0.5.minutes()).build(),
    );

    /* == Propagate the NRHO vehicle == */
    let prop_time = 1.1 * state_luna.period().unwrap();

    let (nrho_final, mut tx_traj) = setup
        .with(tx_nrho_sc, almanac.clone())
        .for_duration_with_traj(prop_time)
        .unwrap();

    tx_traj.name = Some("NRHO Tx SC".to_string());

    println!("{tx_traj}");

    /* == Propagate an LLO vehicle == */
    let llo_orbit =
        Orbit::try_keplerian_altitude(110.0, 1e-4, 90.0, 0.0, 0.0, 0.0, epoch, moon_iau).unwrap();

    let llo_sc = Spacecraft::builder().orbit(llo_orbit).build();

    let (_, llo_traj) = setup
        .with(llo_sc, almanac.clone())
        .until_epoch_with_traj(nrho_final.epoch())
        .unwrap();

    // Export the subset of the first two hours.
    llo_traj
        .clone()
        .filter_by_offset(..2.hours())
        .to_parquet_simple(out.join("05_caps_llo_truth.pq"))?;

    /* == Setup the interlink == */

    let mut measurement_types = IndexSet::new();
    measurement_types.insert(MeasurementType::Range);
    measurement_types.insert(MeasurementType::Doppler);

    let mut stochastics = IndexMap::new();

    let sa45_csac_allan_dev = 1e-11;

    stochastics.insert(
        MeasurementType::Range,
        StochasticNoise::from_hardware_range_km(
            sa45_csac_allan_dev,
            10.0.seconds(),
            link_specific::ChipRate::StandardT4B,
            link_specific::SN0::Average,
        ),
    );

    stochastics.insert(
        MeasurementType::Doppler,
        StochasticNoise::from_hardware_doppler_km_s(
            sa45_csac_allan_dev,
            10.0.seconds(),
            link_specific::CarrierFreq::SBand,
            link_specific::CN0::Average,
        ),
    );

    let interlink = InterlinkTxSpacecraft {
        traj: tx_traj,
        measurement_types,
        integration_time: None,
        timestamp_noise_s: None,
        ab_corr: Aberration::LT,
        stochastic_noises: Some(stochastics),
    };

    // Devices are the transmitter, which is our NRHO vehicle.
    let mut devices = BTreeMap::new();
    devices.insert("NRHO Tx SC".to_string(), interlink);

    let mut configs = BTreeMap::new();
    configs.insert(
        "NRHO Tx SC".to_string(),
        TrkConfig::builder()
            .strands(vec![Strand {
                start: epoch,
                end: nrho_final.epoch(),
            }])
            .build(),
    );

    let mut trk_sim =
        TrackingArcSim::with_seed(devices.clone(), llo_traj.clone(), configs, 0).unwrap();
    println!("{trk_sim}");

    let trk_data = trk_sim.generate_measurements(almanac.clone()).unwrap();
    println!("{trk_data}");

    trk_data
        .to_parquet_simple(out.clone().join("nrho_interlink_msr.pq"))
        .unwrap();

    // Run a truth OD where we estimate the LLO position
    let llo_uncertainty = SpacecraftUncertainty::builder()
        .nominal(llo_sc)
        .x_km(1.0)
        .y_km(1.0)
        .z_km(1.0)
        .vx_km_s(1e-3)
        .vy_km_s(1e-3)
        .vz_km_s(1e-3)
        .build();

    let mut proc_devices = devices.clone();

    // Define the initial estimate, randomized, seed for reproducibility
    let mut initial_estimate = llo_uncertainty.to_estimate_randomized(Some(0)).unwrap();
    // Inflate the covariance -- https://github.com/nyx-space/nyx/issues/339
    initial_estimate.covar *= 2.5;

    // Increase the noise in the devices to accept more measurements.

    for link in proc_devices.values_mut() {
        for noise in &mut link.stochastic_noises.as_mut().unwrap().values_mut() {
            *noise.white_noise.as_mut().unwrap() *= 3.0;
        }
    }

    let init_err = initial_estimate
        .orbital_state()
        .ric_difference(&llo_orbit)
        .unwrap();

    println!("initial estimate:\n{initial_estimate}");
    println!("RIC errors = {init_err}",);

    let odp = InterlinkKalmanOD::new(
        setup.clone(),
        KalmanVariant::ReferenceUpdate,
        Some(ResidRejectCrit::default()),
        proc_devices,
        almanac.clone(),
    );

    // Shrink the data to process.
    let arc = trk_data.filter_by_offset(..2.hours());

    let od_sol = odp.process_arc(initial_estimate, &arc).unwrap();

    println!("{od_sol}");

    od_sol
        .to_parquet(
            out.join("05_caps_interlink_od_sol.pq"),
            ExportCfg::default(),
        )
        .unwrap();

    let od_traj = od_sol.to_traj().unwrap();

    od_traj
        .ric_diff_to_parquet(
            &llo_traj,
            out.join("05_caps_interlink_llo_est_error.pq"),
            ExportCfg::default(),
        )
        .unwrap();

    let final_est = od_sol.estimates.last().unwrap();
    assert!(final_est.within_3sigma(), "should be within 3 sigma");

    println!("ESTIMATE\n{final_est:x}\n");
    let truth = llo_traj.at(final_est.epoch()).unwrap();
    println!("TRUTH\n{truth:x}");

    let final_err = truth
        .orbit
        .ric_difference(&final_est.orbital_state())
        .unwrap();
    println!("ERROR {final_err}");

    // Build the residuals versus reference plot.
    let rvr_sol = odp
        .process_arc(initial_estimate, &arc.resid_vs_ref_check())
        .unwrap();

    rvr_sol
        .to_parquet(
            out.join("05_caps_interlink_resid_v_ref.pq"),
            ExportCfg::default(),
        )
        .unwrap();

    let final_rvr = rvr_sol.estimates.last().unwrap();

    println!("RMAG error {:.3} m", final_err.rmag_km() * 1e3);
    println!(
        "Pure prop error {:.3} m",
        final_rvr
            .orbital_state()
            .ric_difference(&final_est.orbital_state())
            .unwrap()
            .rmag_km()
            * 1e3
    );

    Ok(())
}

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)>

Examples found in repository

nyx-core/examples/03_geo_analysis/raise_optim.rs (line 185)

fn evaluate_weights(
    weights: &[f32],
    prop_time_days: f64,
    state: Arc<SharedState>,
) -> Result<(f64, f64), Box<dyn Error>> {
    let ηthresholds: Vec<f64> = weights.iter().map(|w| *w as f64).collect();

    let eme2k = state.almanac.frame_info(EARTH_J2000).unwrap();
    let epoch = Epoch::from_gregorian_utc_hms(2024, 2, 29, 12, 13, 14);

    let orbit = Orbit::keplerian(24505.9, 0.725, 7.05, 0.0, 0.0, 0.0, epoch, eme2k);

    let sc = Spacecraft::builder()
        .orbit(orbit)
        .mass(Mass::from_dry_and_prop_masses(1000.0, 1000.0))
        .srp(SRPData::from_area(3.0 * 6.0))
        .thruster(Thruster {
            isp_s: 4435.0,
            thrust_N: 0.472,
        })
        .mode(GuidanceMode::Thrust)
        .build();

    let prop_time = prop_time_days * Unit::Day;

    let objectives = &[
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::SemiMajorAxis),
            30_000.0,
            20.0,
        ),
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::Eccentricity),
            0.001,
            5e-5,
        ),
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::Inclination),
            0.05,
            1e-2,
        ),
    ];

    // let kluever_ctrl = Kluever::from_max_eclipse(objectives, &weights_f64, 0.2);
    let ctrl = Ruggiero::from_ηthresholds(objectives, &ηthresholds, sc)?;

    let mut orbital_dyn = OrbitalDynamics::point_masses(vec![MOON, SUN]);
    orbital_dyn.accel_models.push(state.harmonics.clone());

    let sc_dynamics = SpacecraftDynamics::from_model(orbital_dyn, state.srp_dyn.clone())
        .with_guidance_law(ctrl.clone());

    let (final_state, _traj) = Propagator::rk89(
        sc_dynamics.clone(),
        IntegratorOptions::builder()
            .min_step(10.0_f64.seconds())
            .tolerance(1e-8)
            .error_ctrl(ErrorControl::RSSCartesianStep)
            .build(),
    )
    .with(sc, state.almanac.clone())
    .for_duration_with_traj(prop_time)?;

    let prop_usage = sc.mass.prop_mass_kg - final_state.mass.prop_mass_kg;

    let mut penalty = 0.0;
    for obj in objectives {
        let (achieved, error) = obj.assess(&final_state)?;
        if !achieved {
            penalty += error.abs();
        }
        info!("{obj} error: {error:.3}, achieved? {achieved}");
    }

    info!("{ηthresholds:?} -> {prop_usage:.3} kg\tpenalty = {penalty:.3}");

    Ok((prop_usage, penalty * 1000.0))
}

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)>

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)>

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)>

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

nyx-core/examples/03_geo_analysis/raise_optim.rs (line 186)

fn evaluate_weights(
    weights: &[f32],
    prop_time_days: f64,
    state: Arc<SharedState>,
) -> Result<(f64, f64), Box<dyn Error>> {
    let ηthresholds: Vec<f64> = weights.iter().map(|w| *w as f64).collect();

    let eme2k = state.almanac.frame_info(EARTH_J2000).unwrap();
    let epoch = Epoch::from_gregorian_utc_hms(2024, 2, 29, 12, 13, 14);

    let orbit = Orbit::keplerian(24505.9, 0.725, 7.05, 0.0, 0.0, 0.0, epoch, eme2k);

    let sc = Spacecraft::builder()
        .orbit(orbit)
        .mass(Mass::from_dry_and_prop_masses(1000.0, 1000.0))
        .srp(SRPData::from_area(3.0 * 6.0))
        .thruster(Thruster {
            isp_s: 4435.0,
            thrust_N: 0.472,
        })
        .mode(GuidanceMode::Thrust)
        .build();

    let prop_time = prop_time_days * Unit::Day;

    let objectives = &[
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::SemiMajorAxis),
            30_000.0,
            20.0,
        ),
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::Eccentricity),
            0.001,
            5e-5,
        ),
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::Inclination),
            0.05,
            1e-2,
        ),
    ];

    // let kluever_ctrl = Kluever::from_max_eclipse(objectives, &weights_f64, 0.2);
    let ctrl = Ruggiero::from_ηthresholds(objectives, &ηthresholds, sc)?;

    let mut orbital_dyn = OrbitalDynamics::point_masses(vec![MOON, SUN]);
    orbital_dyn.accel_models.push(state.harmonics.clone());

    let sc_dynamics = SpacecraftDynamics::from_model(orbital_dyn, state.srp_dyn.clone())
        .with_guidance_law(ctrl.clone());

    let (final_state, _traj) = Propagator::rk89(
        sc_dynamics.clone(),
        IntegratorOptions::builder()
            .min_step(10.0_f64.seconds())
            .tolerance(1e-8)
            .error_ctrl(ErrorControl::RSSCartesianStep)
            .build(),
    )
    .with(sc, state.almanac.clone())
    .for_duration_with_traj(prop_time)?;

    let prop_usage = sc.mass.prop_mass_kg - final_state.mass.prop_mass_kg;

    let mut penalty = 0.0;
    for obj in objectives {
        let (achieved, error) = obj.assess(&final_state)?;
        if !achieved {
            penalty += error.abs();
        }
        info!("{obj} error: {error:.3}, achieved? {achieved}");
    }

    info!("{ηthresholds:?} -> {prop_usage:.3} kg\tpenalty = {penalty:.3}");

    Ok((prop_usage, penalty * 1000.0))
}

nyx-core/examples/03_geo_analysis/stationkeeping.rs (line 122)

fn main() -> Result<(), Box<dyn Error>> {
    pel::init();
    // Set up the dynamics like in the orbit raise.
    let almanac = Arc::new(MetaAlmanac::latest().map_err(Box::new)?);
    let epoch = Epoch::from_gregorian_utc_hms(2024, 2, 29, 12, 13, 14);

    // Define the GEO orbit, and we're just going to maintain it very tightly.
    let earth_j2000 = almanac.frame_info(EARTH_J2000)?;
    let orbit = Orbit::try_keplerian(42164.0, 1e-5, 0., 163.0, 75.0, 0.0, epoch, earth_j2000)?;
    println!("{orbit:x}");

    let sc = Spacecraft::builder()
        .orbit(orbit)
        .mass(Mass::from_dry_and_prop_masses(1000.0, 1000.0)) // 1000 kg of dry mass and prop, totalling 2.0 tons
        .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
        .thruster(Thruster {
            // "NEXT-STEP" row in Table 2
            isp_s: 4435.0,
            thrust_N: 0.472,
        })
        .mode(GuidanceMode::Thrust) // Start thrusting immediately.
        .build();

    // Set up the spacecraft dynamics like in the orbit raise example.

    let prop_time = 30.0 * Unit::Day;

    // Define the guidance law -- we're just using a Ruggiero controller as demonstrated in AAS-2004-5089.
    let objectives = &[
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::SemiMajorAxis),
            42_165.0,
            20.0,
        ),
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::Eccentricity),
            0.001,
            5e-5,
        ),
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::Inclination),
            0.05,
            1e-2,
        ),
    ];

    let ruggiero_ctrl = Ruggiero::from_max_eclipse(objectives, sc, 0.2)?;
    println!("{ruggiero_ctrl}");

    let mut orbital_dyn = OrbitalDynamics::point_masses(vec![MOON, SUN]);

    let mut jgm3_meta = MetaFile {
        uri: "http://public-data.nyxspace.com/nyx/models/JGM3.cof.gz".to_string(),
        crc32: Some(0xF446F027), // Specifying the CRC32 avoids redownloading it if it's cached.
    };
    jgm3_meta.process(true)?;

    let harmonics = GravityField::from_stor(
        almanac.frame_info(IAU_EARTH_FRAME)?,
        GravityFieldData::from_cof(&jgm3_meta.uri, 8, 8, true)?,
    );
    orbital_dyn.accel_models.push(harmonics);

    let srp_dyn = SolarPressure::default_flux(EARTH_J2000, almanac.clone())?;
    let sc_dynamics = SpacecraftDynamics::from_model(orbital_dyn, srp_dyn)
        .with_guidance_law(ruggiero_ctrl.clone());

    println!("{sc_dynamics}");

    // Finally, let's use the Monte Carlo framework built into Nyx to propagate spacecraft.

    // Let's start by defining the dispersion.
    // 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.
    // Note that additional validation on the MVN is in progress -- https://github.com/nyx-space/nyx/issues/339.
    let mc_rv = MvnSpacecraft::new(
        sc,
        vec![StateDispersion::zero_mean(
            StateParameter::Element(OrbitalElement::SemiMajorAxis),
            3.0,
        )],
    )?;

    let my_mc = MonteCarlo::new(
        sc, // Nominal state
        mc_rv,
        "03_geo_sk".to_string(), // Scenario name
        None, // No specific seed specified, so one will be drawn from the computer's entropy.
    );

    // Build the propagator setup.
    let setup = Propagator::rk89(
        sc_dynamics.clone(),
        IntegratorOptions::builder()
            .min_step(10.0_f64.seconds())
            .error_ctrl(ErrorControl::RSSCartesianStep)
            .build(),
    );

    let num_runs = 25;
    let rslts = my_mc.run_until_epoch(setup, almanac.clone(), sc.epoch() + prop_time, num_runs);

    assert_eq!(rslts.runs.len(), num_runs);

    rslts.to_parquet("03_geo_sk.parquet", ExportCfg::default())?;

    Ok(())
}

nyx-core/examples/03_geo_analysis/raise.rs (line 133)

fn main() -> Result<(), Box<dyn Error>> {
    pel::init();

    // Dynamics models require planetary constants and ephemerides to be defined.
    // Let's start by grabbing those by using ANISE's latest MetaAlmanac.
    // This will automatically download the DE440s planetary ephemeris,
    // the daily-updated Earth Orientation Parameters, the high fidelity Moon orientation
    // parameters (for the Moon Mean Earth and Moon Principal Axes frames), and the PCK11
    // planetary constants kernels.
    // For details, refer to https://github.com/nyx-space/anise/blob/master/data/latest.dhall.
    // Note that we place the Almanac into an Arc so we can clone it cheaply and provide read-only
    // references to many functions.
    let almanac = Arc::new(MetaAlmanac::latest().map_err(Box::new)?);
    // Fetch the EME2000 frame from the Almabac
    let eme2k = almanac.frame_info(EARTH_J2000).unwrap();
    // Define the orbit epoch
    let epoch = Epoch::from_gregorian_utc_hms(2024, 2, 29, 12, 13, 14);

    // Build the spacecraft itself.
    // Using slide 6 of https://aerospace.org/sites/default/files/2018-11/Davis-Mayberry_HPSEP_11212018.pdf
    // for the "next gen" SEP characteristics.

    // GTO start
    let orbit = Orbit::keplerian(24505.9, 0.725, 7.05, 0.0, 0.0, 0.0, epoch, eme2k);

    let sc = Spacecraft::builder()
        .orbit(orbit)
        .mass(Mass::from_dry_and_prop_masses(1000.0, 1000.0)) // 1000 kg of dry mass and prop, totalling 2.0 tons
        .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
        .thruster(Thruster {
            // "NEXT-STEP" row in Table 2
            isp_s: 4435.0,
            thrust_N: 0.472,
        })
        .mode(GuidanceMode::Thrust) // Start thrusting immediately.
        .build();

    let prop_time = 180.0 * Unit::Day;

    // Define the guidance law -- we're just using a Ruggiero controller as demonstrated in AAS-2004-5089.
    let objectives = &[
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::SemiMajorAxis),
            42_165.0,
            20.0,
        ),
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::Eccentricity),
            0.001,
            5e-5,
        ),
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::Inclination),
            0.05,
            1e-2,
        ),
    ];

    // Ensure that we only thrust if we have more than 20% illumination.
    let ruggiero_ctrl = Ruggiero::from_max_eclipse(objectives, sc, 0.2).unwrap();
    println!("{ruggiero_ctrl}");

    // Define the high fidelity dynamics

    // Set up the spacecraft dynamics.

    // Specify that the orbital dynamics must account for the graviational pull of the Moon and the Sun.
    // The gravity of the Earth will also be accounted for since the spaceraft in an Earth orbit.
    let mut orbital_dyn = OrbitalDynamics::point_masses(vec![MOON, SUN]);

    // We want to include the spherical harmonics, so let's download the gravitational data from the Nyx Cloud.
    // We're using the JGM3 model here, which is the default in GMAT.
    let mut jgm3_meta = MetaFile {
        uri: "http://public-data.nyxspace.com/nyx/models/JGM3.cof.gz".to_string(),
        crc32: Some(0xF446F027), // Specifying the CRC32 avoids redownloading it if it's cached.
    };
    // And let's download it if we don't have it yet.
    jgm3_meta.process(true)?;

    // Build the spherical harmonics.
    // The harmonics must be computed in the body fixed frame.
    // We're using the long term prediction of the Earth centered Earth fixed frame, IAU Earth.
    let harmonics = GravityField::from_stor(
        almanac.frame_info(IAU_EARTH_FRAME)?,
        GravityFieldData::from_cof(&jgm3_meta.uri, 8, 8, true).unwrap(),
    );

    // Include the spherical harmonics into the orbital dynamics.
    orbital_dyn.accel_models.push(harmonics);

    // We define the solar radiation pressure, using the default solar flux and accounting only
    // for the eclipsing caused by the Earth.
    let srp_dyn = SolarPressure::default_flux(EARTH_J2000, almanac.clone())?;

    // Finalize setting up the dynamics, specifying the force models (orbital_dyn) separately from the
    // acceleration models (SRP in this case). Use `from_models` to specify multiple accel models.
    let sc_dynamics = SpacecraftDynamics::from_model(orbital_dyn, srp_dyn)
        .with_guidance_law(ruggiero_ctrl.clone());

    println!("{orbit:x}");

    // We specify a minimum step in the propagator because the Ruggiero control would otherwise drive this step very low.
    let (final_state, traj) = Propagator::rk89(
        sc_dynamics.clone(),
        IntegratorOptions::builder()
            .min_step(10.0_f64.seconds())
            .error_ctrl(ErrorControl::RSSCartesianStep)
            .build(),
    )
    .with(sc, almanac.clone())
    .for_duration_with_traj(prop_time)?;

    let prop_usage = sc.mass.prop_mass_kg - final_state.mass.prop_mass_kg;
    println!("{:x}", final_state.orbit);
    println!("prop usage: {prop_usage:.3} kg");

    // Finally, export the results for analysis, including the penumbra percentage throughout the orbit raise.
    traj.to_parquet("./03_geo_raise.parquet", ExportCfg::default())?;

    for status_line in ruggiero_ctrl.status(&final_state) {
        println!("{status_line}");
    }

    ruggiero_ctrl
        .achieved(&final_state)
        .expect("objective not achieved");

    Ok(())
}

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, ())>

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.

impl<__init_step, __min_step, __max_step, __tolerance, __attempts, __fixed_step, __error_ctrl, __integration_frame> IntegratorOptionsBuilder<(__init_step, __min_step, __max_step, __tolerance, __attempts, __fixed_step, __error_ctrl, __integration_frame)>

where IntegratorOptions: for<'__typed_builder_lifetime_for_default> NextFieldDefault<(__init_step,), Output = Duration> + for<'__typed_builder_lifetime_for_default> NextFieldDefault<(&'__typed_builder_lifetime_for_default Duration, __min_step), Output = Duration> + for<'__typed_builder_lifetime_for_default> NextFieldDefault<(&'__typed_builder_lifetime_for_default Duration, &'__typed_builder_lifetime_for_default Duration, __max_step), Output = Duration> + for<'__typed_builder_lifetime_for_default> NextFieldDefault<(&'__typed_builder_lifetime_for_default Duration, &'__typed_builder_lifetime_for_default Duration, &'__typed_builder_lifetime_for_default Duration, __tolerance), Output = f64> + for<'__typed_builder_lifetime_for_default> NextFieldDefault<(&'__typed_builder_lifetime_for_default Duration, &'__typed_builder_lifetime_for_default Duration, &'__typed_builder_lifetime_for_default Duration, &'__typed_builder_lifetime_for_default f64, __attempts), Output = u8> + for<'__typed_builder_lifetime_for_default> NextFieldDefault<(&'__typed_builder_lifetime_for_default Duration, &'__typed_builder_lifetime_for_default Duration, &'__typed_builder_lifetime_for_default Duration, &'__typed_builder_lifetime_for_default f64, &'__typed_builder_lifetime_for_default u8, __fixed_step), Output = bool> + for<'__typed_builder_lifetime_for_default> NextFieldDefault<(&'__typed_builder_lifetime_for_default Duration, &'__typed_builder_lifetime_for_default Duration, &'__typed_builder_lifetime_for_default Duration, &'__typed_builder_lifetime_for_default f64, &'__typed_builder_lifetime_for_default u8, &'__typed_builder_lifetime_for_default bool, __error_ctrl), Output = ErrorControl> + for<'__typed_builder_lifetime_for_default> NextFieldDefault<(&'__typed_builder_lifetime_for_default Duration, &'__typed_builder_lifetime_for_default Duration, &'__typed_builder_lifetime_for_default Duration, &'__typed_builder_lifetime_for_default f64, &'__typed_builder_lifetime_for_default u8, &'__typed_builder_lifetime_for_default bool, &'__typed_builder_lifetime_for_default ErrorControl, __integration_frame), Output = Option<Frame>>,

pub fn build(self) -> IntegratorOptions

Finalise the builder and create its IntegratorOptions instance

Examples found in repository

nyx-core/examples/03_geo_analysis/raise_optim.rs (line 187)

fn evaluate_weights(
    weights: &[f32],
    prop_time_days: f64,
    state: Arc<SharedState>,
) -> Result<(f64, f64), Box<dyn Error>> {
    let ηthresholds: Vec<f64> = weights.iter().map(|w| *w as f64).collect();

    let eme2k = state.almanac.frame_info(EARTH_J2000).unwrap();
    let epoch = Epoch::from_gregorian_utc_hms(2024, 2, 29, 12, 13, 14);

    let orbit = Orbit::keplerian(24505.9, 0.725, 7.05, 0.0, 0.0, 0.0, epoch, eme2k);

    let sc = Spacecraft::builder()
        .orbit(orbit)
        .mass(Mass::from_dry_and_prop_masses(1000.0, 1000.0))
        .srp(SRPData::from_area(3.0 * 6.0))
        .thruster(Thruster {
            isp_s: 4435.0,
            thrust_N: 0.472,
        })
        .mode(GuidanceMode::Thrust)
        .build();

    let prop_time = prop_time_days * Unit::Day;

    let objectives = &[
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::SemiMajorAxis),
            30_000.0,
            20.0,
        ),
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::Eccentricity),
            0.001,
            5e-5,
        ),
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::Inclination),
            0.05,
            1e-2,
        ),
    ];

    // let kluever_ctrl = Kluever::from_max_eclipse(objectives, &weights_f64, 0.2);
    let ctrl = Ruggiero::from_ηthresholds(objectives, &ηthresholds, sc)?;

    let mut orbital_dyn = OrbitalDynamics::point_masses(vec![MOON, SUN]);
    orbital_dyn.accel_models.push(state.harmonics.clone());

    let sc_dynamics = SpacecraftDynamics::from_model(orbital_dyn, state.srp_dyn.clone())
        .with_guidance_law(ctrl.clone());

    let (final_state, _traj) = Propagator::rk89(
        sc_dynamics.clone(),
        IntegratorOptions::builder()
            .min_step(10.0_f64.seconds())
            .tolerance(1e-8)
            .error_ctrl(ErrorControl::RSSCartesianStep)
            .build(),
    )
    .with(sc, state.almanac.clone())
    .for_duration_with_traj(prop_time)?;

    let prop_usage = sc.mass.prop_mass_kg - final_state.mass.prop_mass_kg;

    let mut penalty = 0.0;
    for obj in objectives {
        let (achieved, error) = obj.assess(&final_state)?;
        if !achieved {
            penalty += error.abs();
        }
        info!("{obj} error: {error:.3}, achieved? {achieved}");
    }

    info!("{ηthresholds:?} -> {prop_usage:.3} kg\tpenalty = {penalty:.3}");

    Ok((prop_usage, penalty * 1000.0))
}

nyx-core/examples/03_geo_analysis/stationkeeping.rs (line 123)

fn main() -> Result<(), Box<dyn Error>> {
    pel::init();
    // Set up the dynamics like in the orbit raise.
    let almanac = Arc::new(MetaAlmanac::latest().map_err(Box::new)?);
    let epoch = Epoch::from_gregorian_utc_hms(2024, 2, 29, 12, 13, 14);

    // Define the GEO orbit, and we're just going to maintain it very tightly.
    let earth_j2000 = almanac.frame_info(EARTH_J2000)?;
    let orbit = Orbit::try_keplerian(42164.0, 1e-5, 0., 163.0, 75.0, 0.0, epoch, earth_j2000)?;
    println!("{orbit:x}");

    let sc = Spacecraft::builder()
        .orbit(orbit)
        .mass(Mass::from_dry_and_prop_masses(1000.0, 1000.0)) // 1000 kg of dry mass and prop, totalling 2.0 tons
        .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
        .thruster(Thruster {
            // "NEXT-STEP" row in Table 2
            isp_s: 4435.0,
            thrust_N: 0.472,
        })
        .mode(GuidanceMode::Thrust) // Start thrusting immediately.
        .build();

    // Set up the spacecraft dynamics like in the orbit raise example.

    let prop_time = 30.0 * Unit::Day;

    // Define the guidance law -- we're just using a Ruggiero controller as demonstrated in AAS-2004-5089.
    let objectives = &[
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::SemiMajorAxis),
            42_165.0,
            20.0,
        ),
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::Eccentricity),
            0.001,
            5e-5,
        ),
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::Inclination),
            0.05,
            1e-2,
        ),
    ];

    let ruggiero_ctrl = Ruggiero::from_max_eclipse(objectives, sc, 0.2)?;
    println!("{ruggiero_ctrl}");

    let mut orbital_dyn = OrbitalDynamics::point_masses(vec![MOON, SUN]);

    let mut jgm3_meta = MetaFile {
        uri: "http://public-data.nyxspace.com/nyx/models/JGM3.cof.gz".to_string(),
        crc32: Some(0xF446F027), // Specifying the CRC32 avoids redownloading it if it's cached.
    };
    jgm3_meta.process(true)?;

    let harmonics = GravityField::from_stor(
        almanac.frame_info(IAU_EARTH_FRAME)?,
        GravityFieldData::from_cof(&jgm3_meta.uri, 8, 8, true)?,
    );
    orbital_dyn.accel_models.push(harmonics);

    let srp_dyn = SolarPressure::default_flux(EARTH_J2000, almanac.clone())?;
    let sc_dynamics = SpacecraftDynamics::from_model(orbital_dyn, srp_dyn)
        .with_guidance_law(ruggiero_ctrl.clone());

    println!("{sc_dynamics}");

    // Finally, let's use the Monte Carlo framework built into Nyx to propagate spacecraft.

    // Let's start by defining the dispersion.
    // 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.
    // Note that additional validation on the MVN is in progress -- https://github.com/nyx-space/nyx/issues/339.
    let mc_rv = MvnSpacecraft::new(
        sc,
        vec![StateDispersion::zero_mean(
            StateParameter::Element(OrbitalElement::SemiMajorAxis),
            3.0,
        )],
    )?;

    let my_mc = MonteCarlo::new(
        sc, // Nominal state
        mc_rv,
        "03_geo_sk".to_string(), // Scenario name
        None, // No specific seed specified, so one will be drawn from the computer's entropy.
    );

    // Build the propagator setup.
    let setup = Propagator::rk89(
        sc_dynamics.clone(),
        IntegratorOptions::builder()
            .min_step(10.0_f64.seconds())
            .error_ctrl(ErrorControl::RSSCartesianStep)
            .build(),
    );

    let num_runs = 25;
    let rslts = my_mc.run_until_epoch(setup, almanac.clone(), sc.epoch() + prop_time, num_runs);

    assert_eq!(rslts.runs.len(), num_runs);

    rslts.to_parquet("03_geo_sk.parquet", ExportCfg::default())?;

    Ok(())
}

nyx-core/examples/03_geo_analysis/raise.rs (line 134)

fn main() -> Result<(), Box<dyn Error>> {
    pel::init();

    // Dynamics models require planetary constants and ephemerides to be defined.
    // Let's start by grabbing those by using ANISE's latest MetaAlmanac.
    // This will automatically download the DE440s planetary ephemeris,
    // the daily-updated Earth Orientation Parameters, the high fidelity Moon orientation
    // parameters (for the Moon Mean Earth and Moon Principal Axes frames), and the PCK11
    // planetary constants kernels.
    // For details, refer to https://github.com/nyx-space/anise/blob/master/data/latest.dhall.
    // Note that we place the Almanac into an Arc so we can clone it cheaply and provide read-only
    // references to many functions.
    let almanac = Arc::new(MetaAlmanac::latest().map_err(Box::new)?);
    // Fetch the EME2000 frame from the Almabac
    let eme2k = almanac.frame_info(EARTH_J2000).unwrap();
    // Define the orbit epoch
    let epoch = Epoch::from_gregorian_utc_hms(2024, 2, 29, 12, 13, 14);

    // Build the spacecraft itself.
    // Using slide 6 of https://aerospace.org/sites/default/files/2018-11/Davis-Mayberry_HPSEP_11212018.pdf
    // for the "next gen" SEP characteristics.

    // GTO start
    let orbit = Orbit::keplerian(24505.9, 0.725, 7.05, 0.0, 0.0, 0.0, epoch, eme2k);

    let sc = Spacecraft::builder()
        .orbit(orbit)
        .mass(Mass::from_dry_and_prop_masses(1000.0, 1000.0)) // 1000 kg of dry mass and prop, totalling 2.0 tons
        .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
        .thruster(Thruster {
            // "NEXT-STEP" row in Table 2
            isp_s: 4435.0,
            thrust_N: 0.472,
        })
        .mode(GuidanceMode::Thrust) // Start thrusting immediately.
        .build();

    let prop_time = 180.0 * Unit::Day;

    // Define the guidance law -- we're just using a Ruggiero controller as demonstrated in AAS-2004-5089.
    let objectives = &[
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::SemiMajorAxis),
            42_165.0,
            20.0,
        ),
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::Eccentricity),
            0.001,
            5e-5,
        ),
        Objective::within_tolerance(
            StateParameter::Element(OrbitalElement::Inclination),
            0.05,
            1e-2,
        ),
    ];

    // Ensure that we only thrust if we have more than 20% illumination.
    let ruggiero_ctrl = Ruggiero::from_max_eclipse(objectives, sc, 0.2).unwrap();
    println!("{ruggiero_ctrl}");

    // Define the high fidelity dynamics

    // Set up the spacecraft dynamics.

    // Specify that the orbital dynamics must account for the graviational pull of the Moon and the Sun.
    // The gravity of the Earth will also be accounted for since the spaceraft in an Earth orbit.
    let mut orbital_dyn = OrbitalDynamics::point_masses(vec![MOON, SUN]);

    // We want to include the spherical harmonics, so let's download the gravitational data from the Nyx Cloud.
    // We're using the JGM3 model here, which is the default in GMAT.
    let mut jgm3_meta = MetaFile {
        uri: "http://public-data.nyxspace.com/nyx/models/JGM3.cof.gz".to_string(),
        crc32: Some(0xF446F027), // Specifying the CRC32 avoids redownloading it if it's cached.
    };
    // And let's download it if we don't have it yet.
    jgm3_meta.process(true)?;

    // Build the spherical harmonics.
    // The harmonics must be computed in the body fixed frame.
    // We're using the long term prediction of the Earth centered Earth fixed frame, IAU Earth.
    let harmonics = GravityField::from_stor(
        almanac.frame_info(IAU_EARTH_FRAME)?,
        GravityFieldData::from_cof(&jgm3_meta.uri, 8, 8, true).unwrap(),
    );

    // Include the spherical harmonics into the orbital dynamics.
    orbital_dyn.accel_models.push(harmonics);

    // We define the solar radiation pressure, using the default solar flux and accounting only
    // for the eclipsing caused by the Earth.
    let srp_dyn = SolarPressure::default_flux(EARTH_J2000, almanac.clone())?;

    // Finalize setting up the dynamics, specifying the force models (orbital_dyn) separately from the
    // acceleration models (SRP in this case). Use `from_models` to specify multiple accel models.
    let sc_dynamics = SpacecraftDynamics::from_model(orbital_dyn, srp_dyn)
        .with_guidance_law(ruggiero_ctrl.clone());

    println!("{orbit:x}");

    // We specify a minimum step in the propagator because the Ruggiero control would otherwise drive this step very low.
    let (final_state, traj) = Propagator::rk89(
        sc_dynamics.clone(),
        IntegratorOptions::builder()
            .min_step(10.0_f64.seconds())
            .error_ctrl(ErrorControl::RSSCartesianStep)
            .build(),
    )
    .with(sc, almanac.clone())
    .for_duration_with_traj(prop_time)?;

    let prop_usage = sc.mass.prop_mass_kg - final_state.mass.prop_mass_kg;
    println!("{:x}", final_state.orbit);
    println!("prop usage: {prop_usage:.3} kg");

    // Finally, export the results for analysis, including the penumbra percentage throughout the orbit raise.
    traj.to_parquet("./03_geo_raise.parquet", ExportCfg::default())?;

    for status_line in ruggiero_ctrl.status(&final_state) {
        println!("{status_line}");
    }

    ruggiero_ctrl
        .achieved(&final_state)
        .expect("objective not achieved");

    Ok(())
}

nyx-core/examples/05_cislunar_spacecraft_link_od/main.rs (line 85)

fn main() -> Result<(), Box<dyn Error>> {
    pel::init();

    // ====================== //
    // === ALMANAC SET UP === //
    // ====================== //

    let manifest_dir = PathBuf::from(env!("CARGO_MANIFEST_DIR"));

    let out = manifest_dir.join("data/04_output/");

    let almanac = Arc::new(
        Almanac::new(
            &manifest_dir
                .join("data/01_planetary/pck08.pca")
                .to_string_lossy(),
        )
        .unwrap()
        .load(
            &manifest_dir
                .join("data/01_planetary/de440s.bsp")
                .to_string_lossy(),
        )
        .unwrap(),
    );

    let eme2k = almanac.frame_info(EARTH_J2000).unwrap();
    let moon_iau = almanac.frame_info(IAU_MOON_FRAME).unwrap();

    let epoch = Epoch::from_gregorian_tai(2021, 5, 29, 19, 51, 16, 852_000);
    let nrho = Orbit::cartesian(
        166_473.631_302_239_7,
        -274_715.487_253_382_7,
        -211_233.210_176_686_7,
        0.933_451_604_520_018_4,
        0.436_775_046_841_900_9,
        -0.082_211_021_250_348_95,
        epoch,
        eme2k,
    );

    let tx_nrho_sc = Spacecraft::from(nrho);

    let state_luna = almanac.transform_to(nrho, MOON_J2000, None).unwrap();
    println!("Start state (dynamics: Earth, Moon, Sun gravity):\n{state_luna}");

    let bodies = vec![EARTH, SUN];
    let dynamics = SpacecraftDynamics::new(OrbitalDynamics::point_masses(bodies));

    let setup = Propagator::rk89(
        dynamics,
        IntegratorOptions::builder().max_step(0.5.minutes()).build(),
    );

    /* == Propagate the NRHO vehicle == */
    let prop_time = 1.1 * state_luna.period().unwrap();

    let (nrho_final, mut tx_traj) = setup
        .with(tx_nrho_sc, almanac.clone())
        .for_duration_with_traj(prop_time)
        .unwrap();

    tx_traj.name = Some("NRHO Tx SC".to_string());

    println!("{tx_traj}");

    /* == Propagate an LLO vehicle == */
    let llo_orbit =
        Orbit::try_keplerian_altitude(110.0, 1e-4, 90.0, 0.0, 0.0, 0.0, epoch, moon_iau).unwrap();

    let llo_sc = Spacecraft::builder().orbit(llo_orbit).build();

    let (_, llo_traj) = setup
        .with(llo_sc, almanac.clone())
        .until_epoch_with_traj(nrho_final.epoch())
        .unwrap();

    // Export the subset of the first two hours.
    llo_traj
        .clone()
        .filter_by_offset(..2.hours())
        .to_parquet_simple(out.join("05_caps_llo_truth.pq"))?;

    /* == Setup the interlink == */

    let mut measurement_types = IndexSet::new();
    measurement_types.insert(MeasurementType::Range);
    measurement_types.insert(MeasurementType::Doppler);

    let mut stochastics = IndexMap::new();

    let sa45_csac_allan_dev = 1e-11;

    stochastics.insert(
        MeasurementType::Range,
        StochasticNoise::from_hardware_range_km(
            sa45_csac_allan_dev,
            10.0.seconds(),
            link_specific::ChipRate::StandardT4B,
            link_specific::SN0::Average,
        ),
    );

    stochastics.insert(
        MeasurementType::Doppler,
        StochasticNoise::from_hardware_doppler_km_s(
            sa45_csac_allan_dev,
            10.0.seconds(),
            link_specific::CarrierFreq::SBand,
            link_specific::CN0::Average,
        ),
    );

    let interlink = InterlinkTxSpacecraft {
        traj: tx_traj,
        measurement_types,
        integration_time: None,
        timestamp_noise_s: None,
        ab_corr: Aberration::LT,
        stochastic_noises: Some(stochastics),
    };

    // Devices are the transmitter, which is our NRHO vehicle.
    let mut devices = BTreeMap::new();
    devices.insert("NRHO Tx SC".to_string(), interlink);

    let mut configs = BTreeMap::new();
    configs.insert(
        "NRHO Tx SC".to_string(),
        TrkConfig::builder()
            .strands(vec![Strand {
                start: epoch,
                end: nrho_final.epoch(),
            }])
            .build(),
    );

    let mut trk_sim =
        TrackingArcSim::with_seed(devices.clone(), llo_traj.clone(), configs, 0).unwrap();
    println!("{trk_sim}");

    let trk_data = trk_sim.generate_measurements(almanac.clone()).unwrap();
    println!("{trk_data}");

    trk_data
        .to_parquet_simple(out.clone().join("nrho_interlink_msr.pq"))
        .unwrap();

    // Run a truth OD where we estimate the LLO position
    let llo_uncertainty = SpacecraftUncertainty::builder()
        .nominal(llo_sc)
        .x_km(1.0)
        .y_km(1.0)
        .z_km(1.0)
        .vx_km_s(1e-3)
        .vy_km_s(1e-3)
        .vz_km_s(1e-3)
        .build();

    let mut proc_devices = devices.clone();

    // Define the initial estimate, randomized, seed for reproducibility
    let mut initial_estimate = llo_uncertainty.to_estimate_randomized(Some(0)).unwrap();
    // Inflate the covariance -- https://github.com/nyx-space/nyx/issues/339
    initial_estimate.covar *= 2.5;

    // Increase the noise in the devices to accept more measurements.

    for link in proc_devices.values_mut() {
        for noise in &mut link.stochastic_noises.as_mut().unwrap().values_mut() {
            *noise.white_noise.as_mut().unwrap() *= 3.0;
        }
    }

    let init_err = initial_estimate
        .orbital_state()
        .ric_difference(&llo_orbit)
        .unwrap();

    println!("initial estimate:\n{initial_estimate}");
    println!("RIC errors = {init_err}",);

    let odp = InterlinkKalmanOD::new(
        setup.clone(),
        KalmanVariant::ReferenceUpdate,
        Some(ResidRejectCrit::default()),
        proc_devices,
        almanac.clone(),
    );

    // Shrink the data to process.
    let arc = trk_data.filter_by_offset(..2.hours());

    let od_sol = odp.process_arc(initial_estimate, &arc).unwrap();

    println!("{od_sol}");

    od_sol
        .to_parquet(
            out.join("05_caps_interlink_od_sol.pq"),
            ExportCfg::default(),
        )
        .unwrap();

    let od_traj = od_sol.to_traj().unwrap();

    od_traj
        .ric_diff_to_parquet(
            &llo_traj,
            out.join("05_caps_interlink_llo_est_error.pq"),
            ExportCfg::default(),
        )
        .unwrap();

    let final_est = od_sol.estimates.last().unwrap();
    assert!(final_est.within_3sigma(), "should be within 3 sigma");

    println!("ESTIMATE\n{final_est:x}\n");
    let truth = llo_traj.at(final_est.epoch()).unwrap();
    println!("TRUTH\n{truth:x}");

    let final_err = truth
        .orbit
        .ric_difference(&final_est.orbital_state())
        .unwrap();
    println!("ERROR {final_err}");

    // Build the residuals versus reference plot.
    let rvr_sol = odp
        .process_arc(initial_estimate, &arc.resid_vs_ref_check())
        .unwrap();

    rvr_sol
        .to_parquet(
            out.join("05_caps_interlink_resid_v_ref.pq"),
            ExportCfg::default(),
        )
        .unwrap();

    let final_rvr = rvr_sol.estimates.last().unwrap();

    println!("RMAG error {:.3} m", final_err.rmag_km() * 1e3);
    println!(
        "Pure prop error {:.3} m",
        final_rvr
            .orbital_state()
            .ric_difference(&final_est.orbital_state())
            .unwrap()
            .rmag_km()
            * 1e3
    );

    Ok(())
}

Trait Implementations

impl<TypedBuilderFields> Clone for IntegratorOptionsBuilder<TypedBuilderFields>

where TypedBuilderFields: Clone,

fn clone(&self) -> Self

Returns a duplicate of the value.

fn clone_from(&mut self, source: &Self)

Performs copy-assignment from source.