1 /*****************************************************************************
2 * Copyright (C) 2004-2018 The pagmo development team, *
3 * Advanced Concepts Team (ACT), European Space Agency (ESA) *
4 * http://apps.sourceforge.net/mediawiki/pagmo *
5 * http://apps.sourceforge.net/mediawiki/pagmo/index.php?title=Developers *
6 * http://apps.sourceforge.net/mediawiki/pagmo/index.php?title=Credits *
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 3 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 // NOTE: the order of inclusion in the first two items here is forced by these
26 // two issues:
27 // http://mail.python.org/pipermail/python-list/2004-March/907592.html
28 // http://mail.python.org/pipermail/new-bugs-announce/2011-March/010395.html
29 #if defined(_WIN32)
30 #include <Python.h>
31 #include <cmath>
32 #else
33 #include <Python.h>
34 #include <cmath>
35 #endif
36
37 #if PY_MAJOR_VERSION < 3
38 #error Minimum supported Python version is 2.6.
39 #endif
40
41 #include <boost/math/constants/constants.hpp>
42 #include <boost/python.hpp>
43 #include <boost/python/class.hpp>
44 #include <boost/python/converter/registry.hpp>
45 #include <boost/python/copy_const_reference.hpp>
46 #include <boost/python/def.hpp>
47 #include <boost/python/docstring_options.hpp>
48 #include <boost/python/enum.hpp>
49 #include <boost/python/module.hpp>
50 #include <boost/python/operators.hpp>
51 #include <boost/python/overloads.hpp>
52 #include <boost/python/register_ptr_to_python.hpp>
53 #include <boost/python/self.hpp>
54 #include <boost/utility.hpp>
55
56 #include "../boost_python_container_conversions.h"
57 #include "../utils.h"
58 #include "core_docstrings.hpp"
59 #include <keplerian_toolbox/keplerian_toolbox.hpp>
60
61 using namespace boost::python;
62
par2ic_wrapper(const kep_toolbox::array6D & E,const double & mu)63 static inline tuple par2ic_wrapper(const kep_toolbox::array6D &E, const double &mu)
64 {
65 kep_toolbox::array3D r0, v0;
66 kep_toolbox::par2ic(E, mu, r0, v0);
67 return boost::python::make_tuple(r0, v0);
68 }
69
ic2par_wrapper(const kep_toolbox::array3D & r0,const kep_toolbox::array3D & v0,const double & mu)70 static inline kep_toolbox::array6D ic2par_wrapper(const kep_toolbox::array3D &r0, const kep_toolbox::array3D &v0,
71 const double &mu)
72 {
73 kep_toolbox::array6D E;
74 kep_toolbox::ic2par(r0, v0, mu, E);
75 return E;
76 }
77
par2eq_wrapper(const kep_toolbox::array6D & E,const bool retrograde=false)78 static inline kep_toolbox::array6D par2eq_wrapper(const kep_toolbox::array6D &E, const bool retrograde = false)
79 {
80 kep_toolbox::array6D EQ;
81 kep_toolbox::par2eq(EQ, E, retrograde);
82 return EQ;
83 }
84
eq2par_wrapper(const kep_toolbox::array6D & EQ,const bool retrograde=false)85 static inline kep_toolbox::array6D eq2par_wrapper(const kep_toolbox::array6D &EQ, const bool retrograde = false)
86 {
87 kep_toolbox::array6D E;
88 kep_toolbox::eq2par(E, EQ, retrograde);
89 return E;
90 }
91
eq2ic_wrapper(const kep_toolbox::array6D & EQ,const double & mu,const bool retrograde=false)92 static inline tuple eq2ic_wrapper(const kep_toolbox::array6D &EQ, const double &mu, const bool retrograde = false)
93 {
94 kep_toolbox::array3D r0, v0;
95 kep_toolbox::eq2ic(EQ, mu, r0, v0, retrograde);
96 return boost::python::make_tuple(r0, v0);
97 }
98
ic2eq_wrapper(const kep_toolbox::array3D & r0,const kep_toolbox::array3D & v0,const double & mu,const bool retrograde=false)99 static inline kep_toolbox::array6D ic2eq_wrapper(const kep_toolbox::array3D &r0, const kep_toolbox::array3D &v0,
100 const double &mu, const bool retrograde = false)
101 {
102 kep_toolbox::array6D EQ;
103 kep_toolbox::ic2eq(r0, v0, mu, EQ, retrograde);
104 return EQ;
105 }
106
closest_distance_wrapper(const kep_toolbox::array3D & r0,const kep_toolbox::array3D & v0,const kep_toolbox::array3D & r1,const kep_toolbox::array3D & v1,const double & mu)107 static inline tuple closest_distance_wrapper(const kep_toolbox::array3D &r0, const kep_toolbox::array3D &v0,
108 const kep_toolbox::array3D &r1, const kep_toolbox::array3D &v1,
109 const double &mu)
110 {
111 double min_d, ra;
112 kep_toolbox::closest_distance(min_d, ra, r0, v0, r1, v1, mu);
113 return boost::python::make_tuple(min_d, ra);
114 }
115
propagate_lagrangian_wrapper(const kep_toolbox::array3D & r0,const kep_toolbox::array3D & v0,const double & t,const double & mu)116 static inline tuple propagate_lagrangian_wrapper(const kep_toolbox::array3D &r0, const kep_toolbox::array3D &v0,
117 const double &t, const double &mu)
118 {
119 kep_toolbox::array3D r(r0), v(v0);
120 kep_toolbox::propagate_lagrangian_u(r, v, t, mu);
121 return boost::python::make_tuple(r, v);
122 }
123
propagate_taylor_wrapper(const kep_toolbox::array3D & r0,const kep_toolbox::array3D & v0,const double & m0,const kep_toolbox::array3D & u,const double & t,const double & mu,const double & veff,const int & log10tolerance,const int & log10rtolerance)124 static inline tuple propagate_taylor_wrapper(const kep_toolbox::array3D &r0, const kep_toolbox::array3D &v0,
125 const double &m0, const kep_toolbox::array3D &u, const double &t,
126 const double &mu, const double &veff, const int &log10tolerance,
127 const int &log10rtolerance)
128 {
129 kep_toolbox::array3D r(r0), v(v0);
130 double m(m0);
131 kep_toolbox::propagate_taylor(r, v, m, u, t, mu, veff, log10tolerance, log10rtolerance);
132 return boost::python::make_tuple(kep_toolbox::array3D(r), kep_toolbox::array3D(v), double(m));
133 }
134
propagate_taylor_disturbance_wrapper(const kep_toolbox::array3D & r0,const kep_toolbox::array3D & v0,const double & m0,const kep_toolbox::array3D & thrust,const kep_toolbox::array3D & disturbance,const double & t,const double & mu,const double & veff,const int & log10tolerance,const int & log10rtolerance)135 static inline tuple propagate_taylor_disturbance_wrapper(const kep_toolbox::array3D &r0, const kep_toolbox::array3D &v0,
136 const double &m0, const kep_toolbox::array3D &thrust,
137 const kep_toolbox::array3D &disturbance, const double &t,
138 const double &mu, const double &veff,
139 const int &log10tolerance, const int &log10rtolerance)
140 {
141 kep_toolbox::array3D r(r0), v(v0);
142 double m(m0);
143 kep_toolbox::propagate_taylor_disturbance(r, v, m, thrust, disturbance, t, mu, veff, log10tolerance,
144 log10rtolerance);
145 return boost::python::make_tuple(kep_toolbox::array3D(r), kep_toolbox::array3D(v), double(m));
146 }
147
propagate_taylor_J2_wrapper(const kep_toolbox::array3D & r0,const kep_toolbox::array3D & v0,const double & m0,const kep_toolbox::array3D & u,const double & t,const double & mu,const double & veff,const double & J2RG2,const int & log10tolerance,const int & log10rtolerance)148 static inline tuple propagate_taylor_J2_wrapper(const kep_toolbox::array3D &r0, const kep_toolbox::array3D &v0,
149 const double &m0, const kep_toolbox::array3D &u, const double &t,
150 const double &mu, const double &veff, const double &J2RG2,
151 const int &log10tolerance, const int &log10rtolerance)
152 {
153 kep_toolbox::array3D r(r0), v(v0);
154 double m(m0);
155 kep_toolbox::propagate_taylor_J2(r, v, m, u, t, mu, veff, J2RG2, log10tolerance, log10rtolerance);
156 return boost::python::make_tuple(kep_toolbox::array3D(r), kep_toolbox::array3D(v), double(m));
157 }
158
propagate_taylor_s_wrapper(const kep_toolbox::array3D & r0,const kep_toolbox::array3D & v0,const double & m0,const kep_toolbox::array3D & u,const double & s,const double & mu=1.0,const double & veff=1.0,const double & c=1.0,const double & alpha=1.5,const int & log10tolerance=-10,const int & log10rtolerance=-10)159 static inline tuple propagate_taylor_s_wrapper(const kep_toolbox::array3D &r0, const kep_toolbox::array3D &v0,
160 const double &m0, const kep_toolbox::array3D &u, const double &s,
161 const double &mu = 1.0, const double &veff = 1.0, const double &c = 1.0,
162 const double &alpha = 1.5, const int &log10tolerance = -10,
163 const int &log10rtolerance = -10)
164 {
165 kep_toolbox::array3D r(r0), v(v0);
166 double m(m0);
167 double dt(0);
168 kep_toolbox::propagate_taylor_s(r, v, m, dt, u, s, mu, veff, c, alpha, log10tolerance, log10rtolerance);
169 return boost::python::make_tuple(kep_toolbox::array3D(r), kep_toolbox::array3D(v), double(m), double(dt));
170 }
171
fb_con_wrapper(const kep_toolbox::array3D & vin_rel,const kep_toolbox::array3D & vout_rel,const kep_toolbox::planet::base & pl)172 static inline tuple fb_con_wrapper(const kep_toolbox::array3D &vin_rel, const kep_toolbox::array3D &vout_rel,
173 const kep_toolbox::planet::base &pl)
174 {
175 double eq(0), ineq(0);
176 kep_toolbox::fb_con(eq, ineq, vin_rel, vout_rel, pl);
177 return boost::python::make_tuple(double(eq), double(ineq));
178 }
179
fb_prop_wrapper(const kep_toolbox::array3D & v_in,const kep_toolbox::array3D & v_pla,const double & rp,const double & beta,const double & mu)180 static inline kep_toolbox::array3D fb_prop_wrapper(const kep_toolbox::array3D &v_in, const kep_toolbox::array3D &v_pla,
181 const double &rp, const double &beta, const double &mu)
182 {
183 kep_toolbox::array3D retval;
184 kep_toolbox::fb_prop(retval, v_in, v_pla, rp, beta, mu);
185 return retval;
186 }
187
fb_vel_wrapper(const kep_toolbox::array3D & vin_rel,const kep_toolbox::array3D & vout_rel,const kep_toolbox::planet::base & pl)188 static inline double fb_vel_wrapper(const kep_toolbox::array3D &vin_rel, const kep_toolbox::array3D &vout_rel,
189 const kep_toolbox::planet::base &pl)
190 {
191 double dV(0);
192 kep_toolbox::fb_vel(dV, vin_rel, vout_rel, pl);
193 return dV;
194 }
195
damon_wrapper(const kep_toolbox::array3D & v1,const kep_toolbox::array3D & v2,double tof)196 static inline tuple damon_wrapper(const kep_toolbox::array3D &v1, const kep_toolbox::array3D &v2, double tof)
197 {
198 kep_toolbox::array3D a1, a2;
199 double tau, dv;
200 kep_toolbox::damon_approx(v1, v2, tof, a1, a2, tau, dv);
201 return boost::python::make_tuple(kep_toolbox::array3D(a1), kep_toolbox::array3D(a2), double(tau), double(dv));
202 }
203
204 #define get_constant(arg) \
205 static inline double get_##arg() \
206 { \
207 return ASTRO_##arg; \
208 }
209
210 get_constant(AU);
211 get_constant(JR);
212 get_constant(MU_SUN);
213 get_constant(MU_EARTH);
214 get_constant(EARTH_VELOCITY);
215 get_constant(EARTH_J2);
216 get_constant(EARTH_RADIUS);
217 get_constant(DEG2RAD);
218 get_constant(RAD2DEG);
219 get_constant(DAY2SEC);
220 get_constant(SEC2DAY);
221 get_constant(DAY2YEAR);
222 get_constant(G0);
223
224 #define PYKEP_REGISTER_CONVERTER(T, policy) \
225 { \
226 boost::python::type_info info = boost::python::type_id<T>(); \
227 const boost::python::converter::registration *reg = boost::python::converter::registry::query(info); \
228 if (reg == NULL) { \
229 to_tuple_mapping<T>(); \
230 from_python_sequence<T, policy>(); \
231 } else if ((*reg).m_to_python == NULL) { \
232 to_tuple_mapping<T>(); \
233 from_python_sequence<T, policy>(); \
234 } \
235 }
236
237 #define PYKEP_EXPOSE_CONSTANT(arg) def("_get_" #arg, &get_##arg);
238
239 typedef std::array<double, 8> array8D;
240 typedef std::array<double, 11> array11D;
241
BOOST_PYTHON_MODULE(core)242 BOOST_PYTHON_MODULE(core)
243 {
244 // Disable docstring c++ signature to allow sphinx autodoc to work properly
245 docstring_options doc_options;
246 doc_options.disable_signatures();
247
248 // Register std converters to python lists if not already registered by some
249 // other module
250 PYKEP_REGISTER_CONVERTER(std::vector<double>, variable_capacity_policy)
251 PYKEP_REGISTER_CONVERTER(kep_toolbox::array3D, fixed_size_policy)
252 PYKEP_REGISTER_CONVERTER(kep_toolbox::array6D, fixed_size_policy)
253 PYKEP_REGISTER_CONVERTER(kep_toolbox::array7D, fixed_size_policy)
254 PYKEP_REGISTER_CONVERTER(array8D, fixed_size_policy)
255 PYKEP_REGISTER_CONVERTER(array11D, fixed_size_policy)
256 PYKEP_REGISTER_CONVERTER(std::vector<kep_toolbox::array3D>, variable_capacity_policy)
257 PYKEP_REGISTER_CONVERTER(std::vector<array8D>, variable_capacity_policy)
258 PYKEP_REGISTER_CONVERTER(std::vector<array11D>, variable_capacity_policy)
259
260 // Expose the astrodynamical constants.
261 PYKEP_EXPOSE_CONSTANT(AU);
262 PYKEP_EXPOSE_CONSTANT(JR);
263 PYKEP_EXPOSE_CONSTANT(MU_SUN);
264 PYKEP_EXPOSE_CONSTANT(MU_EARTH);
265 PYKEP_EXPOSE_CONSTANT(EARTH_VELOCITY);
266 PYKEP_EXPOSE_CONSTANT(EARTH_J2);
267 PYKEP_EXPOSE_CONSTANT(EARTH_RADIUS);
268 PYKEP_EXPOSE_CONSTANT(DEG2RAD);
269 PYKEP_EXPOSE_CONSTANT(RAD2DEG);
270 PYKEP_EXPOSE_CONSTANT(DAY2SEC);
271 PYKEP_EXPOSE_CONSTANT(SEC2DAY);
272 PYKEP_EXPOSE_CONSTANT(DAY2YEAR);
273 PYKEP_EXPOSE_CONSTANT(G0);
274
275 // Exposing enums
276 enum_<kep_toolbox::epoch::type>("_epoch_type", "Defines a julian date type exposing the corresponding c++ enum\n\n"
277 "One of\n\n"
278 "* pykep.epoch.epoch_type.JD\n"
279 "* pykep.epoch.epoch_type.MJD\n"
280 "* pykep.epoch.epoch_type.MJD2000\n")
281 .value("MJD", kep_toolbox::epoch::MJD)
282 .value("MJD2000", kep_toolbox::epoch::MJD2000)
283 .value("JD", kep_toolbox::epoch::JD);
284
285 // Epoch class
286 class_<kep_toolbox::epoch>("epoch", pykep::epoch_doc().c_str(),
287 init<optional<const double &, kep_toolbox::epoch::type>>())
288 .add_property("jd", &kep_toolbox::epoch::jd,
289 "Returns the Julian Date\n\n"
290 "Example::\n\n"
291 " jd = e.jd")
292 .add_property("mjd", &kep_toolbox::epoch::mjd,
293 "Returns the Modified Julian Date\n\n"
294 "Example::\n\n"
295 " jd = e.mjd")
296 .add_property("mjd2000", &kep_toolbox::epoch::mjd2000,
297 "Returns the Modified Julian Date 2000\n\n"
298 "Example::\n\n"
299 " jd = e.mjd2000")
300 .def(repr(self))
301 .def_pickle(python_class_pickle_suite<kep_toolbox::epoch>());
302
303 // Epoch constructors helpers
304 def("epoch_from_string", &kep_toolbox::epoch_from_string, pykep::epoch_from_string_doc().c_str());
305 def("epoch_from_iso_string", &kep_toolbox::epoch_from_iso_string, pykep::epoch_from_iso_string_doc().c_str());
306
307 // Lambert.
308 class_<kep_toolbox::lambert_problem>(
309 "lambert_problem", "Represents a multiple revolution Lambert's problem",
310 init<const kep_toolbox::array3D &, const kep_toolbox::array3D &, const double &, const double &, const int &,
311 const int &>(pykep::lambert_problem_doc().c_str(),
312 (arg("r1") = kep_toolbox::array3D{1, 0, 0}, arg("r2") = kep_toolbox::array3D{0, 1, 0},
313 arg("tof") = boost::math::constants::pi<double>() / 2., arg("mu") = 1., arg("cw") = false,
314 arg("max_revs") = 0)))
315 .def("get_v1", &kep_toolbox::lambert_problem::get_v1, return_value_policy<copy_const_reference>(),
316 "Returns a sequence of vectors containing the velocities at r1 of "
317 "all computed solutions to the Lambert's Problem\n\n"
318 "Solutions are stored in order 0 rev, 1rev, 1rev, 2rev, 2rev, ...\n\n"
319 "Example (extracts v1 for the 0 revs solution)::\n\n"
320 " v10 = l.get_v1()[0]")
321 .def("get_v2", &kep_toolbox::lambert_problem::get_v2, return_value_policy<copy_const_reference>(),
322 "Returns a sequence of vectors containing the velocities at r2 of "
323 "all computed solutions to the Lambert's Problem\n\n"
324 "Solutions are stored in order 0 rev, 1rev, 1rev, 2rev, 2rev, ...\n\n"
325 "Example (extracts v2 for the 0 revs solution)::\n\n"
326 " v20 = l.get_v2()[0]")
327 .def("get_r1", &kep_toolbox::lambert_problem::get_r1, return_value_policy<copy_const_reference>(),
328 "Returns a vector containing the r1 defining the Lambert's "
329 "Problem\n\n"
330 "Example ::\n\n"
331 " r1 = l.get_r1()")
332 .def("get_r2", &kep_toolbox::lambert_problem::get_r2, return_value_policy<copy_const_reference>(),
333 "Returns a vector containing the r2 defining the Lambert's "
334 "Problem\n\n"
335 "Example ::\n\n"
336 " r2 = l.get_r2()")
337 .def("get_tof", &kep_toolbox::lambert_problem::get_tof, return_value_policy<copy_const_reference>(),
338 "Returns the time of flight defining the Lambert's Problem\n\n"
339 "Example::\n\n"
340 " t = l.get_tof()")
341 .def("get_mu", &kep_toolbox::lambert_problem::get_mu, return_value_policy<copy_const_reference>(),
342 "Returns the gravitational parameter defining the Lambert's "
343 "Problem\n\n"
344 "Example::\n\n"
345 " mu = l.get_mu()")
346 .def("get_x", &kep_toolbox::lambert_problem::get_x, return_value_policy<copy_const_reference>(),
347 "Returns a sequence containing the x values of all computed "
348 "solutions to the Lambert's Problem\n\n"
349 "Solutions are stored in order 0 rev, 1rev, 1rev, 2rev, 2rev, ...\n\n"
350 "Example (extracts x for the 0 revs solution)::\n\n"
351 " x0 = l.get_x()[0]")
352 .def("get_iters", &kep_toolbox::lambert_problem::get_iters, return_value_policy<copy_const_reference>(),
353 "Returns a sequence containing the number of iterations employed to "
354 "compute each solution to the Lambert's Problem\n\n"
355 "Solutions are stored in order 0 rev, 1rev, 1rev, 2rev, 2rev, ...\n\n"
356 "Example (extracts the number of iterations employed for the 0 revs "
357 "solution)::\n\n"
358 " p0 = l.get_iters()[0]")
359 .def("get_Nmax", &kep_toolbox::lambert_problem::get_Nmax,
360 "Returns the maximum number of revolutions allowed. The total "
361 "number of solution to the Lambert's problem will thus be n_sol = "
362 "Nmax*2 + 1\n\n"
363 "Example::\n\n"
364 " Nmax = l.get_Nmax()\n"
365 " n_sol = Nmax*2+1")
366 .def(repr(self))
367 .def_pickle(python_class_pickle_suite<kep_toolbox::lambert_problem>());
368
369 // Lagrangian propagator for keplerian orbits
370 def("propagate_lagrangian", &propagate_lagrangian_wrapper, pykep::propagate_lagrangian_doc().c_str(),
371 (arg("r0") = kep_toolbox::array3D{1, 0, 0}, arg("v0") = kep_toolbox::array3D{0, 1, 0},
372 arg("tof") = boost::math::constants::pi<double>() / 2, arg("mu") = 1));
373
374 // Taylor propagation of inertially constant thrust arcs
375 def("propagate_taylor", &propagate_taylor_wrapper, pykep::propagate_taylor_doc().c_str(),
376 (arg("r0") = kep_toolbox::array3D{1, 0, 0}, arg("v0") = arg("v0") = kep_toolbox::array3D{0, 1, 0},
377 arg("m0") = 100, arg("thrust") = kep_toolbox::array3D{0, 0, 0},
378 arg("tof") = boost::math::constants::pi<double>() / 2, arg("mu") = 1, arg("veff") = 1, arg("log10tol") = 1e-15,
379 arg("log10rtol") = 1e-15));
380
381 // Taylor propagation of inertially constant thrust arcs with an inertially constant disturbance
382 def("propagate_taylor_disturbance", &propagate_taylor_disturbance_wrapper,
383 "pykep.propagate_taylor(r,v,m,thrust,disturbance,t,mu,veff,log10tol,log10rtol)\n\n"
384 "- r: start position, x,y,z\n"
385 "- v: start velocity, vx,vy,vz\n"
386 "- m: starting mass\n"
387 "- thrust: fixed inertial thrust, cartesian components\n"
388 "- disturbance: fixed inertial disturbance force, cartesian components\n"
389 "- tof: propagation time\n"
390 "- mu: central body gravity constant\n\n"
391 "- veff: the product (Isp g0) defining the engine efficiency \n\n"
392 "- log10tol: the logarithm of the absolute tolerance passed to taylor propagator \n\n"
393 "- log10rtol: the logarithm of the relative tolerance passed to taylor propagator \n\n"
394 "Returns a tuple containing r, v, and m the final position, velocity and mass after the propagation.\n\n"
395 "Example::\n\n"
396 " r,v,m = propagate_taylor_disturbance([1,0,0],[0,1,0],100,[0,0,0],[1e-2,1e-3,1e-4],pi/2,1,1,-15,-15)",
397 (arg("r"), arg("v"), arg("m"), arg("thrust"), arg("disturbance"), arg("tof"), arg("mu"), arg("veff"),
398 arg("log10tol") = 1e-15, arg("log10rtol") = 1e-15));
399 // Taylor propagation of inertially constant thrust arcs under the J2
400 // influence
401 def("propagate_taylor_J2", &propagate_taylor_J2_wrapper,
402 "pykep.propagate_taylor_J2(r,v,m,thrust,t,mu,veff,J2RG2,log10tol,log10rtol)"
403 "\n\n"
404 "- r: start position, x,y,z\n"
405 "- v: start velocity, vx,vy,vz\n"
406 "- m: starting mass\n"
407 "- thrust: fixed inertial thrust, ux,uy,uz\n"
408 "- tof: propagation time\n"
409 "- mu: central body gravity constant\n\n"
410 "- veff: the product (Isp g0) defining the engine efficiency \n\n"
411 "- J2RG2: the product (J2 RE^2) defining the gravity J2 perturbance \n\n"
412 "- log10tol: the logarithm of the absolute tolerance passed to taylor "
413 "propagator \n\n"
414 "- log10rtol: the logarithm of the relative tolerance passed to taylor "
415 "propagator \n\n"
416 "Returns a tuple containing r, v, and m the final position, velocity and "
417 "mass after the propagation.\n\n"
418 "Example::\n\n"
419 " r,v,m = "
420 "propagate_taylor([1,0,0],[0,1,0],100,[0,0,0],pi/2,1,1,1e-3,-15,-15)",
421 (arg("r"), arg("v"), arg("m"), arg("thrust"), arg("tof"), arg("mu"), arg("veff"), arg("log10tol") = 1e-15,
422 arg("log10rtol") = 1e-15));
423 // Taylor propagation of inertially constant thrust arcs (using the
424 // generalized sundmann variable)
425 def("propagate_taylor_s", &propagate_taylor_s_wrapper,
426 "pykep.propagate_taylor_s(r,v,m,u,s,mu,veff,c=1,alpha,log10tol,log10rtol]"
427 ")\n\n"
428 "- r: start position, x,y,z\n"
429 "- v: start velocity, vx,vy,vz\n"
430 "- m: starting mass\n"
431 "- u: fixed inertial thrust, ux,uy,uz\n"
432 "- s: propagation pseudo-time\n"
433 "- mu: central body gravity constant\n\n"
434 "- veff: the product (Isp g0) defining the engine efficiency \n\n"
435 "- c: constant coefficient in the generalized Sundmann transform \n\n"
436 "- alpha: r exponent in the generalized Sundmann transform \n\n"
437 "- log10tol: the logarithm of the absolute tolerance passed to taylor "
438 "propagator \n\n"
439 "- log10rtol: the logarithm of the relative tolerance passed to taylor "
440 "propagator \n\n"
441 "Returns a tuple containing r, v, m and t the final position, velocity, "
442 "mass and time after the propagation pseudo-time s\n\n"
443 "Example::\n\n"
444 " r,v,m,t = "
445 "propagate_taylor_s([1,0,0],[0,1,0],100,[0,0,0],pi/"
446 "2,1.0,1.0,1.0,1.0,-10,-10)");
447
448 // Fly-by helper functions
449 def("fb_con", &fb_con_wrapper,
450 "pykep.fb_con(vin,vout,pl)\n\n"
451 "- vin: cartesian coordinates of the relative hyperbolic velocity before "
452 "the fly-by\n"
453 "- vout: vout, cartesian coordinates of the relative hyperbolic velocity "
454 "after the fly-by\n"
455 "- pl: fly-by planet\n\n"
456 "Returns a tuple containing (eq, ineq). \n"
457 " eq represents the violation of the equality constraint norm(vin)**2 = "
458 "norm(vout)**2.\n"
459 " ineq represents the violation of the inequality constraint on the "
460 "hyperbola asymptote maximum deflection\n\n"
461 "Example::\n\n"
462 " v2_eq, delta_ineq = fb_con(vin, vout, planet_ss('earth'))\n"
463 " v2_eq = v2_eq / EARTH_VELOCITY**2\n");
464 def("fb_prop", &fb_prop_wrapper,
465 "pykep.fb_prop(v,v_pla,rp,beta,mu)\n\n"
466 "- v: spacecarft velocity before the encounter (cartesian, absolute)\n"
467 "- v-pla: planet inertial velocity at encounter (cartesian, absolute)\n"
468 "- rp: fly-by radius\n"
469 "- beta: fly-by plane orientation\n"
470 "- mu: planet gravitational constant\n\n"
471 "Returns the cartesian components of the spacecarft velocity after the "
472 "encounter. \n"
473 "Example::\n\n"
474 " vout = fb_prop([1,0,0],[0,1,0],2,3.1415/2,1)\n");
475 def("fb_vel", &fb_vel_wrapper,
476 "pykep.fb_vel(vin,vout,pl)\n\n"
477 "- vin: cartesian coordinates of the relative hyperbolic velocity before "
478 "the fly-by\n"
479 "- vout: vout, cartesian coordinates of the relative hyperbolic velocity "
480 "after the fly-by\n"
481 "- pl: fly-by planet\n\n"
482 "Returns a number dV. \n"
483 " dV represents the norm of the delta-V needed to make a given fly-by "
484 "possible. dV is necessary to fix the magnitude and the direction of "
485 "vout.\n\n"
486 "Example::\n\n"
487 " dV = fb_vel(vin, vout, planet_ss('earth'))\n");
488
489 // Basic Astrodynamics
490 def("closest_distance", &closest_distance_wrapper,
491 "pykep.closest_distance(r1,v1,r2,v2,mu = 1.0)\n\n"
492 "- r1: initial position (cartesian)\n"
493 "- v1: initial velocity (cartesian)\n"
494 "- r2: final position (cartesian)\n"
495 "- v2: final velocity (cartesian)\n"
496 "- mu: planet gravitational constant\n\n"
497 "Returns a tuple containing the closest distance along a keplerian orbit "
498 "defined by r1,v1,r2,v2 (it assumes the orbit actually exists) and the "
499 "apoapsis radius of the resulting orbit. \n"
500 "Example::\n\n"
501 " d,ra = closest_distance([1,0,0],[0,1,0],[0,1,0],[-1,0,0],1.0)\n",
502 (arg("r1"), arg("v1"), arg("r2"), arg("v2"), arg("mu") = 1.0));
503
504 def("barker", &kep_toolbox::barker,
505 "pykep.barker(r1,r2,mu = 1.0)\n\n"
506 "- r1: initial position (cartesian)\n"
507 "- r2: final position (cartesian)\n"
508 "- mu: gravity parameter\n\n"
509 "Returns the time of flight as evaluated by the Barker equation. \n"
510 "Example:: \n\n"
511 " t = barker([1,0,0],[0,1,0],1.0)",
512 (arg("r1"), arg("r2"), arg("mu") = 1.0));
513
514 def("ic2par", &ic2par_wrapper,
515 "pykep.ic2par(r,v,mu = 1.0)\n\n"
516 "- r: position (cartesian)\n"
517 "- v: velocity (cartesian)\n"
518 "- mu: gravity parameter\n\n"
519 "Returns the osculating keplerian elements a,e,i,W,w,E\n"
520 "E is the eccentric anomaly for e<1, the Gudermannian for e>1\n"
521 "a is the semi-major axis always a positive quantity.\n"
522 "NOTE: The routine is singular when the elements are not defined.\n"
523 "Example:: \n\n"
524 " el = ic2par([1,0,0],[0,1,0],1.0)",
525 (arg("r"), arg("v"), arg("mu") = 1.0));
526
527 def("par2ic", &par2ic_wrapper,
528 "pykep.par2ic(E,mu = 1.0)\n\n"
529 "- E: osculating keplerian elements a,e,i,W,w,E ( l, ND, rad, rad, rad, rad)\n"
530 "- mu: gravity parameter (l^3/s^2)\n\n"
531 "Returns cartesian elements from Keplerian elements\n"
532 "E is the eccentric anomaly for e<1, the Gudermannian for e>1\n"
533 "a is the semi-major axis always a positive quantity.\n"
534 "Example:: \n\n"
535 " r, v = pk.par2ic([1,0.3,0.1,0.1,0.2,0.2], 1)",
536 (arg("E"), arg("mu") = 1.0));
537
538 def("ic2eq", &ic2eq_wrapper,
539 "pykep.ic2eq(r,v,mu = 1.0, retrogade = False)\n\n"
540 "- r: position (cartesian)\n"
541 "- v: velocity (cartesian)\n"
542 "- mu: gravity parameter\n\n"
543 "- retrogade: uses the retrograde parameters. Default value is False.\n\n"
544 "Returns the modified equinoctial elements a(1-e^2),f,g,h,k,L\n"
545 "L is the true mean longitude\n"
546 "Example:: \n\n"
547 " E = ic2eq(r = [1,0,0], v = [0,1,0], mu =1.0)",
548 (arg("r"), arg("v"), arg("mu") = 1.0, arg("retrograde") = false));
549
550 def("eq2ic", &eq2ic_wrapper,
551 "pykep.eq2ic(EQ,mu = 1.0, retrograde = False)\n\n"
552 "- EQ: modified equinoctial elements a(1-e^2),f,g,h,k,L\n"
553 "- mu: gravity parameter (l^3/s^2)\n\n"
554 "Returns cartesian elements from Keplerian elements\n"
555 "L is the true longitude\n"
556 "Example:: \n\n"
557 " r, v = pk.eq2ic(eq = [1,0.3,0.1,0.1,0.2,0.2], mu = 1, retrograde = False)",
558 (arg("eq"), arg("mu") = 1.0, arg("retrograde") = false));
559
560 def("eq2par", &eq2par_wrapper,
561 "pykep.eq2par(EQ, retrograde = False)\n\n"
562 "- EQ: modified equinoctial elements a(1-e^2),f,g,h,k,L\n"
563 "- retrogade: uses the retrograde parameters. Default value is False.\n\n"
564 "Returns the osculating Keplerian elements a,e,i,W,w,E \n"
565 "E is the eccentric anomaly for e<1, the Gudermannian for e>1\n"
566 "a is the semi-major axis, always a positive quantity.\n"
567 "L is the true longitude\n"
568 "Example:: \n\n"
569 " E = pk.eq2par(eq = [1,0.3,0.1,0.1,0.2,0.2], retrograde = False)",
570 (arg("eq"), arg("retrograde") = false));
571
572 def("par2eq", &par2eq_wrapper,
573 "pykep.par2eq(E, retrograde = False)\n\n"
574 "- E: osculating keplerian elements a,e,i,W,w,E ( l, ND, rad, rad, rad, rad)\n"
575 "- retrogade: uses the retrograde parameters. Default value is False.\n\n"
576 "Returns the modified equinoctial elements a(1-e^2),f,g,h,k,L\n"
577 "L is the true mean longitude\n"
578 "Example:: \n\n"
579 " E = pk.par2eq(E = [1.0,0.1,0.2,0.3,0.4,0.5], retrograde = False)",
580 (arg("E"), arg("retrograde") = false));
581
582 def("damon", &damon_wrapper,
583 "pykep.damon(v1,v2,tof)\n\n"
584 "- v1: starting velocity relative to the departure body. This is\n"
585 " computed from the Lambert solution linking initial and\n"
586 " final position in tof\n"
587 "- v2: ending velocity relative to the arrival body. This is\n"
588 " computed from the Lambert solution linking initial and\n"
589 " final position in tof\n"
590 "- tof: time of flight\n\n"
591 "Returns:\n\n"
592 "- a1: acceleration vector in the first part of the transfer\n"
593 "- a2: acceleration vector in the second part of the transfer\n"
594 "- tau: duration of the first part of the transfer\n"
595 "- dv: total estimated DV\n\n"
596 "This function uses the model developed by Damon Landau (JPL)\n"
597 "during GTOC7 to compute an approximation to the low-thrust transfer\n\n"
598 "Example:: \n\n"
599 " a1, a2, tau, dv = pykep.damon(v1,v2,tof)",
600 (arg("v1"), arg("v2"), arg("tof")));
601
602 def("max_start_mass", &kep_toolbox::max_start_mass,
603 "pykep.max_start_mass(a, dv, T_max, Isp)\n\n"
604 "- a: acceleration magnitude as computed from Damon's model\n"
605 "- dv: dv magnitude as computed from Damon's model\n"
606 "- T_max: maximum thrust of the spacecarft NEP propulsion system\n"
607 "- Isp: specific impulse of the spacecarft NEP propulsion system\n\n"
608 "Returns the estimated maximum mass at leg beginning\n"
609 "This function estimates the maximum\n"
610 "mass a NEP spacececraft can have in order for a given\n"
611 "arc to be convertable into a valid low-thrust transfer\n\n"
612 "Example:: \n\n"
613 " a,_,_,dv = pykep.damon(v1,v2,tof)\n"
614 " m = pykep.max_start_mass(norm(a), dv, T_max, Isp)");
615
616 // For the three_impulses_approximation we need to distinguish the overloaded
617 // cases. We thus introduce new function pointers
618 double (*three_imp)(double, double, double, double, double, double, double, double, double)
619 = &kep_toolbox::three_impulses_approx;
620 double (*three_imp_2_pl)(const kep_toolbox::planet::base &, const kep_toolbox::planet::base &)
621 = &kep_toolbox::three_impulses_approx;
622 double (*three_imp_2_pl_2_ep)(const kep_toolbox::planet::base &, const kep_toolbox::planet::base &,
623 const kep_toolbox::epoch &, const kep_toolbox::epoch &)
624 = &kep_toolbox::three_impulses_approx;
625
626 def("_three_impulses_approx", three_imp);
627 def("_three_impulses_approx", three_imp_2_pl);
628 def("_three_impulses_approx", three_imp_2_pl_2_ep);
629 }
630