Merge pull request #6 from cantordust/chrono

Chrono migration

CMakeLists.txt

@@ -228,6 +228,7 @@ endif()
 
 if(kep3_BUILD_BENCHMARKS)
     add_subdirectory(benchmark)
+
 endif()
 
 if(kep3_BUILD_PYTHON_BINDINGS)
@@ -247,3 +248,4 @@ if(kep3_BUILD_PYTHON_BINDINGS)
     # Build directory
     add_subdirectory("${CMAKE_SOURCE_DIR}/pykep")
 endif()
+

include/kep3/core_astro/convert_julian_dates.hpp

@@ -7,16 +7,21 @@
 // Public License v. 2.0. If a copy of the MPL was not distributed
 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
 
-#ifndef kep3_CONVERT_JULIAN_DATES_H
-#define kep3_CONVERT_JULIAN_DATES_H
+#ifndef kep3_CONVERT_JULIAN_DATES_HPP
+#define kep3_CONVERT_JULIAN_DATES_HPP
+
+#include "kep3/epoch.hpp"
+#include <chrono>
+namespace kep3
+{
+    inline double jd2mjd(double in) { return (in - 2400000.5); }
+    inline double jd2mjd2000(double in) { return (in - 2451544.5); }
+    inline double mjd2jd(double in) { return (in + 2400000.5); }
+    inline double mjd2mjd2000(double in) { return (in - 51544); }
+    inline double mjd20002jd(double in) { return (in + 2451544.5); }
+    inline double mjd20002mjd(double in){ return (in + 51544); }
+
 
-namespace kep3 {
-inline double jd2mjd(double in) { return (in - 2400000.5); }
-inline double jd2mjd2000(double in) { return (in - 2451544.5); }
-inline double mjd2jd(double in) { return (in + 2400000.5); }
-inline double mjd2mjd2000(double in) { return (in - 51544); }
-inline double mjd20002jd(double in) { return (in + 2451544.5); }
-inline double mjd20002mjd(double in) { return (in + 51544); }
 } // namespace kep3
 
-#endif // kep3_CONVERT_JULIAN_DATES_H
+#endif // kep3_CONVERT_JULIAN_DATES_HPP

include/kep3/detail/s11n.hpp

@@ -21,6 +21,7 @@
 #include <boost/serialization/base_object.hpp>
 #include <boost/serialization/export.hpp>
 #include <boost/serialization/serialization.hpp>
+#include <boost/serialization/binary_object.hpp>
 #include <boost/serialization/split_free.hpp>
 #include <boost/serialization/split_member.hpp>
 #include <boost/serialization/string.hpp>

include/kep3/epoch.hpp

@@ -7,100 +7,240 @@
 // Public License v. 2.0. If a copy of the MPL was not distributed
 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
 
-#ifndef kep3_EPOCH_H
-#define kep3_EPOCH_H
-
-#include <iostream>
+#ifndef kep3_EPOCH_HPP
+#define kep3_EPOCH_HPP
 
+#include <chrono>
+#include <cstdint>
+#include <ctime>
 #include <fmt/ostream.h>
-
-#include <boost/date_time/gregorian/gregorian.hpp>
-#include <boost/date_time/posix_time/posix_time.hpp>
-#include <boost/lexical_cast.hpp>
+#include <iostream>
 
 #include <kep3/detail/s11n.hpp>
 #include <kep3/detail/visibility.hpp>
+#include <ratio>
+#include <type_traits>
+#include <utility>
 
 /// Keplerian Toolbox
 /**
  * This namespace contains astrodynamics and space flight mechanics routines
  * that are related to keplerian motions or models.
  */
+
 namespace kep3 {
+using namespace std::literals;
+namespace chr = std::chrono;
+using uint = unsigned int;
+
+template <typename T> struct is_duration : std::false_type {};
+
+template <typename Rep, typename Period>
+struct is_duration<chr::duration<Rep, Period>> : std::true_type {};
+
+template <typename T>
+using enable_if_is_duration = std::enable_if_t<is_duration<T>::value>;
+
+struct kep_clock : public chr::system_clock {
+
+  /**
+   * @brief Custom clock.
+   * Used for constructing epochs with a custom reference point (1 Jan 2000).
+   * By defining a custom clock, we avoid the overflow
+   * that std::chrono::system_clock suffers at +/- 292 years.
+   * To do that, we lower the resolution to 1 us (microsecond),
+   * which gives the clock a range of +/- 292 thousand years.
+   *
+   * NOTE: Adding durations of less than 1 us to an epoch (defined below)
+   * would not be registered.
+   *
+   * NOTE: As of C++20, the standard guarantees that std::chrono::system_clock
+   * uses the UNIX time reference point, which is midnight on 1 January 1970
+   * (1970-01-01T00:00:00). We correct for that here in order to bring the
+   * reference point forward to midnight on 1 January 2000
+   * (2000-01-01T00:00:00), which is 0 MJD2000.
+   */
+  using rep = int_fast64_t;
+  // Resolution of (1 / 1'000'000)s = 1 us
+  using period = std::ratio<1, 1'000'000>;
+  using duration = chr::duration<rep, period>;
+  using time_point = chr::time_point<kep_clock, duration>;
+  static constexpr bool is_steady = false;
+  // Number of seconds from midnight on 1 Jan 1970 to midnight on 1 Jan 2000.
+  static constexpr chr::seconds y2k_offset{946684800s};
+
+  static constexpr time_point ref_epoch{kep_clock::time_point{} + y2k_offset};
+
+  static constexpr std::time_t to_time_t(const time_point &t) noexcept {
+    return static_cast<std::time_t>(
+        chr::duration_cast<chr::seconds>(t.time_since_epoch() + y2k_offset)
+            .count());
+  }
+
+  static constexpr time_point from_time_t(std::time_t t) noexcept {
+    return chr::time_point_cast<duration>(
+        time_point(chr::seconds(t) - y2k_offset));
+  }
+};
 
 /// epoch class.
 /**
- * This class defines and contains a non-gregorian date (i.e. a date expressed
- * in julian form). It also provides the user with an interface to boost
- * gregorian dates (see boost documentation at
- * http://www.boost.org/doc/libs/1_42_0/doc/html/date_time.html)
- * using the posix time.
- * The date is defined in MJD2000 (double) as a
- * private member
- *
- * @author Dario Izzo (dario.izzo _AT_ googlemail.com)
+ * This class defines and contains a non-Gregorian date (i.e., a date expressed
+ * in Julian format). The date is defined in MJD2000 format as a
+ * kep_clock::time_point private member. Types of non-Gregorian dates supported:
+ *      - Julian Date (JD): the number of days passed since January 1, 4713 BC
+ * at 12:00 (noon).
+ *      - Modified Julian Date (MJD): the number of days passed since November
+ * 17, 1858 at 00:00 (midnight).
+ *      - Modified Julian Date 2000 (MJD2000): the number of days passed since
+ * Juanuary 1, 2000 at 00:00 (midnight).
  */
 class kep3_DLL_PUBLIC epoch {
 public:
-  /** Types of non gregorian dates supported. Julian Date (JD) is the number of
-   * days passed since January 1, 4713 BC at noon. Modified Julian Date (MJD) is
-   * the number of days passed since November 17, 1858 at 00:00 am. The Modified
-   * Julian Date 2000 (MJD2000) is the number of days passed since Juanuary 1,
-   * 2000 at 00:00am.
-   */
-  enum julian_type { MJD2000, MJD, JD };
+  enum class julian_type { MJD2000, MJD, JD };
 
   /** Constructors */
-  explicit epoch(const double & = 0., julian_type = MJD2000);
-  epoch(const boost::gregorian::greg_day &day,
-        const boost::gregorian::greg_month &month,
-        const boost::gregorian::greg_year &year);
-  explicit epoch(const boost::posix_time::ptime &posix_time);
-
-  /** Computing non-gregorian dates */
-  [[nodiscard]] double mjd2000() const;
-  [[nodiscard]] double jd() const;
-  [[nodiscard]] double mjd() const;
-
-  /** Interface to boost::posix_time::ptime */
-  [[nodiscard]] boost::posix_time::ptime get_posix_time() const;
-  void set_posix_time(const boost::posix_time::ptime &);
-
-  /** operators overloads for sum diff (epoch-days) and the comparison operators
+  // Default constructor
+  epoch();
+
+  // Constructor for days (as a floating-point value)
+  explicit epoch(double epoch_in,
+                 julian_type epoch_type = julian_type::MJD2000);
+
+  /**
+   * Constructs an epoch from a std::chrono::duration.
+   * The reference point is assumed to be MJD2000.
+   * \param[in] time The time as a duration.
+   */
+  template <class Duration, class = enable_if_is_duration<Duration>>
+  explicit epoch(const Duration &duration)
+      : tp{kep_clock::time_point{} + duration} {}
+
+  // Constructor from duration&&)
+  template <class Duration, class = enable_if_is_duration<Duration>>
+  explicit epoch(Duration &&duration)
+      : tp{kep_clock::time_point{} + std::forward(duration)} {}
+
+
+  // Constructor for datetime broken down into its constituents.
+  explicit epoch(const int y, const uint mon, const uint d, const int h = 0,
+                 const int min = 0, const int s = 0, const int ms = 0, const int us = 0);
+
+  /* Computing non-Gregorian dates */
+
+  /**
+   * @return Number of days since 0 JD (including fractional days).
+   */
+  [[nodiscard]] constexpr double jd() const {
+    return chr::duration<double, std::ratio<86400>>(tp.time_since_epoch() - kep_clock::y2k_offset + 211813444800s).count();
+  }
+
+  /**
+   * @return Number of days since 0 MJD (including fractional days).
+   */
+  [[nodiscard]] constexpr double mjd() const {
+    return chr::duration<double, std::ratio<86400>>(tp.time_since_epoch() - kep_clock::y2k_offset + 4453401600s).count();
+  }
+
+  /**
+   * @return Number of days since 0 MJD2000 (including fractional days).
+   */
+  [[nodiscard]] constexpr double mjd2000() const {
+    return chr::duration<double, std::ratio<86400>>(tp.time_since_epoch() - kep_clock::y2k_offset).count();
+  }
+
+  /* Helper functions for constructors */
+  static kep_clock::time_point make_tp(const int y, const uint mon, const uint d, const int h = 0,
+                                       const int min = 0, const int s = 0, const int ms = 0, const int us = 0);
+
+  static kep_clock::time_point make_tp(double epoch_in, julian_type epoch_type);
+
+  // Conversions
+  static constexpr kep_clock::time_point tp_from_days(double days);
+
+  // Duration conversions
+  static constexpr double as_sec(kep_clock::duration &&d) {
+    return std::chrono::duration<double, std::chrono::seconds::period>(d)
+        .count();
+  }
+
+  // Printing
+  std::string as_utc_string() const;
+  static std::string as_utc_string(const kep_clock::time_point& tp);
+
+  /** operator overloads for sum and diff (epoch-days) and comparison
+   * operators
    * **/
-  epoch &operator+=(double rhs);
-  epoch &operator-=(double rhs);
+
+  kep3_DLL_PUBLIC friend std::ostream &operator<<(std::ostream &s,
+                                                  epoch const &epoch_in);
+
+  template <class Duration, class = enable_if_is_duration<Duration>>
+  epoch &operator+=(const Duration &duration) {
+    tp += chr::duration_cast<kep_clock::duration>(duration);
+    return *this;
+  }
+
+  template <class Duration, class = enable_if_is_duration<Duration>>
+  epoch &operator-=(const Duration &duration) {
+    tp -= chr::duration_cast<kep_clock::duration>(duration);
+    return *this;
+  }
+
   kep3_DLL_PUBLIC friend bool operator>(const epoch &c1, const epoch &c2);
   kep3_DLL_PUBLIC friend bool operator<(const epoch &c1, const epoch &c2);
   kep3_DLL_PUBLIC friend bool operator>=(const epoch &c1, const epoch &c2);
   kep3_DLL_PUBLIC friend bool operator<=(const epoch &c1, const epoch &c2);
   kep3_DLL_PUBLIC friend bool operator==(const epoch &c1, const epoch &c2);
   kep3_DLL_PUBLIC friend bool operator!=(const epoch &c1, const epoch &c2);
-  kep3_DLL_PUBLIC friend epoch operator+(epoch lhs, double rhs);
-  kep3_DLL_PUBLIC friend epoch operator-(epoch lhs, double rhs);
-  kep3_DLL_PUBLIC friend double operator-(const epoch &lhs, const epoch &rhs);
 
-private:
+  template <class Duration, class = enable_if_is_duration<Duration>>
+  epoch operator+(const Duration &duration) {
+    return epoch(tp + chr::duration_cast<kep_clock::duration>(duration));
+  }
+
+  template <class Duration, class = enable_if_is_duration<Duration>>
+  epoch operator-(const Duration &duration) {
+    return epoch(tp - chr::duration_cast<kep_clock::duration>(duration));
+  }
+
+  static constexpr auto days(const double value) {
+    return chr::duration_cast<kep_clock::duration>(
+        chr::duration<double, std::ratio<86400>>(value));
+  }
+
+  static constexpr auto sec(const double value) {
+    return chr::duration_cast<kep_clock::duration>(
+        chr::duration<double, std::ratio<1>>(value));
+  }
+
+  kep3_DLL_PUBLIC friend kep_clock::duration operator-(const epoch &lhs,
+                                                       const epoch &rhs);
+
+  private:
+
+  // Constructor for const time_point&)
+  explicit epoch(const kep_clock::time_point &time_point);
+
+  // Constructor for const time_point&&)
+  explicit epoch(kep_clock::time_point &&time_point);
+
   // Serialization code
   friend class boost::serialization::access;
-  template <class Archive> void serialize(Archive &ar, const unsigned int) {
-    ar &m_mjd2000;
+  template <class Archive> void serialize(Archive &ar, const uint) {
+    ar &boost::serialization::make_binary_object(&tp, sizeof(tp));
   }
   // Serialization code (END)
 
-  /// the modified julian date 2000 stored in a double
-  double m_mjd2000;
+  // Time point relative to 1 Jan 2000 (MJD2000)
+  kep_clock::time_point tp;
 };
 
-kep3_DLL_PUBLIC epoch epoch_from_string(const std::string &date);
-kep3_DLL_PUBLIC epoch epoch_from_iso_string(const std::string &date);
-
 kep3_DLL_PUBLIC std::ostream &operator<<(std::ostream &s,
                                          const epoch &epoch_in);
 
-
 } // end of namespace kep3
 
-template <> struct fmt::formatter<kep3::epoch> : ostream_formatter {};
+template <> struct fmt::formatter<kep3::epoch> : fmt::ostream_formatter {};
 
-#endif // kep3_EPOCH_H
+#endif // kep3_EPOCH_HPP

include/kep3/planets/keplerian.hpp

@@ -25,35 +25,38 @@ namespace kep3::udpla {
 class kep3_DLL_PUBLIC keplerian {
 
   kep3::epoch m_ref_epoch;
-  std::array<std::array<double, 3>, 2> m_pos_vel_0;
   std::string m_name;
   double m_mu_central_body;
   double m_mu_self;
   double m_radius;
   double m_safe_radius;
   double m_period;
   bool m_ellipse;
+  std::array<std::array<double, 3>, 2> m_pos_vel_0;
 
   friend class boost::serialization::access;
   template <typename Archive> void serialize(Archive &ar, unsigned) {
     ar &m_ref_epoch;
-    ar &m_pos_vel_0;
     ar &m_name;
     ar &m_mu_central_body;
     ar &m_mu_self;
     ar &m_radius;
     ar &m_safe_radius;
     ar &m_period;
     ar &m_ellipse;
+    ar &m_safe_radius;
+    ar &m_pos_vel_0;
   }
 
 public:
   // NOTE: added_param contains mu_self, radius and safe_radius
   explicit keplerian(const epoch &ref_epoch, const std::array<double, 6> &par,
                      double mu_central_body = 1., std::string name = "Unknown",
+
                      std::array<double, 3> added_params = {-1., -1., -1.}, kep3::elements_type el_t = kep3::elements_type::KEP_F);
   // Constructor from pos_vel
   explicit keplerian(const epoch &ref_epoch = kep3::epoch(0),
+
                      const std::array<std::array<double, 3>, 2> &pos_vel =
                          {{{1.0, 0.0, 0.0}, {0., 1.0, 0.0}}},
                      double mu_central_body = 1., std::string name = "Unknown",