Line data    Source code

       1           0 : // Distributed under the MIT License.
       2             : // See LICENSE.txt for details.
       3             : 
       4             : #pragma once
       5             : 
       6             : #include <algorithm>
       7             : #include <boost/container/static_vector.hpp>
       8             : #include <cstddef>
       9             : #include <optional>
      10             : #include <ostream>
      11             : #include <pup.h>
      12             : #include <pup_stl.h>
      13             : #include <type_traits>
      14             : #include <utility>
      15             : 
      16             : #include "DataStructures/CircularDeque.hpp"
      17             : #include "DataStructures/MathWrapper.hpp"
      18             : #include "DataStructures/StaticDeque.hpp"
      19             : #include "Time/History.hpp"
      20             : #include "Time/TimeStepId.hpp"
      21             : #include "Utilities/Algorithm.hpp"
      22             : #include "Utilities/ErrorHandling/Assert.hpp"
      23             : #include "Utilities/Gsl.hpp"
      24             : #include "Utilities/Literals.hpp"
      25             : #include "Utilities/Serialization/PupBoost.hpp"
      26             : #include "Utilities/StdHelpers.hpp"
      27             : #include "Utilities/StlBoilerplate.hpp"
      28             : #include "Utilities/TMPL.hpp"
      29             : 
      30           1 : namespace TimeSteppers {
      31             : 
      32             : /// \ingroup TimeSteppersGroup
      33             : /// Access to the list of `TimeStepId`s in a `BoundaryHistory`.
      34             : ///
      35             : /// For simplicity of implementation, iterable-container access is not
      36             : /// provided for substeps within a step, but is instead provided
      37             : /// through additional methods on this class.
      38             : /// @{
      39           1 : class ConstBoundaryHistoryTimes
      40             :     : public stl_boilerplate::RandomAccessSequence<ConstBoundaryHistoryTimes,
      41             :                                                    const TimeStepId, false> {
      42             :  protected:
      43           0 :   ~ConstBoundaryHistoryTimes() = default;
      44             : 
      45             :  public:
      46           0 :   virtual size_t size() const = 0;
      47           0 :   virtual const TimeStepId& operator[](size_t n) const = 0;
      48           0 :   virtual const TimeStepId& operator[](
      49             :       const std::pair<size_t, size_t>& step_and_substep) const = 0;
      50           0 :   virtual size_t integration_order(size_t n) const = 0;
      51           0 :   virtual size_t integration_order(const TimeStepId& id) const = 0;
      52           0 :   virtual size_t number_of_substeps(size_t n) const = 0;
      53             :   /// This returns the same value for any substep of the same step.
      54           1 :   virtual size_t number_of_substeps(const TimeStepId& id) const = 0;
      55             : };
      56             : 
      57           0 : class MutableBoundaryHistoryTimes : public ConstBoundaryHistoryTimes {
      58             :  protected:
      59           0 :   ~MutableBoundaryHistoryTimes() = default;
      60             : 
      61             :  public:
      62             :   /// Remove the earliest step and its substeps.
      63           1 :   virtual void pop_front() const = 0;
      64           0 :   virtual void clear() const = 0;
      65             :   /// Remove all substeps for step \p n except for the step itself.
      66           1 :   virtual void clear_substeps(size_t n) const = 0;
      67             : };
      68             : /// @}
      69             : 
      70             : /// \ingroup TimeSteppersGroup
      71             : /// Type erased base class for evaluating BoundaryHistory couplings.
      72             : ///
      73             : /// The results are cached in the `BoundaryHistory` class.
      74             : template <typename UntypedCouplingResult>
      75           1 : class BoundaryHistoryEvaluator {
      76             :  public:
      77           0 :   virtual MathWrapper<const UntypedCouplingResult> operator()(
      78             :       const TimeStepId& local_id, const TimeStepId& remote_id) const = 0;
      79             : 
      80             :  protected:
      81           0 :   ~BoundaryHistoryEvaluator() = default;
      82             : };
      83             : 
      84             : /// \ingroup TimeSteppersGroup
      85             : /// History data used by a TimeStepper for boundary integration.
      86             : ///
      87             : /// \tparam LocalData local data passed to the boundary coupling
      88             : /// \tparam RemoteData remote data passed to the boundary coupling
      89             : /// \tparam CouplingResult type of cached boundary couplings
      90             : template <typename LocalData, typename RemoteData, typename CouplingResult>
      91           1 : class BoundaryHistory {
      92             :  public:
      93           0 :   BoundaryHistory() = default;
      94           0 :   BoundaryHistory(const BoundaryHistory& other) = default;
      95           0 :   BoundaryHistory(BoundaryHistory&&) = default;
      96           0 :   BoundaryHistory& operator=(const BoundaryHistory& other) = default;
      97           0 :   BoundaryHistory& operator=(BoundaryHistory&&) = default;
      98           0 :   ~BoundaryHistory() = default;
      99             : 
     100             :   /// The wrapped types presented by the type-erased history.  One of
     101             :   /// the types in \ref MATH_WRAPPER_TYPES.
     102           1 :   using UntypedCouplingResult = math_wrapper_type<CouplingResult>;
     103             : 
     104             :   // Factored out of ConstSideAccess so that the base classes of
     105             :   // MutableSideAccess can have protected destructors.
     106             :   template <bool Local, bool Mutable>
     107           0 :   class SideAccessCommon
     108             :       : public tmpl::conditional_t<Mutable, MutableBoundaryHistoryTimes,
     109             :                                    ConstBoundaryHistoryTimes> {
     110             :    public:
     111           0 :     using MutableData = tmpl::conditional_t<Local, LocalData, RemoteData>;
     112           0 :     using Data = tmpl::conditional_t<Mutable, MutableData, const MutableData>;
     113             : 
     114           0 :     size_t size() const override { return parent_data().size(); }
     115           0 :     static constexpr size_t max_size() {
     116             :       return decltype(std::declval<ConstSideAccess>()
     117             :                           .parent_data())::max_size();
     118             :     }
     119             : 
     120           0 :     const TimeStepId& operator[](const size_t n) const override {
     121             :       return (*this)[{n, 0}];
     122             :     }
     123           0 :     const TimeStepId& operator[](
     124             :         const std::pair<size_t, size_t>& step_and_substep) const override {
     125             :       return entry(step_and_substep).id;
     126             :     }
     127             : 
     128           0 :     size_t integration_order(const size_t n) const override {
     129             :       return parent_data()[n].integration_order;
     130             :     }
     131           0 :     size_t integration_order(const TimeStepId& id) const override {
     132             :       return step_data(id).integration_order;
     133             :     }
     134             : 
     135           0 :     size_t number_of_substeps(const size_t n) const override {
     136             :       return parent_data()[n].substeps.size();
     137             :     }
     138           0 :     size_t number_of_substeps(const TimeStepId& id) const override {
     139             :       return step_data(id).substeps.size();
     140             :     }
     141             : 
     142             :     /// Access the data stored on the side.  When performed through a
     143             :     /// `MutableSideAccess`, these allow modification of the data.
     144             :     /// Performing such modifications likely invalidates the coupling
     145             :     /// cache for the associated `BoundaryHistory` object, which
     146             :     /// should be cleared.
     147             :     /// @{
     148           1 :     Data& data(const size_t n) const {
     149             :       return parent_data()[n].substeps.front().data;
     150             :     }
     151           1 :     Data& data(const TimeStepId& id) const {
     152             :       return entry(id).data;
     153             :     }
     154             :     /// @}
     155             : 
     156             :     /// Apply \p func to each entry.
     157             :     ///
     158             :     /// The function \p func must accept two arguments, one of type
     159             :     /// `const TimeStepId&` and a second of either type `const Data&`
     160             :     /// or `gsl::not_null<Data*>`.  (Note that `Data` may be a
     161             :     /// const-qualified type.)  If entries are modified, the coupling
     162             :     /// cache must be cleared by calling `clear_coupling_cache()` on
     163             :     /// the parent `BoundaryHistory` object.
     164             :     template <typename Func>
     165           1 :     void for_each(Func&& func) const;
     166             : 
     167             :    protected:
     168           0 :     ~SideAccessCommon() = default;
     169             : 
     170           0 :     auto& parent_data() const {
     171             :       if constexpr (Local) {
     172             :         return parent_->local_data_;
     173             :       } else {
     174             :         return parent_->remote_data_;
     175             :       }
     176             :     }
     177             : 
     178           0 :     auto& step_data(const TimeStepId& id) const {
     179             :       auto entry =
     180             :           std::upper_bound(parent_data().begin(), parent_data().end(), id);
     181             :       ASSERT(entry != parent_data().begin(), "Id " << id << " not present.");
     182             :       --entry;
     183             :       ASSERT(id.substep() < entry->substeps.size() and
     184             :              entry->substeps[id.substep()].id == id,
     185             :              "Id " << id << " not present.");
     186             :       return *entry;
     187             :     }
     188             : 
     189           0 :     auto& entry(const TimeStepId& id) const {
     190             :       // Bounds and consistency are checked in step_data()
     191             :       return step_data(id).substeps[id.substep()];
     192             :     }
     193             : 
     194           0 :     auto& entry(const std::pair<size_t, size_t>& step_and_substep) const {
     195             :       ASSERT(step_and_substep.first < parent_data().size(),
     196             :              "Step out of range: " << step_and_substep.first);
     197             :       auto& substeps = parent_data()[step_and_substep.first].substeps;
     198             :       ASSERT(step_and_substep.second < substeps.size(),
     199             :              "Substep out of range: " << step_and_substep.second);
     200             :       return substeps[step_and_substep.second];
     201             :     }
     202             : 
     203           0 :     using StoredHistory =
     204             :         tmpl::conditional_t<Mutable, BoundaryHistory, const BoundaryHistory>;
     205           0 :     explicit SideAccessCommon(const gsl::not_null<StoredHistory*> parent)
     206             :         : parent_(parent) {}
     207             : 
     208           0 :     gsl::not_null<StoredHistory*> parent_;
     209             :   };
     210             : 
     211             :   /// \cond
     212             :   template <bool Local>
     213             :   class ConstSideAccess;
     214             :   /// \endcond
     215             : 
     216             :   template <bool Local>
     217           0 :   class MutableSideAccess final : public SideAccessCommon<Local, true> {
     218             :    public:
     219           0 :     using Data = tmpl::conditional_t<Local, LocalData, RemoteData>;
     220             : 
     221           0 :     void pop_front() const override;
     222           0 :     void clear() const override;
     223           0 :     void clear_substeps(size_t n) const override;
     224             : 
     225           0 :     void insert(const TimeStepId& id, size_t integration_order,
     226             :                 Data data) const;
     227             : 
     228           0 :     void insert_initial(const TimeStepId& id, size_t integration_order,
     229             :                         Data data) const;
     230             : 
     231             :    private:
     232           0 :     friend class BoundaryHistory;
     233             :     friend class ConstSideAccess<Local>;
     234           0 :     explicit MutableSideAccess(const gsl::not_null<BoundaryHistory*> parent)
     235             :         : SideAccessCommon<Local, true>(parent) {}
     236             :   };
     237             : 
     238             :   template <bool Local>
     239           0 :   class ConstSideAccess final : public SideAccessCommon<Local, false> {
     240             :    private:
     241           0 :     friend class BoundaryHistory;
     242           0 :     explicit ConstSideAccess(const gsl::not_null<const BoundaryHistory*> parent)
     243             :         : SideAccessCommon<Local, false>(parent) {}
     244             :   };
     245             : 
     246           0 :   MutableSideAccess<true> local() { return MutableSideAccess<true>(this); }
     247           0 :   ConstSideAccess<true> local() const { return ConstSideAccess<true>(this); }
     248             : 
     249           0 :   MutableSideAccess<false> remote() { return MutableSideAccess<false>(this); }
     250           0 :   ConstSideAccess<false> remote() const { return ConstSideAccess<false>(this); }
     251             : 
     252             :  private:
     253             :   template <typename Coupling>
     254           0 :   class EvaluatorImpl final
     255             :       : public BoundaryHistoryEvaluator<UntypedCouplingResult> {
     256             :    public:
     257           0 :     MathWrapper<const UntypedCouplingResult> operator()(
     258             :         const TimeStepId& local_id, const TimeStepId& remote_id) const;
     259             : 
     260             :    private:
     261           0 :     friend class BoundaryHistory;
     262             : 
     263           0 :     EvaluatorImpl(const gsl::not_null<const BoundaryHistory*> parent,
     264             :                   Coupling coupling)
     265             :         : parent_(parent), coupling_(std::move(coupling)) {}
     266             : 
     267           0 :     gsl::not_null<const BoundaryHistory*> parent_;
     268           0 :     Coupling coupling_;
     269             :   };
     270             : 
     271             :  public:
     272             :   /// Obtain an object that can evaluate type-erased boundary
     273             :   /// couplings.
     274             :   ///
     275             :   /// The passed functor must take objects of types `LocalData` and
     276             :   /// `RemoteData` and return an object convertible to
     277             :   /// `CouplingResult`.  Results are cached, so different calls to
     278             :   /// this function should pass equivalent couplings.
     279             :   template <typename Coupling>
     280           1 :   auto evaluator(Coupling&& coupling) const {
     281             :     return EvaluatorImpl<Coupling>(this, std::forward<Coupling>(coupling));
     282             :   }
     283             : 
     284             :   /// Clear the cached values.
     285             :   ///
     286             :   /// This is required after existing history entries that have been
     287             :   /// used in coupling calculations are mutated.
     288           1 :   void clear_coupling_cache();
     289             : 
     290           0 :   void pup(PUP::er& p);
     291             : 
     292             :   template <bool IncludeData>
     293           0 :   std::ostream& print(std::ostream& os, size_t padding_size = 0) const;
     294             : 
     295             :  private:
     296             :   template <typename Data>
     297           0 :   struct StepData {
     298           0 :     struct Entry {
     299           0 :       TimeStepId id;
     300           0 :       Data data;
     301             : 
     302           0 :       void pup(PUP::er& p) {
     303             :         p | id;
     304             :         p | data;
     305             :       }
     306             :     };
     307             : 
     308           0 :     size_t integration_order;
     309             :     // Unlike in History, the full step is the first entry, so we need
     310             :     // one more element.
     311           0 :     boost::container::static_vector<Entry, history_max_substeps + 1> substeps;
     312             : 
     313           0 :     void pup(PUP::er& p) {
     314             :       p | integration_order;
     315             :       p | substeps;
     316             :     }
     317             : 
     318           0 :     friend bool operator<(const StepData& a, const StepData& b) {
     319             :       return a.substeps.front().id < b.substeps.front().id;
     320             :     }
     321           0 :     friend bool operator<(const TimeStepId& a, const StepData& b) {
     322             :       return a < b.substeps.front().id;
     323             :     }
     324           0 :     friend bool operator<(const StepData& a, const TimeStepId& b) {
     325             :       return a.substeps.front().id < b;
     326             :     }
     327             :   };
     328             : 
     329           0 :   void insert_local(const TimeStepId& id, size_t integration_order,
     330             :                     LocalData data);
     331           0 :   void insert_remote(const TimeStepId& id, size_t integration_order,
     332             :                      RemoteData data);
     333             : 
     334           0 :   void insert_initial_local(const TimeStepId& id, size_t integration_order,
     335             :                             LocalData data);
     336           0 :   void insert_initial_remote(const TimeStepId& id, size_t integration_order,
     337             :                              RemoteData data);
     338             : 
     339           0 :   void pop_local();
     340           0 :   void pop_remote();
     341             : 
     342           0 :   void clear_substeps_local(size_t n);
     343           0 :   void clear_substeps_remote(size_t n);
     344             : 
     345           0 :   StaticDeque<StepData<LocalData>, history_max_past_steps + 2> local_data_{};
     346           0 :   CircularDeque<StepData<RemoteData>> remote_data_{};
     347             : 
     348             :   template <typename Data>
     349           0 :   using CouplingSubsteps =
     350             :       boost::container::static_vector<Data, history_max_substeps + 1>;
     351             : 
     352             :   // Putting the CircularDeque outermost means that we are inserting
     353             :   // and removing containers that do not allocate, so we don't have to
     354             :   // worry about that.
     355             :   // NOLINTNEXTLINE(spectre-mutable)
     356             :   mutable CircularDeque<CouplingSubsteps<
     357             :       StaticDeque<CouplingSubsteps<std::optional<CouplingResult>>,
     358             :                   decltype(local_data_)::max_size()>>>
     359           0 :       couplings_;
     360             : };
     361             : 
     362             : template <typename LocalData, typename RemoteData, typename CouplingResult>
     363             : template <bool Local, bool Mutable>
     364             : template <typename Func>
     365             : void BoundaryHistory<LocalData, RemoteData, CouplingResult>::SideAccessCommon<
     366             :     Local, Mutable>::for_each(Func&& func) const {
     367             :   for (auto& step : parent_data()) {
     368             :     for (auto& substep : step.substeps) {
     369             :       if constexpr (std::is_invocable_v<Func&, const TimeStepId&,
     370             :                                         const Data&>) {
     371             :         func(std::as_const(substep.id), std::as_const(substep.data));
     372             :       } else {
     373             :         func(std::as_const(substep.id), make_not_null(&substep.data));
     374             :       }
     375             :     }
     376             :   }
     377             : }
     378             : 
     379             : template <typename LocalData, typename RemoteData, typename CouplingResult>
     380             : template <bool Local>
     381             : void BoundaryHistory<LocalData, RemoteData, CouplingResult>::MutableSideAccess<
     382             :     Local>::pop_front() const {
     383             :   if constexpr (Local) {
     384             :     this->parent_->pop_local();
     385             :   } else {
     386             :     this->parent_->pop_remote();
     387             :   }
     388             : }
     389             : 
     390             : template <typename LocalData, typename RemoteData, typename CouplingResult>
     391             : template <bool Local>
     392             : void BoundaryHistory<LocalData, RemoteData,
     393             :                      CouplingResult>::MutableSideAccess<Local>::clear() const {
     394             :   while (not this->empty()) {
     395             :     pop_front();
     396             :   }
     397             : }
     398             : 
     399             : template <typename LocalData, typename RemoteData, typename CouplingResult>
     400             : template <bool Local>
     401             : void BoundaryHistory<LocalData, RemoteData, CouplingResult>::MutableSideAccess<
     402             :     Local>::clear_substeps(const size_t n) const {
     403             :   if constexpr (Local) {
     404             :     this->parent_->clear_substeps_local(n);
     405             :   } else {
     406             :     this->parent_->clear_substeps_remote(n);
     407             :   }
     408             : }
     409             : 
     410             : template <typename LocalData, typename RemoteData, typename CouplingResult>
     411             : template <bool Local>
     412             : void BoundaryHistory<LocalData, RemoteData, CouplingResult>::MutableSideAccess<
     413             :     Local>::insert(const TimeStepId& id, const size_t integration_order,
     414             :                    Data data) const {
     415             :   ASSERT(this->parent_data().empty() or
     416             :              id > this->parent_data().back().substeps.back().id,
     417             :          "New data not newer than current data.");
     418             :   if constexpr (Local) {
     419             :     this->parent_->insert_local(id, integration_order, std::move(data));
     420             :   } else {
     421             :     this->parent_->insert_remote(id, integration_order, std::move(data));
     422             :   }
     423             : }
     424             : 
     425             : template <typename LocalData, typename RemoteData, typename CouplingResult>
     426             : template <bool Local>
     427             : void BoundaryHistory<LocalData, RemoteData, CouplingResult>::MutableSideAccess<
     428             :     Local>::insert_initial(const TimeStepId& id, const size_t integration_order,
     429             :                            Data data) const {
     430             :   ASSERT(id.substep() == 0, "Cannot insert_initial with substeps.");
     431             :   ASSERT(this->parent_data().empty() or
     432             :              id < this->parent_data().front().substeps.front().id,
     433             :          "New data not older than current data.");
     434             :   if constexpr (Local) {
     435             :     this->parent_->insert_initial_local(id, integration_order, std::move(data));
     436             :   } else {
     437             :     this->parent_->insert_initial_remote(id, integration_order,
     438             :                                          std::move(data));
     439             :   }
     440             : }
     441             : 
     442             : template <typename LocalData, typename RemoteData, typename CouplingResult>
     443             : template <typename Coupling>
     444             : auto BoundaryHistory<LocalData, RemoteData, CouplingResult>::EvaluatorImpl<
     445             :     Coupling>::operator()(const TimeStepId& local_id,
     446             :                           const TimeStepId& remote_id) const
     447             :     -> MathWrapper<const UntypedCouplingResult> {
     448             :   auto local_entry = std::upper_bound(
     449             :       parent_->local_data_.begin(), parent_->local_data_.end(), local_id);
     450             :   ASSERT(local_entry != parent_->local_data_.begin(), "local_id not present");
     451             :   --local_entry;
     452             :   ASSERT(local_id.substep() < local_entry->substeps.size() and
     453             :          local_entry->substeps[local_id.substep()].id == local_id,
     454             :          "local_id not present");
     455             :   const auto local_step_offset =
     456             :       static_cast<size_t>(local_entry - parent_->local_data_.begin());
     457             : 
     458             :   auto remote_entry = std::upper_bound(
     459             :       parent_->remote_data_.begin(), parent_->remote_data_.end(), remote_id);
     460             :   ASSERT(remote_entry != parent_->remote_data_.begin(),
     461             :          "remote_id not present");
     462             :   --remote_entry;
     463             :   ASSERT(remote_id.substep() < remote_entry->substeps.size() and
     464             :          remote_entry->substeps[remote_id.substep()].id == remote_id,
     465             :          "remote_id not present");
     466             :   const auto remote_step_offset =
     467             :       static_cast<size_t>(remote_entry - parent_->remote_data_.begin());
     468             : 
     469             :   auto& coupling_entry = parent_->couplings_[remote_step_offset]
     470             :                              [remote_id.substep()][local_step_offset]
     471             :                              [local_id.substep()];
     472             :   if (not coupling_entry.has_value()) {
     473             :     coupling_entry.emplace(coupling_(
     474             :         local_entry->substeps[local_id.substep()].data,
     475             :         remote_entry->substeps[remote_id.substep()].data));
     476             :   }
     477             :   return make_math_wrapper(*coupling_entry);
     478             : }
     479             : 
     480             : template <typename LocalData, typename RemoteData, typename CouplingResult>
     481             : void BoundaryHistory<LocalData, RemoteData,
     482             :                      CouplingResult>::clear_coupling_cache() {
     483             :   for (auto& remote_step : couplings_) {
     484             :     for (auto& remote_substep : remote_step) {
     485             :       for (auto& local_step : remote_substep) {
     486             :         for (auto& local_substep : local_step) {
     487             :           local_substep.reset();
     488             :         }
     489             :       }
     490             :     }
     491             :   }
     492             : }
     493             : 
     494             : template <typename LocalData, typename RemoteData, typename CouplingResult>
     495             : void BoundaryHistory<LocalData, RemoteData, CouplingResult>::pup(PUP::er& p) {
     496             :   p | local_data_;
     497             :   p | remote_data_;
     498             :   p | couplings_;
     499             : }
     500             : 
     501             : template <typename LocalData, typename RemoteData, typename CouplingResult>
     502             : template <bool IncludeData>
     503             : std::ostream& BoundaryHistory<LocalData, RemoteData, CouplingResult>::print(
     504             :     std::ostream& os, const size_t padding_size) const {
     505             :   const std::string pad(padding_size, ' ');
     506             :   using ::operator<<;
     507             :   const auto do_print = [&os, &pad](const auto& times) {
     508             :     for (size_t step = 0; step < times.size(); ++step) {
     509             :       const size_t number_of_substeps = times.number_of_substeps(step);
     510             :       for (size_t substep = 0; substep < number_of_substeps; ++substep) {
     511             :         const auto id = times[{step, substep}];
     512             :         os << pad << " Time: " << id;
     513             :         if (substep == 0) {
     514             :           os << " (order " << times.integration_order(step) << ")";
     515             :         }
     516             :         os << "\n";
     517             :         if constexpr (IncludeData) {
     518             :           os << pad << "  Data: ";
     519             :           // os << times.data(id) fails to compile on gcc-11
     520             :           print_stl(os, times.data(id));
     521             :           os << "\n";
     522             :         }
     523             :       }
     524             :     }
     525             :   };
     526             :   os << pad << "Local Data:\n";
     527             :   do_print(local());
     528             :   os << pad << "Remote Data:\n";
     529             :   do_print(remote());
     530             :   return os;
     531             : }
     532             : 
     533             : template <typename LocalData, typename RemoteData, typename CouplingResult>
     534             : void BoundaryHistory<LocalData, RemoteData, CouplingResult>::insert_local(
     535             :     const TimeStepId& id, const size_t integration_order, LocalData data) {
     536             :   if (id.substep() == 0) {
     537             :     local_data_.push_back({integration_order, {}});
     538             :   } else {
     539             :     ASSERT(integration_order == local_data_.back().integration_order,
     540             :            "Cannot change integration order during a step.");
     541             :   }
     542             :   local_data_.back().substeps.push_back({id, std::move(data)});
     543             :   for (auto& remote_step : couplings_) {
     544             :     for (auto& remote_substep : remote_step) {
     545             :       if (id.substep() == 0) {
     546             :         remote_substep.emplace_back(1_st);
     547             :       } else {
     548             :         remote_substep.back().emplace_back();
     549             :       }
     550             :     }
     551             :   }
     552             : }
     553             : 
     554             : template <typename LocalData, typename RemoteData, typename CouplingResult>
     555             : void BoundaryHistory<LocalData, RemoteData, CouplingResult>::insert_remote(
     556             :     const TimeStepId& id, const size_t integration_order, RemoteData data) {
     557             :   if (id.substep() == 0) {
     558             :     remote_data_.push_back({integration_order, {}});
     559             :   } else {
     560             :     ASSERT(integration_order == remote_data_.back().integration_order,
     561             :            "Cannot change integration order during a step.");
     562             :   }
     563             :   remote_data_.back().substeps.push_back({id, std::move(data)});
     564             :   if (id.substep() == 0) {
     565             :     couplings_.emplace_back(1_st);
     566             :   } else {
     567             :     couplings_.back().emplace_back();
     568             :   }
     569             :   for (const auto& local_step : local_data_) {
     570             :     couplings_.back().back().emplace_back(local_step.substeps.size());
     571             :   }
     572             : }
     573             : 
     574             : template <typename LocalData, typename RemoteData, typename CouplingResult>
     575             : void BoundaryHistory<LocalData, RemoteData, CouplingResult>::
     576             :     insert_initial_local(const TimeStepId& id, const size_t integration_order,
     577             :                          LocalData data) {
     578             :   local_data_.push_front({integration_order, {}});
     579             :   local_data_.front().substeps.push_back({id, std::move(data)});
     580             :   for (auto& remote_step : couplings_) {
     581             :     for (auto& remote_substep : remote_step) {
     582             :       remote_substep.emplace_front(1_st);
     583             :     }
     584             :   }
     585             : }
     586             : 
     587             : template <typename LocalData, typename RemoteData, typename CouplingResult>
     588             : void BoundaryHistory<LocalData, RemoteData, CouplingResult>::
     589             :     insert_initial_remote(const TimeStepId& id, const size_t integration_order,
     590             :                           RemoteData data) {
     591             :   remote_data_.push_front({integration_order, {}});
     592             :   remote_data_.front().substeps.push_back({id, std::move(data)});
     593             :   couplings_.emplace_front(1_st);
     594             :   for (const auto& local_step : local_data_) {
     595             :     couplings_.front().back().emplace_back(local_step.substeps.size());
     596             :   }
     597             : }
     598             : 
     599             : template <typename LocalData, typename RemoteData, typename CouplingResult>
     600             : void BoundaryHistory<LocalData, RemoteData, CouplingResult>::pop_local() {
     601             :   local_data_.pop_front();
     602             :   for (auto& remote_step : couplings_) {
     603             :     for (auto& remote_substep : remote_step) {
     604             :       remote_substep.pop_front();
     605             :     }
     606             :   }
     607             : }
     608             : 
     609             : template <typename LocalData, typename RemoteData, typename CouplingResult>
     610             : void BoundaryHistory<LocalData, RemoteData, CouplingResult>::pop_remote() {
     611             :   remote_data_.pop_front();
     612             :   couplings_.pop_front();
     613             : }
     614             : 
     615             : template <typename LocalData, typename RemoteData, typename CouplingResult>
     616             : void BoundaryHistory<LocalData, RemoteData,
     617             :                      CouplingResult>::clear_substeps_local(const size_t n) {
     618             :   local_data_[n].substeps.erase(local_data_[n].substeps.begin() + 1,
     619             :                                 local_data_[n].substeps.end());
     620             :   for (auto& remote_step : couplings_) {
     621             :     for (auto& remote_substep : remote_step) {
     622             :       auto& local_step = remote_substep[n];
     623             :       local_step.erase(local_step.begin() + 1, local_step.end());
     624             :     }
     625             :   }
     626             : }
     627             : 
     628             : template <typename LocalData, typename RemoteData, typename CouplingResult>
     629             : void BoundaryHistory<LocalData, RemoteData,
     630             :                      CouplingResult>::clear_substeps_remote(const size_t n) {
     631             :   remote_data_[n].substeps.erase(remote_data_[n].substeps.begin() + 1,
     632             :                                  remote_data_[n].substeps.end());
     633             :   auto& remote_step = couplings_[n];
     634             :   remote_step.erase(remote_step.begin() + 1, remote_step.end());
     635             : }
     636             : 
     637             : template <typename LocalData, typename RemoteData, typename CouplingResult>
     638           0 : std::ostream& operator<<(
     639             :     std::ostream& os,
     640             :     const BoundaryHistory<LocalData, RemoteData, CouplingResult>& history) {
     641             :   return history.template print<true>(os);
     642             : }
     643             : }  // namespace TimeSteppers