Line data    Source code

       1           0 : // Distributed under the MIT License.
       2             : // See LICENSE.txt for details.
       3             : 
       4             : #pragma once
       5             : 
       6             : #include <array>
       7             : #include <boost/math/quaternion.hpp>
       8             : #include <cstddef>
       9             : #include <map>
      10             : #include <memory>
      11             : #include <ostream>
      12             : #include <pup.h>
      13             : 
      14             : #include "DataStructures/DataVector.hpp"
      15             : #include "Domain/FunctionsOfTime/FunctionOfTime.hpp"
      16             : #include "Domain/FunctionsOfTime/PiecewisePolynomial.hpp"
      17             : #include "Domain/FunctionsOfTime/ThreadsafeList.hpp"
      18             : #include "Utilities/Gsl.hpp"
      19             : #include "Utilities/Serialization/CharmPupable.hpp"
      20             : 
      21             : namespace domain::FunctionsOfTime {
      22             : /*!
      23             :  * \ingroup ComputationalDomainGroup
      24             :  * \brief A FunctionOfTime that stores quaternions for the rotation map
      25             :  *
      26             :  * \details This FunctionOfTime stores quaternions that will be used in the
      27             :  * time-dependent rotation map as well as the orbital angular velocity that
      28             :  * will be controlled by the rotation control sytem. To get the quaternion, an
      29             :  * ODE is solved of the form \f$ \dot{q} = \frac{1}{2} q \times \omega \f$
      30             :  * where \f$ \omega \f$ is the orbital angular velocity which is stored
      31             :  * internally as the derivative of an angle `PiecewisePolynomial`, and
      32             :  * \f$ \times \f$ here is quaternion multiplication.
      33             :  *
      34             :  * Different from a `PiecewisePolynomial`, only the quaternion
      35             :  * itself is stored, not any of the derivatives because the derivatives must be
      36             :  * calculated from the solved ODE at every function call. Because
      37             :  * derivatives of the quaternion are not stored, the template parameter
      38             :  * `MaxDeriv` refers to both the max derivative of the stored angle
      39             :  * PiecewisePolynomial and the max derivative returned by the
      40             :  * QuaternionFunctionOfTime. The `update` function is then just a wrapper
      41             :  * around the internal `PiecewisePolynomial::update` function with the addition
      42             :  * that it then updates the stored quaternions as well.
      43             :  *
      44             :  * The angle PiecewisePolynomial is accessible through the `angle_func`,
      45             :  * `angle_func_and_deriv`, and `angle_func_and_2_derivs` functions which
      46             :  * correspond to the function calls of a normal PiecewisePolynomial except
      47             :  * without the `angle_` prefix.
      48             :  *
      49             :  * It is encouraged to use `quat_func` and `angle_func` when you want the
      50             :  * specific values of the functions to avoid ambiguity in what you are
      51             :  * calling. However, the original three `func` functions inherited from the
      52             :  * FunctionOfTime base class are necessary because the maps use the generic
      53             :  * `func` functions, thus they return the quaternion and its derivatives (which
      54             :  * are needed for the map). This is all to keep the symmetry of naming
      55             :  * `angle_func` and `quat_func` so that function calls won't be ambiguous.
      56             :  *
      57             :  * \note This class conforms to the requirements of the
      58             :  * `Parallel::GlobalCache` for objects held by mutable global cache tags.
      59             :  */
      60             : template <size_t MaxDeriv>
      61           1 : class QuaternionFunctionOfTime : public FunctionOfTime {
      62             :  public:
      63           0 :   QuaternionFunctionOfTime();
      64           0 :   QuaternionFunctionOfTime(QuaternionFunctionOfTime&&);
      65           0 :   QuaternionFunctionOfTime(const QuaternionFunctionOfTime&);
      66           0 :   QuaternionFunctionOfTime& operator=(QuaternionFunctionOfTime&&);
      67           0 :   QuaternionFunctionOfTime& operator=(const QuaternionFunctionOfTime&);
      68           0 :   ~QuaternionFunctionOfTime() override;
      69             : 
      70           0 :   QuaternionFunctionOfTime(
      71             :       double t, const std::array<DataVector, 1>& initial_quat_func,
      72             :       std::array<DataVector, MaxDeriv + 1> initial_angle_func,
      73             :       double expiration_time);
      74             : 
      75             :   // LCOV_EXCL_START
      76           0 :   explicit QuaternionFunctionOfTime(CkMigrateMessage* /*unused*/);
      77             :   // LCOV_EXCL_STOP
      78             : 
      79           0 :   auto get_clone() const -> std::unique_ptr<FunctionOfTime> override;
      80             : 
      81             :   // clang-tidy: google-runtime-references
      82             :   // clang-tidy: cppcoreguidelines-owning-memory,-warnings-as-errors
      83           0 :   WRAPPED_PUPable_decl_template(QuaternionFunctionOfTime<MaxDeriv>);  // NOLINT
      84             : 
      85             :   /// Returns domain of validity for the function of time
      86           1 :   std::array<double, 2> time_bounds() const override;
      87             : 
      88           1 :   double expiration_after(double time) const override;
      89             : 
      90             :   /// Updates the `MaxDeriv`th derivative of the angle piecewisepolynomial at
      91             :   /// the given time, then updates the stored quaternions.
      92             :   ///
      93             :   /// `updated_max_deriv` is a datavector of the `MaxDeriv`s for each component.
      94             :   /// `next_expiration_time` is the next expiration time.
      95             :   ///
      96             :   /// The \p time_of_update must be the same as the old expiration
      97             :   /// time.  It is passed as a check that the calling code is
      98             :   /// computing the other arguments with the correct value.
      99           1 :   void update(double time_of_update, DataVector updated_max_deriv,
     100             :               double next_expiration_time) override;
     101             : 
     102             :   // NOLINTNEXTLINE(google-runtime-references)
     103           0 :   void pup(PUP::er& p) override;
     104             : 
     105             :   /// Returns the quaternion at an arbitrary time `t`.
     106           1 :   std::array<DataVector, 1> func(const double t) const override {
     107             :     return quat_func(t);
     108             :   }
     109             : 
     110             :   /// Returns the quaternion and its first derivative at an arbitrary time `t`.
     111           1 :   std::array<DataVector, 2> func_and_deriv(const double t) const override {
     112             :     return quat_func_and_deriv(t);
     113             :   }
     114             : 
     115             :   /// Returns the quaternion and the first two derivatives at an arbitrary
     116             :   /// time `t`.
     117           1 :   std::array<DataVector, 3> func_and_2_derivs(const double t) const override {
     118             :     return quat_func_and_2_derivs(t);
     119             :   }
     120             : 
     121             :   /// Returns the quaternion at an arbitrary time `t`.
     122           1 :   std::array<DataVector, 1> quat_func(double t) const;
     123             : 
     124             :   /// Returns the quaternion and its first derivative at an arbitrary time `t`.
     125           1 :   std::array<DataVector, 2> quat_func_and_deriv(double t) const;
     126             : 
     127             :   /// Returns the quaternion and the first two derivatives at an arbitrary
     128             :   /// time `t`.
     129           1 :   std::array<DataVector, 3> quat_func_and_2_derivs(double t) const;
     130             : 
     131             :   /// Returns stored angle at an arbitrary time `t`.
     132           1 :   std::array<DataVector, 1> angle_func(const double t) const {
     133             :     return angle_f_of_t_.func(t);
     134             :   }
     135             : 
     136             :   /// Returns stored angle and its first derivative (omega) at an arbitrary time
     137             :   /// `t`.
     138           1 :   std::array<DataVector, 2> angle_func_and_deriv(const double t) const {
     139             :     return angle_f_of_t_.func_and_deriv(t);
     140             :   }
     141             : 
     142             :   /// Returns stored angle and the first two derivatives at an arbitrary
     143             :   /// time `t`.
     144           1 :   std::array<DataVector, 3> angle_func_and_2_derivs(const double t) const {
     145             :     return angle_f_of_t_.func_and_2_derivs(t);
     146             :   }
     147             : 
     148             :  private:
     149             :   template <size_t LocalMaxDeriv>
     150           0 :   friend bool operator==(  // NOLINT(readability-redundant-declaration)
     151             :       const QuaternionFunctionOfTime<LocalMaxDeriv>& lhs,
     152             :       const QuaternionFunctionOfTime<LocalMaxDeriv>& rhs);
     153             : 
     154             :   template <size_t LocalMaxDeriv>
     155           0 :   friend std::ostream& operator<<(  // NOLINT(readability-redundant-declaration)
     156             :       std::ostream& os,
     157             :       const QuaternionFunctionOfTime<LocalMaxDeriv>& quaternion_f_of_t);
     158             : 
     159             :   FunctionOfTimeHelpers::ThreadsafeList<boost::math::quaternion<double>>
     160           0 :       stored_quaternions_and_times_{};
     161           0 :   domain::FunctionsOfTime::PiecewisePolynomial<MaxDeriv> angle_f_of_t_{};
     162           0 :   std::map<double, double> update_backlog_{};
     163             : 
     164           0 :   void unpack_old_version(PUP::er& p, size_t version);
     165             : 
     166             :   /// Integrates the ODE \f$ \dot{q} = \frac{1}{2} q \times \omega \f$ from time
     167             :   /// `t0` to time `t`. On input, `quaternion_to_integrate` is the initial
     168             :   /// quaternion at time `t0` and on output, it stores the result at time `t`
     169           1 :   void solve_quaternion_ode(
     170             :       gsl::not_null<boost::math::quaternion<double>*> quaternion_to_integrate,
     171             :       double t0, double t) const;
     172             : 
     173             :   /// Does common operations to all the `func` functions such as updating stored
     174             :   /// info, solving the ODE, and returning the normalized quaternion as a boost
     175             :   /// quaternion for easy calculations
     176           1 :   boost::math::quaternion<double> setup_func(double t) const;
     177             : };
     178             : 
     179             : template <size_t MaxDeriv>
     180           0 : bool operator!=(const QuaternionFunctionOfTime<MaxDeriv>& lhs,
     181             :                 const QuaternionFunctionOfTime<MaxDeriv>& rhs);
     182             : 
     183             : /// \cond
     184             : template <size_t MaxDeriv>
     185             : PUP::able::PUP_ID QuaternionFunctionOfTime<MaxDeriv>::my_PUP_ID = 0;  // NOLINT
     186             : /// \endcond
     187             : }  // namespace domain::FunctionsOfTime