Skip to content

Commit

Permalink
Merge pull request #1 from esa/planet_class
Browse files Browse the repository at this point in the history
Planet class
  • Loading branch information
darioizzo authored Aug 25, 2023
2 parents be44ba1 + cbe7690 commit d022e0e
Show file tree
Hide file tree
Showing 27 changed files with 2,320 additions and 1,158 deletions.
394 changes: 394 additions & 0 deletions .clang-tidy

Large diffs are not rendered by default.

2 changes: 2 additions & 0 deletions .clangd
Original file line number Diff line number Diff line change
@@ -0,0 +1,2 @@
CompileFlags:
Add: '-std=c++17'
11 changes: 3 additions & 8 deletions .gitignore
Original file line number Diff line number Diff line change
@@ -1,16 +1,11 @@
# Build dirs.
build*

# Docs.
doc/_build
doc/examples/out
doc/generated

# Notebooks bits.
.ipynb_checkpoints

# clangd bits.
.cache

# vscode bits.
.vscode

# config files
include/kep3/config.hpp
7 changes: 6 additions & 1 deletion CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -130,6 +130,8 @@ endif()
# List of source files.
set(kep3_SRC_FILES
"${CMAKE_CURRENT_SOURCE_DIR}/src/epoch.cpp"
"${CMAKE_CURRENT_SOURCE_DIR}/src/planet.cpp"
"${CMAKE_CURRENT_SOURCE_DIR}/src/detail/type_name.cpp"
)

# Setup of the kep3 shared library.
Expand Down Expand Up @@ -176,9 +178,12 @@ endif()
find_package(spdlog CONFIG REQUIRED)
target_link_libraries(kep3 PRIVATE spdlog::spdlog)

# Configure config.hpp.
configure_file("${CMAKE_CURRENT_SOURCE_DIR}/config.hpp.in" "${CMAKE_CURRENT_SOURCE_DIR}/include/kep3/config.hpp" @ONLY)

# Installation of the header files.
install(DIRECTORY "${CMAKE_CURRENT_SOURCE_DIR}/include/kep3" DESTINATION include)
#install(FILES "${CMAKE_CURRENT_BINARY_DIR}/include/kep3/config.hpp" DESTINATION include/kep3)
install(FILES "${CMAKE_CURRENT_SOURCE_DIR}/include/kep3/config.hpp" DESTINATION include/kep3)

# Installation of the library.
install(TARGETS kep3
Expand Down
373 changes: 373 additions & 0 deletions COPYING

Large diffs are not rendered by default.

674 changes: 0 additions & 674 deletions LICENSE

This file was deleted.

151 changes: 75 additions & 76 deletions benchmark/convert_anomalies_benchmark.cpp
Original file line number Diff line number Diff line change
@@ -1,6 +1,7 @@
// Copyright 2020, 2021, 2022 Francesco Biscani ([email protected]), Dario Izzo ([email protected])
// Copyright 2023, 2024 Dario Izzo ([email protected]), Francesco Biscani
// ([email protected])
//
// This file is part of the dsyre library.
// This file is part of the kep3 library.
//
// This Source Code Form is subject to the terms of the Mozilla
// Public License v. 2.0. If a copy of the MPL was not distributed
Expand All @@ -12,91 +13,89 @@

#include <fmt/core.h>

#include<kep3/core_astro/convert_anomalies.hpp>

#include <kep3/core_astro/convert_anomalies.hpp>

using namespace kep3;
using namespace std::chrono;

// In this benchmark we test the speed and accuracy of the Kepler's equation solvers
// In this benchmark we test the speed and accuracy of the Kepler's equation
// solvers

void perform_test_speed(double min_ecc, double max_ecc, unsigned N) {
//
// Engines
//
std::mt19937 rng_engine(122012203u);
//
// Distribtuions
//
std::uniform_real_distribution<double> ecc_d(min_ecc, max_ecc);
std::uniform_real_distribution<double> M_d(-100, 100.);

// We generate the random dataset
std::vector<double> eccenricities(N);
std::vector<double> mean_anomalies(N);

void perform_test_speed(double min_ecc, double max_ecc, unsigned N)
{
//
// Engines
//
std::mt19937 rng_engine(122012203u);
//
// Distribtuions
//
std::uniform_real_distribution<double> ecc_d(min_ecc, max_ecc);
std::uniform_real_distribution<double> M_d(-100, 100.);

// We generate the random dataset
std::vector<double> eccenricities(N);
std::vector<double> mean_anomalies(N);
for (auto i = 0u; i < N; ++i) {
mean_anomalies[i] = M_d(rng_engine);
eccenricities[i] = ecc_d(rng_engine);
}

for (auto i = 0u; i < N; ++i) {
mean_anomalies[i] = M_d(rng_engine);
eccenricities[i] = ecc_d(rng_engine);
}
// We log progress
fmt::print("{:.2f} min_ecc, {:.2f} max_ecc, on {} data points: ", min_ecc,
max_ecc, N);

// We log progress
fmt::print("{:.2f} min_ecc, {:.2f} max_ecc, on {} data points: ", min_ecc, max_ecc, N);

auto start = high_resolution_clock::now();
for (auto i = 0u; i < N; ++i) {
m2e(mean_anomalies[i], eccenricities[i]);
}
auto stop = high_resolution_clock::now();
auto duration = duration_cast<microseconds>(stop - start);
fmt::print("{:.3f}s\n", duration.count() / 1e6);
auto start = high_resolution_clock::now();
for (auto i = 0u; i < N; ++i) {
m2e(mean_anomalies[i], eccenricities[i]);
}
auto stop = high_resolution_clock::now();
auto duration = duration_cast<microseconds>(stop - start);
fmt::print("{:.3f}s\n", duration.count() / 1e6);
}

void perform_test_accuracy(double min_ecc, double max_ecc, unsigned N)
{
//
// Engines
//
std::mt19937 rng_engine(122012203u);
//
// Distribtuions
//
std::uniform_real_distribution<double> ecc_d(min_ecc, max_ecc);
std::uniform_real_distribution<double> M_d(-100, 100.);

// We generate the random dataset
std::vector<double> eccenricities(N);
std::vector<double> mean_anomalies(N);
void perform_test_accuracy(double min_ecc, double max_ecc, unsigned N) {
//
// Engines
//
std::mt19937 rng_engine(122012203u);
//
// Distribtuions
//
std::uniform_real_distribution<double> ecc_d(min_ecc, max_ecc);
std::uniform_real_distribution<double> M_d(-100, 100.);

for (auto i = 0u; i < N; ++i) {
mean_anomalies[i] = M_d(rng_engine);
eccenricities[i] = ecc_d(rng_engine);
}
// We generate the random dataset
std::vector<double> eccenricities(N);
std::vector<double> mean_anomalies(N);

// We log progress
fmt::print("{:.2f} min_ecc, {:.2f} max_ecc, on {} data points: ", min_ecc, max_ecc, N);
std::vector<double> err(N);
for (auto i = 0u; i < N; ++i) {
auto res = e2m(m2e(mean_anomalies[i], eccenricities[i]), eccenricities[i]);
err[i] = std::abs(res - mean_anomalies[i]);
}
auto max_it = max_element(std::begin(err), std::end(err));
auto min_it = min_element(std::begin(err), std::end(err));
auto avg = std::accumulate(err.begin(), err.end(), 0.0) / err.size();
fmt::print("{:.3e} avg, {:.3e} min, {:.3e} max\n", avg, *min_it, *max_it );
}
for (auto i = 0u; i < N; ++i) {
mean_anomalies[i] = M_d(rng_engine);
eccenricities[i] = ecc_d(rng_engine);
}

// We log progress
fmt::print("{:.2f} min_ecc, {:.2f} max_ecc, on {} data points: ", min_ecc,
max_ecc, N);
std::vector<double> err(N);
for (auto i = 0u; i < N; ++i) {
auto res = e2m(m2e(mean_anomalies[i], eccenricities[i]), eccenricities[i]);
err[i] = std::abs(res - mean_anomalies[i]);
}
auto max_it = max_element(std::begin(err), std::end(err));
auto min_it = min_element(std::begin(err), std::end(err));
auto avg = std::accumulate(err.begin(), err.end(), 0.0) / err.size();
fmt::print("{:.3e} avg, {:.3e} min, {:.3e} max\n", avg, *min_it, *max_it);
}

int main()
{
unsigned seed = 7898935u;
fmt::print("\nComputes speed different eccentricity ranges:\n");
perform_test_speed(0, 0.5, 1000000);
perform_test_speed(0.5, 0.9, 1000000);
perform_test_speed(0.9, 0.99, 1000000);
fmt::print("\nComputes error at different eccentricity ranges:\n");
perform_test_accuracy(0, 0.5, 100000);
perform_test_accuracy(0.5, 0.9, 100000);
perform_test_accuracy(0.9, 0.99, 100000);
int main() {
unsigned seed = 7898935u;
fmt::print("\nComputes speed different eccentricity ranges:\n");
perform_test_speed(0, 0.5, 1000000);
perform_test_speed(0.5, 0.9, 1000000);
perform_test_speed(0.9, 0.99, 1000000);
fmt::print("\nComputes error at different eccentricity ranges:\n");
perform_test_accuracy(0, 0.5, 100000);
perform_test_accuracy(0.5, 0.9, 100000);
perform_test_accuracy(0.9, 0.99, 100000);
}
2 changes: 1 addition & 1 deletion clang-format
Original file line number Diff line number Diff line change
Expand Up @@ -9,7 +9,7 @@ AllowAllParametersOfDeclarationOnNextLine: true
AllowShortBlocksOnASingleLine: false
AllowShortCaseLabelsOnASingleLine: false
AllowShortFunctionsOnASingleLine: Empty
AllowShortIfStatementsOnASingleLine: true
AllowShortIfStatementsOnASingleLine: false
AllowShortLoopsOnASingleLine: false
AlwaysBreakAfterDefinitionReturnType: None
AlwaysBreakBeforeMultilineStrings: false
Expand Down
34 changes: 34 additions & 0 deletions config.hpp.in
Original file line number Diff line number Diff line change
@@ -0,0 +1,34 @@
// Copyright 2023, 2024 Dario Izzo ([email protected]), Francesco Biscani
// ([email protected])
//
// This file is part of the kep3 library.
//
// This Source Code Form is subject to the terms of the Mozilla
// 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_CONFIG_HPP
#define kep3_CONFIG_HPP

// NOTE: include this so that we can
// detect _LIBCPP_VERSION below.
#include <ciso646>

// Start of defines instantiated by CMake.
#define kep3_VERSION "@kep3_VERSION@"
#define kep3_VERSION_MAJOR @kep3_VERSION_MAJOR@
#define kep3_VERSION_MINOR @kep3_VERSION_MINOR@
#define kep3_VERSION_PATCH @kep3_VERSION_PATCH@

// End of defines instantiated by CMake.

#if defined(__clang__) && defined(_LIBCPP_VERSION)

// When using clang + libc++, prefer the name-based
// extract() implementation for UDx classes. See
// the explanation in typeid_name_extract.hpp.
#define kep3_PREFER_TYPEID_NAME_EXTRACT

#endif

#endif
41 changes: 14 additions & 27 deletions include/kep3/core_astro/convert_anomalies.hpp
Original file line number Diff line number Diff line change
@@ -1,27 +1,11 @@
/*****************************************************************************
* Copyright (C) 2023 The pykep development team, *
* Advanced Concepts Team (ACT), European Space Agency (ESA) *
* *
* https://gitter.im/esa/pykep *
* https://github.com/esa/pykep *
* *
* [email protected] *
* *
* This program is free software; you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation; either version 2 of the License, or *
* (at your option) any later version. *
* *
* This program is distributed in the hope that it will be useful, *
* but WITHOUT ANY WARRANTY; without even the implied warranty of *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
* GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with this program; if not, write to the *
* Free Software Foundation, Inc., *
* 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. *
*****************************************************************************/
// Copyright 2023, 2024 Dario Izzo ([email protected]), Francesco Biscani
// ([email protected])
//
// This file is part of the kep3 library.
//
// This Source Code Form is subject to the terms of the Mozilla
// 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_ANOMALIES_H
#define kep3_CONVERT_ANOMALIES_H
Expand All @@ -35,17 +19,20 @@

namespace kep3 {

// mean to eccentric (only ellipses) e<1. Preserves the sign and integer number of revolutions.
// mean to eccentric (only ellipses) e<1. Preserves the sign and integer number
// of revolutions.
inline double m2e(double M, double ecc) {
// We use as intial guess the Mean Anomaly
// (tests indicated that any higher order expansion does not really improve)
double IG = M;
const int digits = std::numeric_limits<double>::digits;
double sol = boost::math::tools::halley_iterate(
[M, ecc](double E) {
return std::make_tuple(kepE(E, M, ecc), d_kepE(E, ecc), dd_kepE(E, ecc));
return std::make_tuple(kepE(E, M, ecc), d_kepE(E, ecc),
dd_kepE(E, ecc));
},
IG, IG - boost::math::constants::pi<double>(), IG + boost::math::constants::pi<double>(), digits);
IG, IG - boost::math::constants::pi<double>(),
IG + boost::math::constants::pi<double>(), digits);
return sol;
}
// eccentric to mean (only ellipses) e<1
Expand Down
67 changes: 16 additions & 51 deletions include/kep3/core_astro/convert_julian_dates.hpp
Original file line number Diff line number Diff line change
@@ -1,57 +1,22 @@
/*****************************************************************************
* Copyright (C) 2023 The pykep development team, *
* Advanced Concepts Team (ACT), European Space Agency (ESA) *
* *
* https://gitter.im/esa/pykep *
* https://github.com/esa/pykep *
* *
* [email protected] *
* *
* This program is free software; you can redistribute it and/or modify *
* it under the terms of the GNU General Public License as published by *
* the Free Software Foundation; either version 2 of the License, or *
* (at your option) any later version. *
* *
* This program is distributed in the hope that it will be useful, *
* but WITHOUT ANY WARRANTY; without even the implied warranty of *
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
* GNU General Public License for more details. *
* *
* You should have received a copy of the GNU General Public License *
* along with this program; if not, write to the *
* Free Software Foundation, Inc., *
* 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. *
*****************************************************************************/
// Copyright 2023, 2024 Dario Izzo ([email protected]), Francesco Biscani
// ([email protected])
//
// This file is part of the kep3 library.
//
// This Source Code Form is subject to the terms of the Mozilla
// 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

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
Loading

0 comments on commit d022e0e

Please sign in to comment.