nyx_space/propagators/
propagator.rs

1/*
2    Nyx, blazing fast astrodynamics
3    Copyright (C) 2018-onwards Christopher Rabotin <christopher.rabotin@gmail.com>
4
5    This program is free software: you can redistribute it and/or modify
6    it under the terms of the GNU Affero General Public License as published
7    by the Free Software Foundation, either version 3 of the License, or
8    (at your option) any later version.
9
10    This program is distributed in the hope that it will be useful,
11    but WITHOUT ANY WARRANTY; without even the implied warranty of
12    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
13    GNU Affero General Public License for more details.
14
15    You should have received a copy of the GNU Affero General Public License
16    along with this program.  If not, see <https://www.gnu.org/licenses/>.
17*/
18
19use std::sync::Arc;
20
21use anise::almanac::Almanac;
22
23use super::{IntegrationDetails, IntegratorMethod, IntegratorOptions, PropInstance};
24use crate::dynamics::Dynamics;
25use crate::linalg::allocator::Allocator;
26use crate::linalg::{DefaultAllocator, OVector};
27use crate::time::Duration;
28use crate::State;
29
30/// A Propagator allows propagating a set of dynamics forward or backward in time.
31/// It is an EventTracker, without any event tracking. It includes the options, the integrator
32/// details of the previous step, and the set of coefficients used for the monomorphic instance.
33#[derive(Clone, Debug)]
34pub struct Propagator<D: Dynamics>
35where
36    DefaultAllocator: Allocator<<D::StateType as State>::Size>
37        + Allocator<<D::StateType as State>::Size, <D::StateType as State>::Size>
38        + Allocator<<D::StateType as State>::Size, <D::StateType as State>::Size>
39        + Allocator<<D::StateType as State>::VecLength>,
40{
41    pub dynamics: D, // Stores the dynamics used. *Must* use this to get the latest values
42    pub opts: IntegratorOptions, // Stores the integration options (tolerance, min/max step, init step, etc.)
43    pub method: IntegratorMethod,
44}
45
46/// The `Propagator` trait defines the functions of a propagator and of an event tracker.
47impl<D: Dynamics> Propagator<D>
48where
49    DefaultAllocator: Allocator<<D::StateType as State>::Size>
50        + Allocator<<D::StateType as State>::Size, <D::StateType as State>::Size>
51        + Allocator<<D::StateType as State>::Size, <D::StateType as State>::Size>
52        + Allocator<<D::StateType as State>::VecLength>,
53{
54    /// Each propagator must be initialized with `new` which stores propagator information.
55    pub fn new(dynamics: D, method: IntegratorMethod, opts: IntegratorOptions) -> Self {
56        Self {
57            dynamics,
58            opts,
59            method,
60        }
61    }
62
63    /// Set the tolerance for the propagator
64    pub fn set_tolerance(&mut self, tol: f64) {
65        self.opts.tolerance = tol;
66    }
67
68    /// Set the maximum step size for the propagator and sets the initial step to that value if currently greater
69    pub fn set_max_step(&mut self, step: Duration) {
70        self.opts.set_max_step(step);
71    }
72
73    pub fn set_min_step(&mut self, step: Duration) {
74        self.opts.set_min_step(step);
75    }
76
77    /// An RK89 propagator (the default) with custom propagator options.
78    pub fn rk89(dynamics: D, opts: IntegratorOptions) -> Self {
79        Self::new(dynamics, IntegratorMethod::RungeKutta89, opts)
80    }
81
82    /// A Dormand Prince 7-8 propagator with custom propagator options: it's about 20% faster than an RK98, and more stable in two body dynamics.
83    /// WARNINGS: Dormand Prince may have issues with generating proper trajectories, leading to glitches in event finding.
84    pub fn dp78(dynamics: D, opts: IntegratorOptions) -> Self {
85        Self::new(dynamics, IntegratorMethod::DormandPrince78, opts)
86    }
87
88    pub fn with(&self, state: D::StateType, almanac: Arc<Almanac>) -> PropInstance<D> {
89        // Pre-allocate the k used in the propagator
90        let mut k = Vec::with_capacity(self.method.stages() + 1);
91        for _ in 0..self.method.stages() {
92            k.push(OVector::<f64, <D::StateType as State>::VecLength>::zeros());
93        }
94        PropInstance {
95            state,
96            prop: self,
97            details: IntegrationDetails {
98                step: self.opts.init_step,
99                error: 0.0,
100                attempts: 1,
101            },
102            log_progress: true,
103            almanac,
104            step_size: self.opts.init_step,
105            fixed_step: self.opts.fixed_step,
106            k,
107        }
108    }
109
110    /// Default propagator is an RK89 with the default PropOpts.
111    pub fn default(dynamics: D) -> Self {
112        Self::rk89(dynamics, IntegratorOptions::default())
113    }
114
115    /// A default Dormand Prince 78 propagator with the default PropOpts.
116    /// Faster and more stable than an RK89 (`default`) but seems to cause issues for event finding.
117    /// WARNINGS: Dormand Prince may have issues with generating proper trajectories, leading to glitches in event finding.
118    pub fn default_dp78(dynamics: D) -> Self {
119        Self::dp78(dynamics, IntegratorOptions::default())
120    }
121}