diff --git a/.clang-format b/.clang-format index 3132f91d1..5ca562b50 100644 --- a/.clang-format +++ b/.clang-format @@ -56,6 +56,6 @@ SpacesInContainerLiterals: false SpacesInCStyleCastParentheses: false SpacesInParentheses: false SpacesInSquareBrackets: false -Standard: Cpp11 +Standard: Auto TabWidth: 4 UseTab: Never diff --git a/.clang-tidy b/.clang-tidy index 254b32a44..41de96a23 100644 --- a/.clang-tidy +++ b/.clang-tidy @@ -1,5 +1,5 @@ --- -Checks: 'clang-diagnostic-*,clang-analyzer-*,*,-hicpp-uppercase-literal-suffix,-readability-uppercase-literal-suffix,-modernize-use-trailing-return-type,-readability-named-parameter,-hicpp-named-parameter,-cppcoreguidelines-pro-bounds-array-to-pointer-decay,-hicpp-no-array-decay,-llvm-header-guard,-cppcoreguidelines-macro-usage,-google-runtime-references,-readability-isolate-declaration,-fuchsia-default-arguments-calls,-fuchsia-overloaded-operator,-fuchsia-default-arguments-declarations,-readability-else-after-return,-google-runtime-int,-hicpp-signed-bitwise,-cert-dcl21-cpp,-cppcoreguidelines-avoid-magic-numbers,-readability-magic-numbers,-cppcoreguidelines-pro-bounds-pointer-arithmetic,-cppcoreguidelines-pro-type-reinterpret-cast,-cppcoreguidelines-avoid-c-arrays,-hicpp-avoid-c-arrays,-modernize-avoid-c-arrays,-modernize-use-transparent-functors,-cert-dcl16-c,-cppcoreguidelines-pro-type-union-access,-bugprone-branch-clone,-fuchsia-statically-constructed-objects,-cppcoreguidelines-pro-bounds-constant-array-index,-readability-static-accessed-through-instance,-cppcoreguidelines-pro-type-vararg,-hicpp-vararg,-llvmlibc-restrict-system-libc-headers,-llvmlibc-callee-namespace,-llvmlibc-implementation-in-namespace,-llvm-else-after-return,-fuchsia-trailing-return,-readability-identifier-length,-altera-unroll-loops,-altera-id-dependent-backward-branch,-altera-struct-pack-align,-performance-no-int-to-ptr,-readability-function-cognitive-complexity' +Checks: 'clang-diagnostic-*,clang-analyzer-*,*,-hicpp-uppercase-literal-suffix,-readability-uppercase-literal-suffix,-modernize-use-trailing-return-type,-readability-named-parameter,-hicpp-named-parameter,-cppcoreguidelines-pro-bounds-array-to-pointer-decay,-hicpp-no-array-decay,-llvm-header-guard,-cppcoreguidelines-macro-usage,-google-runtime-references,-readability-isolate-declaration,-fuchsia-default-arguments-calls,-fuchsia-overloaded-operator,-fuchsia-default-arguments-declarations,-readability-else-after-return,-google-runtime-int,-hicpp-signed-bitwise,-cert-dcl21-cpp,-cppcoreguidelines-avoid-magic-numbers,-readability-magic-numbers,-cppcoreguidelines-pro-bounds-pointer-arithmetic,-cppcoreguidelines-pro-type-reinterpret-cast,-cppcoreguidelines-avoid-c-arrays,-hicpp-avoid-c-arrays,-modernize-avoid-c-arrays,-modernize-use-transparent-functors,-cert-dcl16-c,-cppcoreguidelines-pro-type-union-access,-bugprone-branch-clone,-fuchsia-statically-constructed-objects,-cppcoreguidelines-pro-bounds-constant-array-index,-readability-static-accessed-through-instance,-cppcoreguidelines-pro-type-vararg,-hicpp-vararg,-llvmlibc-restrict-system-libc-headers,-llvmlibc-callee-namespace,-llvmlibc-implementation-in-namespace,-llvm-else-after-return,-fuchsia-trailing-return,-readability-identifier-length,-altera-unroll-loops,-altera-id-dependent-backward-branch,-altera-struct-pack-align,-performance-no-int-to-ptr,-readability-function-cognitive-complexity,-misc-include-cleaner,-llvmlibc-inline-function-decl' WarningsAsErrors: '*' AnalyzeTemporaryDtors: false FormatStyle: none diff --git a/.github/workflows/gha_ci.yml b/.github/workflows/gha_ci.yml index a10a63f29..278c2c74e 100644 --- a/.github/workflows/gha_ci.yml +++ b/.github/workflows/gha_ci.yml @@ -114,6 +114,12 @@ jobs: fetch-depth: 2 - name: Build run: bash tools/gha_conda_coverage.sh + conda_llvm16_asan: + runs-on: ubuntu-latest + steps: + - uses: actions/checkout@v3 + - name: Build + run: bash tools/gha_llvm16_conda_asan.sh conda_llvm15_asan: runs-on: ubuntu-latest steps: diff --git a/CMakeLists.txt b/CMakeLists.txt index 1e5d8a6c3..d4e2d98ec 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -11,7 +11,7 @@ if(NOT CMAKE_BUILD_TYPE) FORCE) endif() -project(heyoka VERSION 2.0.0 LANGUAGES CXX C) +project(heyoka VERSION 2.0.1 LANGUAGES CXX C) list(APPEND CMAKE_MODULE_PATH "${CMAKE_CURRENT_SOURCE_DIR}/cmake" "${CMAKE_CURRENT_SOURCE_DIR}/cmake/yacma") diff --git a/doc/changelog.rst b/doc/changelog.rst index fc803b57e..e92bf3839 100644 --- a/doc/changelog.rst +++ b/doc/changelog.rst @@ -1,6 +1,16 @@ Changelog ========= +2.0.1 (unreleased) +------------------ + +Fix +~~~ + +- Fix orbital elements singularity when using the VSOP2013 + theory at low precision + (`#348 `__). + 2.0.0 (2023-09-22) ------------------ diff --git a/include/heyoka/expression.hpp b/include/heyoka/expression.hpp index 6c983fecb..2ec8059d4 100644 --- a/include/heyoka/expression.hpp +++ b/include/heyoka/expression.hpp @@ -20,11 +20,9 @@ #include #include #include -#include #include #include #include -#include #include #include #include @@ -79,7 +77,7 @@ class HEYOKA_DLL_PUBLIC expression template void serialize(Archive &ar, unsigned) { - ar &m_value; + ar & m_value; } public: @@ -138,9 +136,9 @@ HEYOKA_DLL_PUBLIC expression operator""_ldbl(unsigned long long); #if defined(HEYOKA_HAVE_REAL128) template -inline expression operator"" _f128() +inline expression operator""_f128() { - return expression{mppp::literals::operator"" _rq()}; + return expression{mppp::literals::operator""_rq < Chars... > ()}; } #endif diff --git a/src/llvm_state.cpp b/src/llvm_state.cpp index e7003123e..d84652112 100644 --- a/src/llvm_state.cpp +++ b/src/llvm_state.cpp @@ -364,7 +364,15 @@ struct llvm_state::jit { auto lljit = lljit_builder.create(); // LCOV_EXCL_START if (!lljit) { - throw std::invalid_argument("Error creating an LLJIT object"); + auto err = lljit.takeError(); + + std::string err_report; + llvm::raw_string_ostream ostr(err_report); + + ostr << err; + + throw std::invalid_argument( + fmt::format("Could not create an LLJIT object. The full error message is:\n{}", ostr.str())); } // LCOV_EXCL_STOP m_lljit = std::move(*lljit); @@ -579,6 +587,7 @@ auto llvm_state_bc_to_module(const std::string &module_name, const std::string & } // namespace detail +// NOLINTNEXTLINE(cppcoreguidelines-rvalue-reference-param-not-moved) llvm_state::llvm_state(std::tuple &&tup) : m_jitter(std::make_unique()), m_opt_level(std::get<1>(tup)), m_fast_math(std::get<2>(tup)), m_force_avx512(std::get<3>(tup)), m_slp_vectorize(std::get<4>(tup)), m_module_name(std::move(std::get<0>(tup))) @@ -1071,11 +1080,10 @@ void llvm_state::optimise() // NOTE: the reason for this inconsistency is that opt uses PB.parsePassPipeline() // (instead of PB.buildPerModuleDefaultPipeline()) to set up the optimisation // pipeline. Indeed, if we replace PB.buildPerModuleDefaultPipeline(ol) with - // PB.buildPerModuleDefaultPipeline(MPM, "default") (which - // corresponds to invoking "opt -passes='default'"), we do NOT need to set - // SLP vectorization on here to get the SLP vectorizer. Not sure if we should consider - // switching to this alternative way of setting up the optimisation pipeline - // in the future. + // PB.parsePassPipeline(MPM, "default") (which corresponds to invoking + // "opt -passes='default'"), we do NOT need to set SLP vectorization on + // here to get the SLP vectorizer. Not sure if we should consider switching to this + // alternative way of setting up the optimisation pipeline in the future. llvm::PipelineTuningOptions pto; pto.SLPVectorization = m_slp_vectorize; llvm::PassBuilder PB(m_jitter->m_tm.get(), pto); diff --git a/src/model/vsop2013.cpp b/src/model/vsop2013.cpp index 5ec31564a..96179289b 100644 --- a/src/model/vsop2013.cpp +++ b/src/model/vsop2013.cpp @@ -324,6 +324,12 @@ std::vector vsop2013_cartesian_impl(std::uint32_t pl_idx, expression qp = sqrt(qp_2); }); + // Check if qp is zero. This indicates a zero inclination orbit. + bool zero_inc = false; + if (const auto *num_ptr = std::get_if(&qp.value())) { + zero_inc = is_zero(*num_ptr); + } + // E, e, sqrt(1 - e**2), cos(i), sin(i), cos(Om), sin(Om), sin(E), cos(E) expression E, e, sqrt_1me2, ci, si, cOm, sOm, sin_E, cos_E; tbb::parallel_invoke( @@ -337,12 +343,26 @@ std::vector vsop2013_cartesian_impl(std::uint32_t pl_idx, expression ci = 1_dbl - 2_dbl * qp_2; si = sqrt(1_dbl - ci * ci); }, - [&]() { cOm = q / qp; }, [&]() { sOm = p / qp; }); + // NOTE: for zero inclination, set conventionally Om to zero. + [&]() { cOm = zero_inc ? 1_dbl : q / qp; }, [&]() { sOm = zero_inc ? 0_dbl : p / qp; }); + + // Check for zero eccentricity. + bool zero_ecc = false; + // NOTE: not sure if this is possible at all, but better safe than sorry. + // Let's exclude it from code coverage checks for the time being. + // LCOV_EXCL_START + if (const auto *num_ptr = std::get_if(&e.value())) { + zero_ecc = is_zero(*num_ptr); + } + // LCOV_EXCL_STOP // cos(om), sin(om), q1/a, q2/a. expression com, som, q1_a, q2_a; - tbb::parallel_invoke([&]() { com = (k * cOm + h * sOm) / e; }, [&]() { som = (h * cOm - k * sOm) / e; }, - [&]() { q1_a = cos_E - e; }, [&]() { q2_a = sqrt_1me2 * sin_E; }); + tbb::parallel_invoke( + // NOTE: for zero eccentricity, set conventionally om to zero. + [&]() { com = zero_ecc ? 1_dbl : (k * cOm + h * sOm) / e; }, + [&]() { som = zero_ecc ? 0_dbl : (h * cOm - k * sOm) / e; }, [&]() { q1_a = cos_E - e; }, + [&]() { q2_a = sqrt_1me2 * sin_E; }); // Prepare the entries of the rotation matrix, and a few auxiliary quantities. expression R00, R01, R10, R11, R20, R21, v_num, v_den; diff --git a/src/taylor_02.cpp b/src/taylor_02.cpp index e1f0e0764..61ae7c790 100644 --- a/src/taylor_02.cpp +++ b/src/taylor_02.cpp @@ -975,6 +975,7 @@ std::pair taylor_compute_jet_compact_mode( // // to pass the data necessary to the parallel workers. par_data_t = llvm::StructType::get(context, {builder.getInt32Ty(), ext_fp_ptr_t, ext_fp_ptr_t}); + // NOLINTNEXTLINE(cppcoreguidelines-owning-memory) gl_par_data = new llvm::GlobalVariable(md, par_data_t, false, llvm::GlobalVariable::InternalLinkage, llvm::ConstantAggregateZero::get(par_data_t)); diff --git a/test/opt_checks.cpp b/test/opt_checks.cpp index cd09b55e3..d248bbcf7 100644 --- a/test/opt_checks.cpp +++ b/test/opt_checks.cpp @@ -63,9 +63,10 @@ TEST_CASE("pow vect") REQUIRE(!boost::algorithm::contains(md_ir, "llvm.pow")); } -#if LLVM_VERSION_MAJOR >= 16 + // NOTE: no idea why, but in LLVM 17 it seems like there's no autovec + // happening at all in this specific case. This needs to be investigated. +#if LLVM_VERSION_MAJOR == 16 - // NOTE: LLVM16 is currently the version tested in the CI on arm64. if (tf.aarch64) { REQUIRE(!boost::algorithm::contains(md_ir, "llvm.pow")); } diff --git a/test/vsop2013.cpp b/test/vsop2013.cpp index 7d9fdae9e..6ce2d1721 100644 --- a/test/vsop2013.cpp +++ b/test/vsop2013.cpp @@ -14,6 +14,7 @@ #include #include +#include #include #include @@ -780,3 +781,26 @@ TEST_CASE("vsop2013 mus") REQUIRE(mus[3] == 8.9970116036316091182e-10); REQUIRE(mus[9] == 2.1886997654259696800e-12); } + +// Bug: low-precision solutions may have zero inclination, +// which leads to singularities when converting to Cartesian. +TEST_CASE("vsop2013 low prec zero inc") +{ + auto ex = vsop2013_cartesian(1, kw::time = "tm"_var, kw::thresh = 1e-1)[0]; + + llvm_state s; + + add_cfunc(s, "f", {ex}); + + s.compile(); + + auto *cf_ptr + = reinterpret_cast(s.jit_lookup("f")); + + double out = 0, in = 0; + + cf_ptr(&out, &in, nullptr, nullptr); + + REQUIRE(!std::isnan(out)); + REQUIRE(out != 0); +} diff --git a/tools/circleci_ubuntu_arm64.sh b/tools/circleci_ubuntu_arm64.sh index 9bffa9b96..390b19719 100644 --- a/tools/circleci_ubuntu_arm64.sh +++ b/tools/circleci_ubuntu_arm64.sh @@ -14,7 +14,7 @@ wget https://github.com/conda-forge/miniforge/releases/latest/download/Mambaforg export deps_dir=$HOME/local export PATH="$HOME/miniconda/bin:$PATH" bash miniconda.sh -b -p $HOME/miniconda -mamba create -y -q -p $deps_dir cxx-compiler c-compiler cmake llvmdev tbb-devel tbb boost-cpp sleef xtensor xtensor-blas blas blas-devel fmt spdlog 'mppp>=0.27' +mamba create -y -p $deps_dir cxx-compiler c-compiler cmake llvmdev tbb-devel tbb boost-cpp sleef xtensor xtensor-blas blas blas-devel fmt spdlog 'mppp>=0.27' source activate $deps_dir # Create the build dir and cd into it. diff --git a/tools/gha_conda_release_static.sh b/tools/gha_conda_release_static.sh index 1e8f79ab0..60d864d57 100644 --- a/tools/gha_conda_release_static.sh +++ b/tools/gha_conda_release_static.sh @@ -16,7 +16,7 @@ export PATH="$HOME/miniconda/bin:$PATH" bash miniconda.sh -b -p $HOME/miniconda conda config --add channels conda-forge conda config --set channel_priority strict -conda create -y -q -p $deps_dir c-compiler cxx-compiler cmake llvmdev tbb-devel tbb boost-cpp 'mppp>=0.27' sleef xtensor xtensor-blas blas blas-devel fmt spdlog zlib +conda create -y -q -p $deps_dir c-compiler cxx-compiler cmake 'llvmdev=16.*' tbb-devel tbb boost-cpp 'mppp>=0.27' sleef xtensor xtensor-blas blas blas-devel fmt spdlog zlib source activate $deps_dir # Create the build dir and cd into it. diff --git a/tools/gha_llvm16_conda_asan.sh b/tools/gha_llvm16_conda_asan.sh new file mode 100755 index 000000000..24c930669 --- /dev/null +++ b/tools/gha_llvm16_conda_asan.sh @@ -0,0 +1,30 @@ +#!/usr/bin/env bash + +# Echo each command +set -x + +# Exit on error. +set -e + +# Core deps. +sudo apt-get install wget + +# Install conda+deps. +wget https://github.com/conda-forge/miniforge/releases/latest/download/Mambaforge-Linux-x86_64.sh -O mambaforge.sh +export deps_dir=$HOME/local +export PATH="$HOME/mambaforge/bin:$PATH" +bash mambaforge.sh -b -p $HOME/mambaforge +mamba create -y -q -p $deps_dir c-compiler cxx-compiler cmake 'llvmdev=16.*' tbb-devel tbb boost-cpp 'mppp>=0.27' sleef xtensor xtensor-blas blas blas-devel 'fmt=9.*' spdlog +source activate $deps_dir + +# Create the build dir and cd into it. +mkdir build +cd build + +# GCC build. +cmake ../ -DCMAKE_PREFIX_PATH=$deps_dir -DCMAKE_BUILD_TYPE=Debug -DHEYOKA_BUILD_TESTS=yes -DHEYOKA_BUILD_TUTORIALS=ON -DHEYOKA_WITH_MPPP=yes -DHEYOKA_WITH_SLEEF=yes -DCMAKE_CXX_FLAGS="-fsanitize=address" -DBoost_NO_BOOST_CMAKE=ON +make -j2 VERBOSE=1 +ctest -V -j2 -E vsop2013 + +set +e +set +x