From 6af39d8726aaa5baf3b1f18e9762713cd076b598 Mon Sep 17 00:00:00 2001 From: Florian Weik Date: Fri, 18 Oct 2019 16:36:52 +0200 Subject: [PATCH 1/7] utils: Implemented tuple protocol for Array --- src/utils/include/utils/Array.hpp | 15 +++++++++++++++ src/utils/include/utils/get.hpp | 9 ++++++++- src/utils/tests/Array_test.cpp | 10 ++++++++++ 3 files changed, 33 insertions(+), 1 deletion(-) diff --git a/src/utils/include/utils/Array.hpp b/src/utils/include/utils/Array.hpp index 94fa5a62854..9074735aae2 100644 --- a/src/utils/include/utils/Array.hpp +++ b/src/utils/include/utils/Array.hpp @@ -20,6 +20,7 @@ #define UTILS_ARRAY_HPP #include "device_qualifier.hpp" +#include "get.hpp" #include @@ -158,5 +159,19 @@ template struct Array { ar &m_storage; } }; + +template +struct tuple_element> { + using type = T; +}; + +template +struct tuple_size> : std::integral_constant {}; + +template +auto get(Array const &a) -> std::enable_if_t<(I < N), const T &> { + return a[I]; +} + } // namespace Utils #endif diff --git a/src/utils/include/utils/get.hpp b/src/utils/include/utils/get.hpp index 683ea4edcbc..b99d8f54667 100644 --- a/src/utils/include/utils/get.hpp +++ b/src/utils/include/utils/get.hpp @@ -20,10 +20,17 @@ #define UTILS_GET_HPP namespace Utils { - template auto get(const T &v) { return std::get(v); } +template struct tuple_size : std::tuple_size {}; + +template +struct tuple_element : std::tuple_element {}; + +template +using tuple_element_t = typename tuple_element::type; + } // namespace Utils #endif diff --git a/src/utils/tests/Array_test.cpp b/src/utils/tests/Array_test.cpp index cab5a7fb299..cb4a855e0c5 100644 --- a/src/utils/tests/Array_test.cpp +++ b/src/utils/tests/Array_test.cpp @@ -103,3 +103,13 @@ BOOST_AUTO_TEST_CASE(zero_size) { Array const a{}; BOOST_CHECK(a.empty()); } + +BOOST_AUTO_TEST_CASE(tuple_protocol) { + using A = Utils::Array; + + static_assert(std::is_same, int>::value, ""); + static_assert(std::is_same, int>::value, ""); + static_assert(A{}.size() == Utils::tuple_size::value, ""); + + BOOST_CHECK_EQUAL(Utils::get<1>(A{1, 2, 3, 4}), 2); +} \ No newline at end of file From 2a01f9bec76161d71ac8b9f6afe5fdd99afa4dcd Mon Sep 17 00:00:00 2001 From: Florian Weik Date: Fri, 18 Oct 2019 16:48:39 +0200 Subject: [PATCH 2/7] utils: Utils::mask + test --- src/utils/include/utils/get.hpp | 1 - src/utils/include/utils/mask.hpp | 37 +++++++++++++++++++++++++ src/utils/tests/CMakeLists.txt | 1 + src/utils/tests/mask_test.cpp | 47 ++++++++++++++++++++++++++++++++ 4 files changed, 85 insertions(+), 1 deletion(-) create mode 100644 src/utils/include/utils/mask.hpp create mode 100644 src/utils/tests/mask_test.cpp diff --git a/src/utils/include/utils/get.hpp b/src/utils/include/utils/get.hpp index b99d8f54667..48a99fa5331 100644 --- a/src/utils/include/utils/get.hpp +++ b/src/utils/include/utils/get.hpp @@ -31,6 +31,5 @@ struct tuple_element : std::tuple_element {}; template using tuple_element_t = typename tuple_element::type; - } // namespace Utils #endif diff --git a/src/utils/include/utils/mask.hpp b/src/utils/include/utils/mask.hpp new file mode 100644 index 00000000000..0e3da15eebc --- /dev/null +++ b/src/utils/include/utils/mask.hpp @@ -0,0 +1,37 @@ +#ifndef ESPRESSO_MASK_HPP +#define ESPRESSO_MASK_HPP + +#include +#include + +#include +#include + +namespace Utils { +namespace detail { +template +auto mask_impl(Integral mask, T t, std::index_sequence) { + return T{((mask & (1u << I)) ? get(t) : tuple_element_t{})...}; +} +} // namespace detail + +/** + * @brief Set elements of a tuple-like to zero by a bit mask. + * + * @tparam T implements the tuple interface(get, tuple_size, ...) + * @tparam Integral An unsigned integral type + * @param mask bit mask, if the i-th bit is set, the i-th element + * in @param t is set to zero. + * @param t + * @return t partially zeroed out according to mask + */ +template +std::enable_if_t::value && + (8 * sizeof(Integral) >= tuple_size::value), + T> +mask(Integral mask, T t) { + return detail::mask_impl(mask, t, std::make_index_sequence::value>{}); +} +} // namespace Utils + +#endif // ESPRESSO_MASK_HPP diff --git a/src/utils/tests/CMakeLists.txt b/src/utils/tests/CMakeLists.txt index a197f71c19d..cf397593f2e 100644 --- a/src/utils/tests/CMakeLists.txt +++ b/src/utils/tests/CMakeLists.txt @@ -37,6 +37,7 @@ unit_test(NAME sampling_test SRC sampling_test.cpp DEPENDS utils) unit_test(NAME coordinate_transformation_test SRC coordinate_transformation.cpp DEPENDS utils) unit_test(NAME rotation_matrix_test SRC rotation_matrix_test.cpp DEPENDS utils) unit_test(NAME quaternion_test SRC quaternion_test.cpp DEPENDS utils) +unit_test(NAME mask_test SRC mask_test.cpp DEPENDS utils) unit_test(NAME gather_buffer_test SRC gather_buffer_test.cpp DEPENDS utils Boost::mpi MPI::MPI_CXX NUM_PROC 4) unit_test(NAME scatter_buffer_test SRC scatter_buffer_test.cpp DEPENDS utils Boost::mpi MPI::MPI_CXX NUM_PROC 4) diff --git a/src/utils/tests/mask_test.cpp b/src/utils/tests/mask_test.cpp new file mode 100644 index 00000000000..f01d80cd65b --- /dev/null +++ b/src/utils/tests/mask_test.cpp @@ -0,0 +1,47 @@ +/* + * Copyright (C) 2019 The ESPResSo project + * + * This file is part of ESPResSo. + * + * ESPResSo 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 3 of the License, or + * (at your option) any later version. + * + * ESPResSo 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, see . + */ + +#define BOOST_TEST_MODULE mask test +#define BOOST_TEST_DYN_LINK +#include + +#include + +#include +#include +#include + +BOOST_AUTO_TEST_CASE(mask_) { + using Utils::get; + using Utils::tuple_element_t; + + const uint8_t mask = 1u | 4u; + + auto const a = std::make_tuple(std::string("a"), 3, 4.5); + using input_type = decltype(a); + + auto const result = Utils::mask(mask, a); + using result_type = decltype(result); + + static_assert(std::is_same::value, ""); + + BOOST_CHECK_EQUAL(get<0>(result), get<0>(a)); + BOOST_CHECK_EQUAL(get<1>(result), 0); + BOOST_CHECK_EQUAL(get<2>(result), get<2>(a)); +} From 03dccc3a01ecf9c050a386b51e61b11aec449247 Mon Sep 17 00:00:00 2001 From: Florian Weik Date: Fri, 18 Oct 2019 18:56:28 +0200 Subject: [PATCH 3/7] utils: tuple protocol for vector --- src/core/rotation.cpp | 10 +++------- src/utils/include/utils/Vector.hpp | 12 ++++++++++++ 2 files changed, 15 insertions(+), 7 deletions(-) diff --git a/src/core/rotation.cpp b/src/core/rotation.cpp index abe1bb4ffa8..3a8c4841c17 100644 --- a/src/core/rotation.cpp +++ b/src/core/rotation.cpp @@ -51,6 +51,7 @@ #include #include +#include /** Calculate the derivatives of the quaternion and angular acceleration * for a given particle @@ -135,19 +136,14 @@ void define_Qdd(Particle const &p, double Qd[4], double Qdd[4], double S[3], void propagate_omega_quat_particle(Particle &p) { // If rotation for the particle is disabled entirely, return early. - if (!p.p.rotation) + if (p.p.rotation == ROTATION_FIXED) return; Utils::Vector4d Qd{}, Qdd{}; Utils::Vector3d S{}, Wd{}; // Clear rotational velocity for blocked rotation axes. - if (!(p.p.rotation & ROTATION_X)) - p.m.omega[0] = 0; - if (!(p.p.rotation & ROTATION_Y)) - p.m.omega[1] = 0; - if (!(p.p.rotation & ROTATION_Z)) - p.m.omega[2] = 0; + p.m.omega = Utils::mask(p.p.rotation, p.m.omega); define_Qdd(p, Qd.data(), Qdd.data(), S.data(), Wd.data()); diff --git a/src/utils/include/utils/Vector.hpp b/src/utils/include/utils/Vector.hpp index 2fc2571276f..dd23081c683 100644 --- a/src/utils/include/utils/Vector.hpp +++ b/src/utils/include/utils/Vector.hpp @@ -368,5 +368,17 @@ template struct decay_to_scalar> { template struct decay_to_scalar> { using type = T; }; +template +struct tuple_element> { + using type = T; +}; + +template +struct tuple_size> : std::integral_constant {}; + +template +auto get(Vector const &a) -> std::enable_if_t<(I < N), const T &> { + return a[I]; +} } // namespace Utils #endif From b44f3a2a2283fc0cf95e98f829a538d751a1bff9 Mon Sep 17 00:00:00 2001 From: Florian Weik Date: Fri, 18 Oct 2019 19:28:30 +0200 Subject: [PATCH 4/7] core: Return torque as value from rotational langevin thermo + modernization of torque trafo --- src/core/rotation.cpp | 23 +++++------------------ src/core/thermostat.hpp | 25 ++++++------------------- 2 files changed, 11 insertions(+), 37 deletions(-) diff --git a/src/core/rotation.cpp b/src/core/rotation.cpp index 3a8c4841c17..3e55725e738 100644 --- a/src/core/rotation.cpp +++ b/src/core/rotation.cpp @@ -168,25 +168,12 @@ void propagate_omega_quat_particle(Particle &p) { } inline void convert_torque_to_body_frame_apply_fix_and_thermostat(Particle &p) { - auto const t = convert_vector_space_to_body(p, p.f.torque); - p.f.torque = Utils::Vector3d{}; + auto const torque = (thermo_switch & THERMO_LANGEVIN) + ? friction_thermo_langevin_rotation(p) + + convert_vector_space_to_body(p, p.f.torque) + : convert_vector_space_to_body(p, p.f.torque); - if (thermo_switch & THERMO_LANGEVIN) { - friction_thermo_langevin_rotation(p); - - p.f.torque += t; - } else { - p.f.torque = t; - } - - if (!(p.p.rotation & ROTATION_X)) - p.f.torque[0] = 0; - - if (!(p.p.rotation & ROTATION_Y)) - p.f.torque[1] = 0; - - if (!(p.p.rotation & ROTATION_Z)) - p.f.torque[2] = 0; + p.f.torque = mask(p.p.rotation, torque); } /** convert the torques to the body-fixed frames and propagate angular diff --git a/src/core/thermostat.hpp b/src/core/thermostat.hpp index ca77f10a04c..82f42a0e202 100644 --- a/src/core/thermostat.hpp +++ b/src/core/thermostat.hpp @@ -261,7 +261,7 @@ inline Utils::Vector3d friction_thermo_langevin(Particle const &p) { \xi_i\f$. The same friction coefficient \f$\gamma\f$ is used as that for translation. */ -inline void friction_thermo_langevin_rotation(Particle &p) { +inline Utils::Vector3d friction_thermo_langevin_rotation(const Particle &p) { extern Thermostat::GammaType langevin_pref2_rotation; Thermostat::GammaType langevin_pref_friction_buf, langevin_pref_noise_buf; @@ -297,28 +297,15 @@ inline void friction_thermo_langevin_rotation(Particle &p) { } #endif /* LANGEVIN_PER_PARTICLE */ - // Rotational degrees of virtual sites are thermostatted, - // so no switching here - // Here the thermostats happens - Utils::Vector3d noise = v_noise(p.p.identity); - for (int j = 0; j < 3; j++) { + auto const noise = v_noise(p.p.identity); #ifdef PARTICLE_ANISOTROPY - if (langevin_pref_noise_buf[j] > 0.0) { - p.f.torque[j] = -langevin_pref_friction_buf[j] * p.m.omega[j] + - langevin_pref_noise_buf[j] * noise[j]; - } else { - p.f.torque[j] = -langevin_pref_friction_buf[j] * p.m.omega[j]; - } + return -hadamard_product(langevin_pref_friction_buf, p.m.omega) + + hadamard_product(langevin_pref_noise_buf, noise); #else - if (langevin_pref_noise_buf > 0.0) { - p.f.torque[j] = -langevin_pref_friction_buf * p.m.omega[j] + - langevin_pref_noise_buf * noise[j]; - } else { - p.f.torque[j] = -langevin_pref_friction_buf * p.m.omega[j]; - } + return -langevin_pref_friction_buf * p.m.omega + + langevin_pref_noise_buf * noise; #endif - } } #endif // ROTATION From d25e32c50feadcedddad53ef280c01823b0bbfc3 Mon Sep 17 00:00:00 2001 From: Florian Weik Date: Fri, 18 Oct 2019 19:41:07 +0200 Subject: [PATCH 5/7] core: rotation: Modernize local_rotate_particle --- src/core/rotation.cpp | 30 +++++++++--------------------- src/utils/include/utils/mask.hpp | 3 ++- src/utils/tests/mask_test.cpp | 2 +- 3 files changed, 12 insertions(+), 23 deletions(-) diff --git a/src/core/rotation.cpp b/src/core/rotation.cpp index 3e55725e738..258806422c9 100644 --- a/src/core/rotation.cpp +++ b/src/core/rotation.cpp @@ -273,31 +273,19 @@ Utils::Vector3d convert_vector_space_to_body(const Particle &p, */ void local_rotate_particle(Particle &p, const Utils::Vector3d &axis_space_frame, const double phi) { - // Convert rotation axis to body-fixed frame - Utils::Vector3d axis = convert_vector_space_to_body(p, axis_space_frame); - // Rotation turned off entirely? if (!p.p.rotation) return; - // Per coordinate fixing - if (!(p.p.rotation & ROTATION_X)) - axis[0] = 0; - if (!(p.p.rotation & ROTATION_Y)) - axis[1] = 0; - if (!(p.p.rotation & ROTATION_Z)) - axis[2] = 0; - // Re-normalize rotation axis - double l = axis.norm(); - // Check, if the rotation axis is nonzero - if (l < std::numeric_limits::epsilon()) - return; - - axis /= l; - - double s = sin(phi / 2); - Utils::Vector4d q = {cos(phi / 2), s * axis[0], s * axis[1], s * axis[2]}; - q.normalize(); + // Convert rotation axis to body-fixed frame + auto const axis = + mask(p.p.rotation, convert_vector_space_to_body(p, axis_space_frame)) + .normalize(); + + auto const s = std::sin(phi / 2); + auto const q = + Utils::Vector4d{cos(phi / 2), s * axis[0], s * axis[1], s * axis[2]} + .normalize(); // Rotate the particle p.r.quat = Utils::multiply_quaternions(p.r.quat, q); diff --git a/src/utils/include/utils/mask.hpp b/src/utils/include/utils/mask.hpp index 0e3da15eebc..eae80f141cd 100644 --- a/src/utils/include/utils/mask.hpp +++ b/src/utils/include/utils/mask.hpp @@ -30,7 +30,8 @@ std::enable_if_t::value && (8 * sizeof(Integral) >= tuple_size::value), T> mask(Integral mask, T t) { - return detail::mask_impl(mask, t, std::make_index_sequence::value>{}); + return detail::mask_impl(mask, t, + std::make_index_sequence::value>{}); } } // namespace Utils diff --git a/src/utils/tests/mask_test.cpp b/src/utils/tests/mask_test.cpp index f01d80cd65b..a44343c4476 100644 --- a/src/utils/tests/mask_test.cpp +++ b/src/utils/tests/mask_test.cpp @@ -23,9 +23,9 @@ #include -#include #include #include +#include BOOST_AUTO_TEST_CASE(mask_) { using Utils::get; From bc07c6432a551811259beb954dadc2033088d1b5 Mon Sep 17 00:00:00 2001 From: Florian Weik Date: Sat, 19 Oct 2019 22:15:12 +0200 Subject: [PATCH 6/7] utils: Improve docstring of mask --- src/utils/include/utils/mask.hpp | 22 +++++++++++++++------- 1 file changed, 15 insertions(+), 7 deletions(-) diff --git a/src/utils/include/utils/mask.hpp b/src/utils/include/utils/mask.hpp index eae80f141cd..30ff3a970eb 100644 --- a/src/utils/include/utils/mask.hpp +++ b/src/utils/include/utils/mask.hpp @@ -16,20 +16,28 @@ auto mask_impl(Integral mask, T t, std::index_sequence) { } // namespace detail /** - * @brief Set elements of a tuple-like to zero by a bit mask. + * @brief Pick elements of a tuple-like by a bit mask. + * + * E.g. every element of the input for which the corresponding + * bit is set in the mask is set is copied to the output unmodified, + * the elements that are not set are set to zero (default constructed + * instance of the type). + * + * Example: + * mask(0b1011, {1, 2, 3, 4}) => {1, 0, 3, 4} * * @tparam T implements the tuple interface(get, tuple_size, ...) * @tparam Integral An unsigned integral type * @param mask bit mask, if the i-th bit is set, the i-th element - * in @param t is set to zero. - * @param t + * in @param t is copied to the output, otherwise it is set to zero. + * @param t input elements * @return t partially zeroed out according to mask */ template -std::enable_if_t::value && - (8 * sizeof(Integral) >= tuple_size::value), - T> -mask(Integral mask, T t) { +auto mask(Integral mask, T t) + -> std::enable_if_t::value && + (8 * sizeof(Integral) >= tuple_size::value), + T> { return detail::mask_impl(mask, t, std::make_index_sequence::value>{}); } From a48b99dfe9fffcf2de70d881079fe71cd06b0acf Mon Sep 17 00:00:00 2001 From: Florian Weik Date: Sun, 20 Oct 2019 12:43:58 +0200 Subject: [PATCH 7/7] utils: size_in_bits + test --- src/utils/include/utils/mask.hpp | 5 ++-- src/utils/include/utils/type_traits.hpp | 8 ++++++ src/utils/tests/CMakeLists.txt | 1 + src/utils/tests/type_traits_test.cpp | 36 +++++++++++++++++++++++++ 4 files changed, 48 insertions(+), 2 deletions(-) create mode 100644 src/utils/tests/type_traits_test.cpp diff --git a/src/utils/include/utils/mask.hpp b/src/utils/include/utils/mask.hpp index 30ff3a970eb..64bc1c78985 100644 --- a/src/utils/include/utils/mask.hpp +++ b/src/utils/include/utils/mask.hpp @@ -3,8 +3,8 @@ #include #include +#include -#include #include namespace Utils { @@ -36,7 +36,8 @@ auto mask_impl(Integral mask, T t, std::index_sequence) { template auto mask(Integral mask, T t) -> std::enable_if_t::value && - (8 * sizeof(Integral) >= tuple_size::value), + (size_in_bits::value >= + tuple_size::value), T> { return detail::mask_impl(mask, t, std::make_index_sequence::value>{}); diff --git a/src/utils/include/utils/type_traits.hpp b/src/utils/include/utils/type_traits.hpp index 5b947a0c6d8..99ebfe5e587 100644 --- a/src/utils/include/utils/type_traits.hpp +++ b/src/utils/include/utils/type_traits.hpp @@ -19,6 +19,7 @@ #ifndef UTILS_TYPE_TRAITS_HPP #define UTILS_TYPE_TRAITS_HPP +#include #include namespace Utils { @@ -44,6 +45,13 @@ template struct conjunction : std::conditional, B1>::type {}; +/** + * @brief Size of a type in bits. + * + * On posix plattforms this is 8 * sizeof(T). + */ +template +struct size_in_bits : std::integral_constant {}; } // namespace Utils #endif diff --git a/src/utils/tests/CMakeLists.txt b/src/utils/tests/CMakeLists.txt index cf397593f2e..a61da92cd0b 100644 --- a/src/utils/tests/CMakeLists.txt +++ b/src/utils/tests/CMakeLists.txt @@ -38,6 +38,7 @@ unit_test(NAME coordinate_transformation_test SRC coordinate_transformation.cpp unit_test(NAME rotation_matrix_test SRC rotation_matrix_test.cpp DEPENDS utils) unit_test(NAME quaternion_test SRC quaternion_test.cpp DEPENDS utils) unit_test(NAME mask_test SRC mask_test.cpp DEPENDS utils) +unit_test(NAME type_traits_test SRC type_traits_test.cpp DEPENDS utils) unit_test(NAME gather_buffer_test SRC gather_buffer_test.cpp DEPENDS utils Boost::mpi MPI::MPI_CXX NUM_PROC 4) unit_test(NAME scatter_buffer_test SRC scatter_buffer_test.cpp DEPENDS utils Boost::mpi MPI::MPI_CXX NUM_PROC 4) diff --git a/src/utils/tests/type_traits_test.cpp b/src/utils/tests/type_traits_test.cpp new file mode 100644 index 00000000000..30a8c5e37cf --- /dev/null +++ b/src/utils/tests/type_traits_test.cpp @@ -0,0 +1,36 @@ +/* + * Copyright (C) 2018-2019 The ESPResSo project + * + * This file is part of ESPResSo. + * + * ESPResSo 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 3 of the License, or + * (at your option) any later version. + * + * ESPResSo 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, see . + */ + +/** \file + * Unit test for Utils Type Traits. + * + */ + +#define BOOST_TEST_MODULE type traits tests +#define BOOST_TEST_DYN_LINK +#include + +#include + +#include + +BOOST_AUTO_TEST_CASE(size_in_bits) { + static_assert(CHAR_BIT == Utils::size_in_bits::value, ""); + static_assert(CHAR_BIT * sizeof(int) == Utils::size_in_bits::value, ""); +}