PaGMO  1.1.5
laplace.cpp
1 /*****************************************************************************
2  * Copyright (C) 2004-2015 The PaGMO development team, *
3  * Advanced Concepts Team (ACT), European Space Agency (ESA) *
4  * *
5  * https://github.com/esa/pagmo *
6  * *
7  * act@esa.int *
8  * *
9  * This program is free software; you can redistribute it and/or modify *
10  * it under the terms of the GNU General Public License as published by *
11  * the Free Software Foundation; either version 2 of the License, or *
12  * (at your option) any later version. *
13  * *
14  * This program is distributed in the hope that it will be useful, *
15  * but WITHOUT ANY WARRANTY; without even the implied warranty of *
16  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
17  * GNU General Public License for more details. *
18  * *
19  * You should have received a copy of the GNU General Public License *
20  * along with this program; if not, write to the *
21  * Free Software Foundation, Inc., *
22  * 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. *
23  *****************************************************************************/
24 
25 #include <boost/numeric/conversion/cast.hpp>
26 #include <string>
27 #include <keplerian_toolbox/epoch.h>
28 
29 #include "../AstroToolbox/mga_dsm.h"
30 #include "../AstroToolbox/misc4Tandem.h"
31 #include "../exceptions.h"
32 
33 #include "laplace.h"
34 
35 namespace pagmo { namespace problem {
36 
37 static inline int laplace_get_dimension(const std::vector<int> &seq)
38 {
39  if (seq.size() < 2)
40  {
41  pagmo_throw(value_error,"fly-by sequence must contain at least two planets!! Earth and Jupiter");
42  }
43  if (seq[0]!=3)
44  {
45  pagmo_throw(value_error,"Starting planet must be the Earth (3) for the Laplace mission");
46  }
47  if (seq.back()!=5)
48  {
49  pagmo_throw(value_error,"Final planet must be Jupiter (5) for the Laplace mission");
50  }
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");
54  }
55  }
56  return boost::numeric_cast<int>(4*seq.size() - 2);
57 }
58 
59 const int laplace::default_sequence[5] = {3,2,3,3,5};
60 const int* laplace::get_default_sequence(){return default_sequence;}
61 
63 
67 laplace::laplace(const std::vector<int> &seq):base(-10000.0,10000.0,laplace_get_dimension(seq)),
68  problem(orbit_insertion,&seq[0],seq.size(),0,0,0,.97,4 * 71492.0)
69 {
70  // Set the bounds.
71  set_lb(0,5475.0);
72  set_ub(0,9132.0);
73  set_lb(1,0.1);
74  set_ub(1,3.5);
75  set_lb(2,0.0);
76  set_ub(2,1.0);
77  set_lb(3,0.0);
78  set_ub(3,1.0);
79  for (size_t i = 0; i < seq.size() - 1; ++i)
80  {
81  set_lb(i+4,20);
82  set_ub(i+4,2500);
83  }
84  for (size_t i = 0; i < seq.size() - 1; ++i)
85  {
86  set_lb(i+3+seq.size(),0.01);
87  set_ub(i+3+seq.size(),0.99);
88  }
89  for (size_t i = 0; i < seq.size() - 2; ++i)
90  {
91  set_lb(i+2*(seq.size()+ 1),1.05);
92  set_ub(i+2*(seq.size()+ 1),100);
93  }
94  for (size_t i = 0; i < seq.size() - 2; ++i)
95  {
96  set_lb(i+3*seq.size(),-M_PI);
97  set_ub(i+3*seq.size(),M_PI);
98  }
99 }
100 
103 {
104  return base_ptr(new laplace(*this));
105 }
106 
109 {
110  f[0] = 0;
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];
116  }
117  const double delta = totaltime - 8 * 365.25;
118  // Penalise trajectory longer than 8 years by 200 meters/s per month.
119  f[0] = std::max<double>(f[0],f[0] + 0.2 / 30 * delta);
120 }
121 
123 
132  std::string laplace::pretty(const std::vector<double> &x) const
133 {
134  double obj = 0;
135  MGA_DSM(x, problem, obj);
136  std::ostringstream s;
137  s.precision(15);
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];
145  }
146  s << '\n';
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] << ' ';
152  }
153  s << '\n';
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];
158  }
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';
162  }
163  for (size_t i = 0; i < seq_size - 2; ++i) {
164  s << "Vinf at flyby: " << std::sqrt(problem.vrelin_vec[i]) << '\n';
165  }
166  for (size_t i = 0; i < seq_size - 1; ++i) {
167  s << "dsm" << i+1 << ": " << problem.DV[i+1] << '\n';
168  }
169  return s.str();
170 }
171 
173 
176 //void laplace::set_sparsity(int &lenG, std::vector<int> &iGfun, std::vector<int> &jGvar) const
177 //{
178 // lenG=(4*problem.sequence.size()-2);
179 // iGfun.resize(lenG);
180 // jGvar.resize(lenG);
181 // for (int i = 0; i<lenG; ++i)
182 // {
183 // iGfun[i] = 0;
184 // jGvar[i] = i;
185 // }
186 //}
187 
188 std::string laplace::get_name() const
189 {
190  return "Laplace";
191 }
192 
193 }} //namespaces
194 
195 BOOST_CLASS_EXPORT_IMPLEMENT(pagmo::problem::laplace)
Root PaGMO namespace.
boost::shared_ptr< base > base_ptr
Alias for shared pointer to base problem.
Definition: problem/base.h:62
std::vector< double > decision_vector
Decision vector type.
Definition: types.h:40
Base problem class.
Definition: problem/base.h:148
std::string pretty(const std::vector< double > &x) const
Outputs a stream with the trajectory data.
Definition: laplace.cpp:132
void set_lb(const decision_vector &)
Set lower bounds from pagmo::decision_vector.
base_ptr clone() const
Clone method.
Definition: laplace.cpp:102
void set_ub(const decision_vector &)
Set upper bounds from pagmo::decision_vector.
Laplace problem.
Definition: laplace.h:65
std::vector< double > fitness_vector
Fitness vector type.
Definition: types.h:42
void objfun_impl(fitness_vector &, const decision_vector &) const
Implementation of the objective function.
Definition: laplace.cpp:108
laplace(const std::vector< int > &=std::vector< int >(get_default_sequence(), get_default_sequence()+5))
Constructor.
Definition: laplace.cpp:67
std::string get_name() const
Implementation of the sparsity structure.
Definition: laplace.cpp:188