Skip to content

Commit

Permalink
WIP on testing.
Browse files Browse the repository at this point in the history
  • Loading branch information
bluescarni committed Jul 4, 2024
1 parent 2dd7363 commit 667b5d9
Show file tree
Hide file tree
Showing 2 changed files with 89 additions and 2 deletions.
4 changes: 4 additions & 0 deletions include/heyoka/model/sgp4.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -69,6 +69,10 @@ class HEYOKA_DLL_PUBLIC_INLINE_CLASS sgp4_propagator
template <typename Input, typename... KwArgs>
static auto parse_ctor_args(const Input &in, const KwArgs &...kw_args)
{
if (in.data_handle() == nullptr) [[unlikely]] {
throw std::invalid_argument("Cannot initialise an sgp4_propagator with a null list of satellites");
}

if (in.extent(1) == 0u) [[unlikely]] {
throw std::invalid_argument("Cannot initialise an sgp4_propagator with an empty list of satellites");
}
Expand Down
87 changes: 85 additions & 2 deletions test/model_sgp4.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -7,6 +7,8 @@
// with this file, You can obtain one at http://mozilla.org/MPL/2.0/.

#include <cstddef>
#include <ranges>
#include <stdexcept>
#include <vector>

#include <boost/math/constants/constants.hpp>
Expand Down Expand Up @@ -150,7 +152,7 @@ TEST_CASE("propagator")
const prop_t::in_1d<double> tm_in{tm.data(), 2};

for (auto cm : {false, true}) {
prop_t prop{md_input_t{ins.data(), 9, 2}, kw::compact_mode = cm};
prop_t prop{md_input_t{ins.data(), 2}, kw::compact_mode = cm};

std::vector<double> outs(12u);
prop_t::out_2d out{outs.data(), 6, 2};
Expand Down Expand Up @@ -222,7 +224,7 @@ TEST_CASE("propagator batch")
const prop_t::in_2d<double> tm_in{tm.data(), 2, 2};

for (auto cm : {false, true}) {
prop_t prop{md_input_t{ins.data(), 9, 2}, kw::compact_mode = cm};
prop_t prop{md_input_t{ins.data(), 2}, kw::compact_mode = cm};

std::vector<double> outs(24u);
prop_t::out_3d out{outs.data(), 2, 6, 2};
Expand Down Expand Up @@ -288,3 +290,84 @@ TEST_CASE("propagator batch")
REQUIRE(out(1, 5, 1) == approximately(1.785277174529374, 10000.));
}
}

TEST_CASE("error handling")
{
detail::edb_disabler ed;

using Catch::Matchers::Message;
using md_input_t = mdspan<const double, extents<std::size_t, 9, std::dynamic_extent>>;
using prop_t = model::sgp4_propagator<double>;

// Propagator with null list or zero satellites.
REQUIRE_THROWS_MATCHES((prop_t{md_input_t{nullptr, 0}}), std::invalid_argument,
Message("Cannot initialise an sgp4_propagator with a null list of satellites"));

std::vector<double> input(9u);

REQUIRE_THROWS_MATCHES((prop_t{md_input_t{input.data(), 0}}), std::invalid_argument,
Message("Cannot initialise an sgp4_propagator with an empty list of satellites"));
}

TEST_CASE("derivatives")
{
detail::edb_disabler ed;

using md_input_t = mdspan<const double, extents<std::size_t, 9, std::dynamic_extent>>;
using prop_t = model::sgp4_propagator<double>;

// First compute the order-2 derivatives of the whole model.
const auto sgp4_func = model::sgp4();
const auto inputs = make_vars("n0", "e0", "i0", "node0", "omega0", "m0", "bstar", "tsince");
const auto dt = diff_tensors(sgp4_func, std::vector(inputs.begin(), inputs.begin() + 6), kw::diff_order = 2);

// Make a compiled function with the derivatives.
auto diff_cf = cfunc<double>(dt | std::views::transform([](const auto &p) { return p.second; }), inputs);

// Create a propagator with derivatives.
const std::vector<double> ins = {revday2radmin(13.75091047972192),
revday2radmin(15.50103472202482),
0.0024963,
0.0007417,
deg2rad(90.2039),
deg2rad(51.6439),
deg2rad(55.5633),
deg2rad(211.2001),
deg2rad(320.5956),
deg2rad(17.6667),
deg2rad(91.4738),
deg2rad(85.6398),
0.75863e-3,
.38792e-4,
2460486.5,
2458826.5,
0.6478633000000116,
0.6933954099999937};

const auto tm = std::array{1440., 0.};
const prop_t::in_1d<double> tm_in{tm.data(), 2};

prop_t prop{md_input_t{ins.data(), 2}, kw::diff_order = 2};

// Prepare the input buffer for the cfunc.
std::vector<double> cf_in(ins.begin(), ins.begin() + 14);
cf_in.insert(cf_in.end(), tm.begin(), tm.end());
cfunc<double>::in_2d cf_in_span(cf_in.data(), 8, 2);

// Prepare the output buffers.
std::vector<double> cf_out(dt.size() * 2u), prop_out(cf_out);
cfunc<double>::out_2d cf_out_span(cf_out.data(), dt.size(), 2);
cfunc<double>::out_2d prop_out_span(prop_out.data(), dt.size(), 2);

// Evaluate the cfunc.
diff_cf(cf_out_span, cf_in_span);

// Evaluate the propagation.
prop(prop_out_span, tm_in);

for (std::size_t i = 0; i < prop_out_span.extent(0); ++i) {
for (std::size_t j = 0; j < prop_out_span.extent(1); ++j) {
REQUIRE(prop_out_span(i, j) == approximately(cf_out_span(i, j)));
}
}
}

0 comments on commit 667b5d9

Please sign in to comment.