Structure for the python module

.github/workflows/main.yml

@@ -1,4 +1,4 @@
-name: GitHub CI
+name: C++ Library
 on:
   push:
     branches:

.github/workflows/pykep.yml

@@ -0,0 +1,17 @@
+name: Python Module
+on:
+  push:
+    branches:
+      - main
+    tags:
+      - 'v[0-9]+.[0-9]+.[0-9]+'
+  pull_request:
+    branches:
+      - main
+jobs:
+  linux_pykep:
+    runs-on: ubuntu-latest
+    steps:
+      - uses: actions/checkout@v3
+      - name: Build
+        run: bash tools/gha_linux_pykep.sh

.gitignore

@@ -9,3 +9,6 @@ build*
 
 # config files
 include/kep3/config.hpp
+
+# python stuff
+*__pycache__*

CMakeLists.txt

@@ -165,9 +165,12 @@ target_include_directories(kep3 PUBLIC
     $<BUILD_INTERFACE:${CMAKE_CURRENT_BINARY_DIR}/include>
     $<INSTALL_INTERFACE:include>)
 # Boost.
-# NOTE: need 1.73 for atomic_ref.
-find_package(Boost 1.73 REQUIRED)
+# NOTE: need 1.69 for safe numerics.
+set(_kep3_MIN_BOOST_VERSION "1.69")
+find_package(Boost ${_kep3_MIN_BOOST_VERSION} REQUIRED)
 target_link_libraries(kep3 PUBLIC Boost::boost)
+# NOTE: quench warnings from Boost when building the library.
+target_compile_definitions(kep3 PRIVATE BOOST_ALLOW_DEPRECATED_HEADERS)
 
 # fmt.
 find_package(fmt CONFIG REQUIRED)
@@ -207,6 +210,17 @@ install(TARGETS kep3
     ARCHIVE DESTINATION "${kep3_INSTALL_LIBDIR}"
     RUNTIME DESTINATION bin
 )
+# Setup of the CMake config file.
+configure_file("${CMAKE_CURRENT_SOURCE_DIR}/kep3-config.cmake.in" "${CMAKE_CURRENT_BINARY_DIR}/kep3-config.cmake" @ONLY)
+install(FILES "${CMAKE_CURRENT_BINARY_DIR}/kep3-config.cmake"
+    DESTINATION "${kep3_INSTALL_LIBDIR}/cmake/kep3")
+install(EXPORT kep3_export NAMESPACE kep3:: DESTINATION "${kep3_INSTALL_LIBDIR}/cmake/kep3")
+# Take care of versioning.
+include(CMakePackageConfigHelpers)
+write_basic_package_version_file("${CMAKE_CURRENT_BINARY_DIR}/kep3-config-version.cmake" COMPATIBILITY SameMinorVersion)
+install(FILES "${CMAKE_CURRENT_BINARY_DIR}/kep3-config-version.cmake" DESTINATION "${kep3_INSTALL_LIBDIR}/cmake/kep3")
+
+
 if(kep3_BUILD_TESTS)
     enable_testing()
     add_subdirectory(test)
@@ -215,3 +229,21 @@ endif()
 if(kep3_BUILD_BENCHMARKS)
     add_subdirectory(benchmark)
 endif()
+
+if(kep3_BUILD_PYTHON_BINDINGS)
+    # Find python
+    find_package(Python3 REQUIRED COMPONENTS Interpreter Development)
+    message(STATUS "Python3 interpreter: ${Python3_EXECUTABLE}")
+    message(STATUS "Python3 installation directory: ${Python3_SITEARCH}")
+    set(PYKEP_INSTALL_PATH "" CACHE STRING "pykep module installation path")
+    mark_as_advanced(PYKEP_INSTALL_PATH)
+
+    # pybind11.
+    find_package(pybind11 REQUIRED)
+
+    # Finding kep3 (cpp) 
+    find_package(kep3 ${kep3_VERSION} EXACT REQUIRED)
+
+    # Build directory
+    add_subdirectory("${CMAKE_SOURCE_DIR}/pykep")
+endif()

include/kep3/core_astro/constants.hpp

@@ -28,7 +28,7 @@ constexpr double pi = boost::math::constants::pi<double>();
 constexpr double half_pi = boost::math::constants::half_pi<double>();
 
 constexpr double AU = 149597870700.0;     // Astronomical Unit (m)
-constexpr double CAVENDISH = 73.6687e-11; // Cavendish constant
+constexpr double CAVENDISH = 73.6687e-11; // Cavendish constant (N M^2 / kg^2)
 constexpr double MU_SUN =
     1.32712440018e20; // Sun's gravitational parameter (m^3/s^2 kg)
 constexpr double MU_EARTH =

include/kep3/core_astro/convert_anomalies.hpp

@@ -10,12 +10,13 @@
 #ifndef kep3_CONVERT_ANOMALIES_H
 #define kep3_CONVERT_ANOMALIES_H
 
-#include "kep3/core_astro/constants.hpp"
 #include <cmath>
+#include <stdexcept>
 
 #include <boost/math/constants/constants.hpp>
 #include <boost/math/tools/roots.hpp>
 
+#include <kep3/core_astro/constants.hpp>
 #include <kep3/core_astro/kepler_equations.hpp>
 
 namespace kep3 {
@@ -24,16 +25,18 @@ namespace kep3 {
 // of revolutions.
 // NOLINTNEXTLINE(bugprone-easily-swappable-parameters)
 inline double m2e(double M, double ecc) {
+  if (ecc >= 1) {
+    throw std::domain_error(
+        "m2e: Eccentric anomaly is not defined for ecc >=1.");
+  }
   // We compute the sin and cos of the mean anomaly which are used also later
   // for the initial guess (3rd order expansion of Kepler's equation).
   double sinM = std::sin(M), cosM = std::cos(M);
   // Here we use the atan2 to recover the mean anomaly in the [0,2pi] range.
   // This makes sure that for high value of M no catastrophic cancellation
   // occurs, as would be the case using std::fmod(M, 2pi)
   double M_cropped = std::atan2(sinM, cosM);
-  if (M_cropped < 0) {
-    M_cropped += 2 * kep3::pi;
-  }
+
   // The Initial guess follows from a third order expansion of Kepler's
   // equation.
   double IG = M_cropped + ecc * sinM + ecc * ecc * sinM * cosM +
@@ -57,38 +60,73 @@ inline double m2e(double M, double ecc) {
   return sol;
 }
 // eccentric to mean (only ellipses) e<1
-inline double e2m(double E, double ecc) { return (E - ecc * std::sin(E)); }
+inline double e2m(double E, double ecc) {
+  if (ecc >= 1) {
+    throw std::domain_error(
+        "e2m: Eccentric anomaly is not defined for ecc >=1.");
+  }
+  return (E - ecc * std::sin(E));
+}
 
 // eccentric to true (only ellipses) e<1 (returns in range [-pi,pi])
 inline double e2f(double E, double ecc) {
+  if (ecc >= 1) {
+    throw std::domain_error(
+        "e2f: Eccentric anomaly is not defined for ecc >=1.");
+  }
   return 2 * std::atan(std::sqrt((1 + ecc) / (1 - ecc)) * std::tan(E / 2));
 }
 
 // true to eccentric (only ellipses) e<1 (returns in range [-pi,pi])
 inline double f2e(double f, double ecc) {
+  if (ecc >= 1) {
+    throw std::domain_error(
+        "f2e: Eccentric anomaly is not defined for ecc >=1.");
+  }
   return 2 * std::atan(std::sqrt((1 - ecc) / (1 + ecc)) * std::tan(f / 2));
 }
 
 // mean to true (only ellipses) e<1 (returns in range [-pi,pi])
-inline double m2f(double M, double ecc) { return e2f(m2e(M, ecc), ecc); }
+inline double m2f(double M, double ecc) {
+  if (ecc >= 1) {
+    throw std::domain_error("m2f: Mean anomaly is not defined for ecc >=1.");
+  };
+  return e2f(m2e(M, ecc), ecc);
+}
 
 // true to mean (only ellipses) e<1 (returns in range [-pi,pi])
-inline double f2m(double f, double ecc) { return e2m(f2e(f, ecc), ecc); }
+inline double f2m(double f, double ecc) {
+  if (ecc >= 1) {
+    throw std::domain_error("f2m: Mean anomaly is not defined for ecc >=1.");
+  };
+  return e2m(f2e(f, ecc), ecc);
+}
 
 // gudermannian to true (only hyperbolas) e>1 (returns in range [-pi,pi])
 inline double zeta2f(double f, double ecc) {
+  if (ecc <= 1) {
+    throw std::domain_error(
+        "zeta2f: zeta (Gudermannian) is not defined for ecc <=1.");
+  };
   return 2 * std::atan(std::sqrt((1 + ecc) / (ecc - 1)) * std::tan(f / 2));
 }
 
 // true to gudermannian (only hyperbolas) e>1 (returns in range [-pi,pi])
 inline double f2zeta(double zeta, double ecc) {
+  if (ecc <= 1) {
+    throw std::domain_error(
+        "f2zeta: zeta (Gudermannian) is not defined for ecc <=1.");
+  };
   return 2 * std::atan(std::sqrt((ecc - 1) / (1 + ecc)) * std::tan(zeta / 2));
 }
 
 // mean to hyperbolic (only hyperbolas) e>1.
 // NOLINTNEXTLINE(bugprone-easily-swappable-parameters)
 inline double n2h(double N, double ecc) {
-
+  if (ecc <= 1) {
+    throw std::domain_error(
+        "n2h: hyperbolic anomaly is not defined for ecc <=1.");
+  };
   // The Initial guess (TODO(darioizo) improve)
   double IG = 1.;
 
@@ -110,23 +148,49 @@ inline double n2h(double N, double ecc) {
 }
 
 // hyperbolic H to hyperbolic mean N (only hyperbolas) e>1
-inline double h2n(double H, double ecc) { return (ecc * std::sinh(H) - H); }
+inline double h2n(double H, double ecc) {
+  if (ecc <= 1) {
+    throw std::domain_error(
+        "h2n: hyperbolic anomaly is not defined for ecc <=1.");
+  };
+  return (ecc * std::sinh(H) - H);
+}
 
 // hyperbolic H to true (only hyperbolas) e>1 (returns in range [-pi,pi])
 inline double h2f(double H, double ecc) {
+  if (ecc <= 1) {
+    throw std::domain_error(
+        "h2f: hyperbolic anomaly is not defined for ecc <=1.");
+  };
   return 2 * std::atan(std::sqrt((1 + ecc) / (ecc - 1)) * std::tanh(H / 2));
 }
 
 // true to hyperbolic H (only hyperbolas) e>1 (returns in range [-pi,pi])
 inline double f2h(double f, double ecc) {
+  if (ecc <= 1) {
+    throw std::domain_error(
+        "f2h: hyperbolic anomaly is not defined for ecc <=1.");
+  };
   return 2 * std::atanh(std::sqrt((ecc - 1) / (1 + ecc)) * std::tan(f / 2));
 }
 
 // mean hyperbolic to true (only hyperbolas) e>1 (returns in range [-pi,pi])
-inline double n2f(double N, double ecc) { return h2f(n2h(N, ecc), ecc); }
+inline double n2f(double N, double ecc) {
+  if (ecc <= 1) {
+    throw std::domain_error(
+        "n2f: mean hyperbolic anomaly is not defined for ecc <=1.");
+  };
+  return h2f(n2h(N, ecc), ecc);
+}
 
 // true to mean hyperbolic (only hyperbolas) e>1 (returns in range [-pi,pi])
-inline double f2n(double f, double ecc) { return h2n(f2h(f, ecc), ecc); }
+inline double f2n(double f, double ecc) {
+  if (ecc <= 1) {
+    throw std::domain_error(
+        "f2n: mean hyperbolic anomaly is not defined for ecc <=1.");
+  };
+  return h2n(f2h(f, ecc), ecc);
+}
 
 } // namespace kep3
 #endif // kep3_TOOLBOX_M2E_H

include/kep3/planets/jpl_lp.hpp

@@ -54,7 +54,7 @@ class kep3_DLL_PUBLIC jpl_lp {
 
 public:
   // Constructor
-  explicit jpl_lp(const std::string & = "earth");
+  explicit jpl_lp(std::string = "earth");
   // Mandatory UDPLA methods
   [[nodiscard]] std::array<std::array<double, 3>, 2> eph(const epoch &) const;
 

kep3-config.cmake.in

@@ -0,0 +1,14 @@
+# Mandatory public dependencies on Boost and fmt.
+find_package(Boost @_kep3_MIN_BOOST_VERSION@ REQUIRED serialization)
+find_package(fmt REQUIRED CONFIG)
+find_package(heyoka REQUIRED CONFIG)
+find_package(xtensor REQUIRED CONFIG)
+find_package(xtensor-blas REQUIRED CONFIG)
+
+# Get current dir.
+get_filename_component(_kep3_CONFIG_SELF_DIR "${CMAKE_CURRENT_LIST_FILE}" PATH)
+
+include(${_kep3_CONFIG_SELF_DIR}/kep3_export.cmake)
+
+# Clean up.
+unset(_kep3_CONFIG_SELF_DIR)

kep3_devel.yml

@@ -12,4 +12,5 @@ dependencies:
   - spdlog
   - xtensor
   - xtensor-blas
+  - pybind11
 

pykep/CMakeLists.txt

@@ -0,0 +1,61 @@
+# Configure the version file.
+configure_file("${CMAKE_CURRENT_SOURCE_DIR}/_version.py.in" "${CMAKE_CURRENT_SOURCE_DIR}/_version.py" @ONLY)   
+
+# The list of pydsyre's Python files.
+set(kep3_PYTHON_FILES __init__.py _version.py test.py)
+
+# Copy the python files in the current binary dir,
+# so that we can import pydesyre from the build dir.
+# NOTE: importing from the build dir will work
+# only on single-configuration generators.
+foreach(kep3_PYTHON_FILE ${kep3_PYTHON_FILES})
+    configure_file("${CMAKE_CURRENT_SOURCE_DIR}/${kep3_PYTHON_FILE}"
+        "${CMAKE_CURRENT_BINARY_DIR}/${kep3_PYTHON_FILE}" COPYONLY)
+endforeach()
+
+# Core module.
+Python3_add_library(core MODULE WITH_SOABI
+    core.cpp
+    docstrings.cpp
+)
+
+target_link_libraries(core PRIVATE kep3::kep3)
+target_link_libraries(core PRIVATE "${pybind11_LIBRARIES}")
+
+target_include_directories(core SYSTEM PRIVATE "${pybind11_INCLUDE_DIR}" "${Python3_INCLUDE_DIRS}")
+target_compile_definitions(core PRIVATE "${pybind11_DEFINITIONS}")
+target_compile_options(core PRIVATE
+    "$<$<CONFIG:Debug>:${KEP3_CXX_FLAGS_DEBUG}>"
+    "$<$<CONFIG:Release>:${KEP3_CXX_FLAGS_RELEASE}>"
+    "$<$<CONFIG:RelWithDebInfo>:${KEP3_CXX_FLAGS_RELEASE}>"
+    "$<$<CONFIG:MinSizeRel>:${KEP3_CXX_FLAGS_RELEASE}>"
+
+)
+target_include_directories(core PUBLIC
+    $<BUILD_INTERFACE:${CMAKE_CURRENT_SOURCE_DIR}/include>
+    $<BUILD_INTERFACE:${CMAKE_CURRENT_BINARY_DIR}/include>
+    $<INSTALL_INTERFACE:include>)
+set_target_properties(core PROPERTIES CXX_VISIBILITY_PRESET hidden)
+set_target_properties(core PROPERTIES VISIBILITY_INLINES_HIDDEN TRUE)
+target_compile_features(core PRIVATE cxx_std_17)
+set_property(TARGET core PROPERTY CXX_EXTENSIONS NO)
+
+# Installation setup.
+if(PYKEP_INSTALL_PATH STREQUAL "")
+    message(STATUS "pykep will be installed in the default location: ${Python3_SITEARCH}")
+    set(_PYKEP_INSTALL_DIR "${Python3_SITEARCH}/pykep")
+else()
+    message(STATUS "pykep will be installed in the custom location: ${PYKEP_INSTALL_PATH}")
+    set(_PYKEP_INSTALL_DIR "${PYKEP_INSTALL_PATH}/pykep")
+endif()
+
+# Install the core module.
+install(TARGETS core
+    RUNTIME DESTINATION ${_PYKEP_INSTALL_DIR}
+    LIBRARY DESTINATION ${_PYKEP_INSTALL_DIR}
+)
+
+# Install the Python files.
+install(FILES ${PYKEP_PYTHON_FILES} DESTINATION ${_PYKEP_INSTALL_DIR})
+
+unset(_PYKEP_INSTALL_DIR)

pykep/__init__.py

@@ -0,0 +1,18 @@
+"""Pykep is a coolbox for interplanetary trajectory design developed by ESA's Advanced Concepts Team. Its main
+ purpose is fast prototyping of reseacrh ideas, and is not intended for operational usage.
+
+ Some important conventions followed:
+ 1 - All units expected are S.I. (m,sec,kg,N) unless explicitly stated.
+ 2 - The default set of osculating orbital parameters is, in this order: [sma, ecc, incl, W, w, f], where f is the true anomaly
+ 3 - The default option to represent epochs as floats is the modified julian date 2000 (MJD2000)."""
+
+# Version setup.
+from ._version import __version__
+del _version
+
+
+# Importing cpp functionalities
+from .core import *
+
+# We import the unit test submodule
+from . import test

pykep/_version.py

@@ -0,0 +1 @@
+__version__ = "0.0.1"

pykep/_version.py.in

@@ -0,0 +1 @@
+__version__ = "@kep3_VERSION@"