nyx_space/md/opti/multipleshooting/
equidistant_heuristic.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 super::ctrlnodes::Node;
20use super::multishoot::MultipleShooting;
21pub use super::CostFunction;
22use crate::errors::TargetingError;
23use crate::md::prelude::*;
24use crate::{Orbit, Spacecraft};
25
26impl<'a> MultipleShooting<'a, Node, 3, 3> {
27    /// Builds a multiple shooting structure assuming that the optimal trajectory is a straight line
28    /// between the start and end points. The position of the nodes will be update at each iteration
29    /// of the outer loop.
30    /// NOTE: this may cause some nodes to be below the surface of a celestial object if in low orbit
31    pub fn equidistant_nodes(
32        x0: Spacecraft,
33        xf: Orbit,
34        node_count: usize,
35        prop: &'a Propagator<SpacecraftDynamics>,
36    ) -> Result<Self, TargetingError> {
37        if node_count < 3 {
38            error!("At least three nodes are needed for a multiple shooting optimization");
39            return Err(TargetingError::UnderdeterminedProblem);
40        }
41
42        // Compute the direction of the objective
43        let mut direction = xf.radius_km - x0.orbit.radius_km;
44        if direction.norm() < f64::EPSILON {
45            return Err(TargetingError::TargetsTooClose);
46        }
47        let distance_increment = direction.norm() / (node_count as f64);
48        let duration_increment = (xf.epoch - x0.epoch()) / (node_count as f64);
49        direction /= direction.norm();
50
51        // Build each node successively (includes xf)
52        let mut nodes = Vec::with_capacity(node_count + 1);
53        let mut prev_node_radius = x0.orbit.radius_km;
54        let mut prev_node_epoch = x0.epoch();
55
56        for _ in 0..node_count {
57            // Compute the position we want.
58            let this_node = prev_node_radius + distance_increment * direction;
59            let this_epoch = prev_node_epoch + duration_increment;
60            nodes.push(Node {
61                x: this_node[0],
62                y: this_node[1],
63                z: this_node[2],
64                vmag: 0.0,
65                frame: x0.orbit.frame,
66                epoch: this_epoch,
67            });
68            prev_node_radius = this_node;
69            prev_node_epoch = this_epoch;
70        }
71        Ok(Self {
72            prop,
73            targets: nodes,
74            x0,
75            xf,
76            current_iteration: 0,
77            max_iterations: 50,
78            improvement_threshold: 0.01,
79            variables: [
80                Vary::VelocityX.into(),
81                Vary::VelocityY.into(),
82                Vary::VelocityZ.into(),
83            ],
84            all_dvs: Vec::with_capacity(node_count),
85        })
86    }
87}