25 #include <boost/numeric/conversion/cast.hpp>
27 #include <keplerian_toolbox/epoch.h>
29 #include "../AstroToolbox/mga_dsm.h"
30 #include "../AstroToolbox/misc4Tandem.h"
31 #include "../exceptions.h"
35 namespace pagmo {
namespace problem {
37 static inline int laplace_get_dimension(
const std::vector<int> &seq)
41 pagmo_throw(value_error,
"fly-by sequence must contain at least two planets!! Earth and Jupiter");
45 pagmo_throw(value_error,
"Starting planet must be the Earth (3) for the Laplace mission");
49 pagmo_throw(value_error,
"Final planet must be Jupiter (5) for the Laplace mission");
51 for (
size_t i=0; i< seq.size(); ++i) {
52 if ((seq[i] >= 7) || (seq[i] <= 0)) {
53 pagmo_throw(value_error,
"Planet sequence must contains numbers from 1 to 6");
56 return boost::numeric_cast<
int>(4*seq.size() - 2);
59 const int laplace::default_sequence[5] = {3,2,3,3,5};
60 const int* laplace::get_default_sequence(){
return default_sequence;}
68 problem(orbit_insertion,&seq[0],seq.size(),0,0,0,.97,4 * 71492.0)
79 for (
size_t i = 0; i < seq.size() - 1; ++i)
84 for (
size_t i = 0; i < seq.size() - 1; ++i)
86 set_lb(i+3+seq.size(),0.01);
87 set_ub(i+3+seq.size(),0.99);
89 for (
size_t i = 0; i < seq.size() - 2; ++i)
91 set_lb(i+2*(seq.size()+ 1),1.05);
92 set_ub(i+2*(seq.size()+ 1),100);
94 for (
size_t i = 0; i < seq.size() - 2; ++i)
96 set_lb(i+3*seq.size(),-M_PI);
97 set_ub(i+3*seq.size(),M_PI);
111 MGA_DSM(x, problem, f[0]);
112 const size_t sequence_size = (x.size() + 2) / 4;
113 double totaltime = 0;
114 for (
size_t i = 0; i < sequence_size - 1 ; ++i) {
115 totaltime += x[i + 4];
117 const double delta = totaltime - 8 * 365.25;
119 f[0] = std::max<double>(f[0],f[0] + 0.2 / 30 * delta);
135 MGA_DSM(x, problem, obj);
136 std::ostringstream s;
138 s << std::scientific;
139 const size_t seq_size = (x.size() + 2) / 4;
140 pagmo_assert((x.size() + 2) % 4 == 0 && seq_size >= 2);
141 pagmo_assert(problem.sequence.size() == seq_size);
142 s <<
"Flyby sequence: ";
143 for (
size_t i = 0; i < seq_size; ++i) {
144 s << problem.sequence[i];
147 s <<
"Departure epoch (mjd2000): " << x[0] <<
'\n';
148 s <<
"Departure epoch: " << ::kep_toolbox::epoch(x[0],::kep_toolbox::epoch::MJD2000) <<
'\n';
149 s <<
"Vinf polar components: ";
150 for (
size_t i = 0; i < 3; ++i) {
151 s << x[i + 1] <<
' ';
154 double totaltime = 0;
155 for (
size_t i = 0; i < seq_size - 1; ++i) {
156 s <<
"Leg time of flight: " << x[i + 4] <<
'\n';
157 totaltime += x[i + 4];
159 s <<
"Total time of flight: " << totaltime <<
'\n';
160 for (
size_t i = 0; i < seq_size - 2; ++i) {
161 s <<
"Flyby radius: " << x[i + 2 * (seq_size + 1)] <<
'\n';
163 for (
size_t i = 0; i < seq_size - 2; ++i) {
164 s <<
"Vinf at flyby: " << std::sqrt(problem.vrelin_vec[i]) <<
'\n';
166 for (
size_t i = 0; i < seq_size - 1; ++i) {
167 s <<
"dsm" << i+1 <<
": " << problem.DV[i+1] <<
'\n';
boost::shared_ptr< base > base_ptr
Alias for shared pointer to base problem.
std::vector< double > decision_vector
Decision vector type.
std::string pretty(const std::vector< double > &x) const
Outputs a stream with the trajectory data.
void set_lb(const decision_vector &)
Set lower bounds from pagmo::decision_vector.
base_ptr clone() const
Clone method.
void set_ub(const decision_vector &)
Set upper bounds from pagmo::decision_vector.
std::vector< double > fitness_vector
Fitness vector type.
void objfun_impl(fitness_vector &, const decision_vector &) const
Implementation of the objective function.
laplace(const std::vector< int > &=std::vector< int >(get_default_sequence(), get_default_sequence()+5))
Constructor.
std::string get_name() const
Implementation of the sparsity structure.