Add IMEX protocols #5618

wthrowe · 2023-11-07T18:42:23Z

The examples for the ImplicitSector protocol that will be added later with the actual solver are

template <bool TestWithAnalyticSolution>
struct ImplicitSector : tt::ConformsTo<imex::protocols::ImplicitSector> {
  using tensors = tmpl::list<Var2, Var3>;
  using initial_guess = tmpl::conditional_t<TestWithAnalyticSolution,
                                            AnalyticSolution, InitialGuess>;
 
  struct SolveAttempt {
    using tags_from_evolution =
        tmpl::list<Var1, NonTensor, VariablesFromEvolution>;
    using simple_tags = tmpl::list<RecordPreparersForTest, VariablesTemporary>;
    using compute_tags = tmpl::list<SomeComputeTag>;
 
    using source_prep =
        tmpl::list<Preparer<PrepId::Shared>, Preparer<PrepId::Source>>;
    using jacobian_prep =
        tmpl::list<Preparer<PrepId::Shared>, Preparer<PrepId::Jacobian>>;
 
    using source = Source;
    using jacobian =
        tmpl::conditional_t<TestWithAnalyticSolution,
                            imex::NoJacobianBecauseSolutionIsAnalytic,
                            Jacobian>;
  };
 
  using solve_attempts = tmpl::list<SolveAttempt>;
};

struct AnalyticSolution {
  using return_tags = tmpl::list<Var2, Var3>;
  using argument_tags = tmpl::list<Var1, NonTensor>;
  static std::vector<imex::GuessResult> apply(
      const gsl::not_null<tnsr::II<DataVector, 2>*> var2,
      const gsl::not_null<tnsr::I<DataVector, 2>*> var3,
      const Scalar<DataVector>& var1, const double non_tensor,
      const Variables<tmpl::list<Var2, Var3>>& inhomogeneous_terms,
      const double implicit_weight) {
    // Solution for source terms
    // S[v2^ij] = v3^i v3^j - nt v2^ij
    // S[v3^i] = -v1 v3^i
 
    // Solving  v3^i = X - w v1 v3^i  gives  v3^i = X / (1 + w v1)
    tenex::evaluate<ti::I>(var3, get<Var3>(inhomogeneous_terms)(ti::I) /
                                     (1.0 + implicit_weight * var1()));
    tenex::evaluate<ti::I, ti::J>(
        var2, (get<Var2>(inhomogeneous_terms)(ti::I, ti::J) +
               implicit_weight * (*var3)(ti::I) * (*var3)(ti::J)) /
                  (1.0 + implicit_weight * non_tensor));
    return {get(var1).size(), imex::GuessResult::ExactSolution};
  }
};

Edit: decided to add two more examples:

struct Source {
  using return_tags = tmpl::list<::Tags::Source<Var2>, ::Tags::Source<Var3>>;
  using argument_tags = tmpl::list<Var1, Var2, Var3, NonTensor,
                                   RecordPreparersForTest, SomeComputeTagBase>;
 
  static void apply(const gsl::not_null<tnsr::II<DataVector, 2>*> source_var2,
                    const gsl::not_null<tnsr::I<DataVector, 2>*> source_var3,
                    const Scalar<DataVector>& var1,
                    const tnsr::II<DataVector, 2>& var2,
                    const tnsr::I<DataVector, 2>& var3, const double non_tensor,
                    const RecordPreparersForTest::type& prep_run_values,
                    const double compute_tag_value) {

struct Jacobian {
  using return_tags =
      tmpl::list<imex::Tags::Jacobian<Var2, ::Tags::Source<Var2>>,
                 imex::Tags::Jacobian<Var3, ::Tags::Source<Var2>>,
                 imex::Tags::Jacobian<Var3, ::Tags::Source<Var3>>>;
  using argument_tags =
      tmpl::list<Var1, Var3, NonTensor, RecordPreparersForTest>;
 
  static void apply(const gsl::not_null<tnsr::iiJJ<DataVector, 2>*> dvar2_dvar2,
                    const gsl::not_null<tnsr::iJJ<DataVector, 2>*> dvar2_dvar3,
                    const gsl::not_null<tnsr::iJ<DataVector, 2>*> dvar3_dvar3,
                    const Scalar<DataVector>& var1,
                    const tnsr::I<DataVector, 2>& var3, const double non_tensor,
                    const RecordPreparersForTest::type& prep_run_values) {

Proposed changes

Upgrade instructions

Code review checklist

The code is documented and the documentation renders correctly. Run
make doc to generate the documentation locally into BUILD_DIR/docs/html.
Then open index.html.
The code follows the stylistic and code quality guidelines listed in the
code review guide.
The PR lists upgrade instructions and is labeled bugfix or
new feature if appropriate.

Further comments

nilsdeppe

You can squash any changes directly

nilsdeppe · 2023-11-11T17:30:51Z

src/Evolution/Imex/CMakeLists.txt

@@ -0,0 +1,28 @@
+# Distributed under the MIT License.


Shouldn't there be a add_subdirectory one level higher in this commit?

nilsdeppe · 2023-11-11T17:30:57Z

tests/Unit/Evolution/Imex/CMakeLists.txt

@@ -0,0 +1,26 @@
+# Distributed under the MIT License.


Shouldn't there be a add_subdirectory one level higher in this commit?

nilsdeppe · 2023-11-11T17:32:17Z

src/Evolution/Imex/Protocols/ImexSystem.hpp

+///
+/// In addition to the usual requirements for an evolution system, an
+/// IMEX system must specify an `implicit_sectors` typelist of structs
+/// conforming to protocols::ImplicitSector, each of which described


described -> describes?

yoonso0-0

Some minor fix requests on documentation as the first customer :)

Please squash everything

yoonso0-0 · 2023-11-15T19:37:41Z

src/Evolution/Imex/GuessResult.hpp

+#include "Utilities/TMPL.hpp"
+
+namespace imex {
+/// Accuracy of a guess returned from an implicit sector's


[nitpick] Accuracy -> Type?

yoonso0-0 · 2023-11-15T19:42:56Z

src/Evolution/Imex/Protocols/ImplicitSector.hpp

+///   `std::vector<GuessResult>` indicating, for each point, whether
+///   the implicit equation has been solved analytically or whether
+///   the numerical solve should continue.  An empty return is
+///   equivalent to `GuessResult::InitialGuess` for every point.  When


imex::GuessResult::InitialGuess

Q. Does this mean that an empty return is exactly same as using imex::GuessExplicitResult ?

No. GuessExplicitResult guesses the explicit result. You are free to guess something else. I added a bit to the documentation of the GuessExplicitResult struct. Do you think something should be added here?

yoonso0-0 · 2023-11-15T19:46:17Z

src/Evolution/Imex/Protocols/ImplicitSector.hpp

+///   equivalent to `GuessResult::InitialGuess` for every point.  When
+///   this is called, the sector variables will have the value from
+///   the explicit part of the time step.  This mutator will not be
+///   called if the implicit weight is zero, as a system-independent


I think the comment after ".. weight is zero" is a little confusing, since implicit weight being zero means no implicit step is required, not that it's analytically invertible.

Sort of. The step becomes explicit and no source evaluations are performed, but there is still a step. Admittedly, calling the result of solving y = constant for y an analytic solution is a bit grandiose, so I changed it to "This mutator will not be called if the implicit weight is zero, as the solution is trivial in that case.".

yoonso0-0 · 2023-11-15T19:49:02Z

src/Evolution/Imex/Protocols/ImplicitSector.hpp

+///     * `tags_from_evolution` for tags to be made available from the
+///        evolution DataBox.  Volume quantities will be reduced to
+///        have one grid point, with the appropriate value for the
+///        point being solved for.


Please mention that this excludes the tags from tensors (that needs to be solved).

Done, and also added a check to the protocol verification below.

wthrowe

Also added a sanity check now assumed after some local changes to the solver.

wthrowe · 2023-11-28T17:58:22Z

src/Evolution/Imex/Protocols/ImplicitSector.hpp

+///   `std::vector<GuessResult>` indicating, for each point, whether
+///   the implicit equation has been solved analytically or whether
+///   the numerical solve should continue.  An empty return is
+///   equivalent to `GuessResult::InitialGuess` for every point.  When


No. GuessExplicitResult guesses the explicit result. You are free to guess something else. I added a bit to the documentation of the GuessExplicitResult struct. Do you think something should be added here?

wthrowe · 2023-11-28T18:17:58Z

src/Evolution/Imex/Protocols/ImplicitSector.hpp

+///   equivalent to `GuessResult::InitialGuess` for every point.  When
+///   this is called, the sector variables will have the value from
+///   the explicit part of the time step.  This mutator will not be
+///   called if the implicit weight is zero, as a system-independent


Sort of. The step becomes explicit and no source evaluations are performed, but there is still a step. Admittedly, calling the result of solving y = constant for y an analytic solution is a bit grandiose, so I changed it to "This mutator will not be called if the implicit weight is zero, as the solution is trivial in that case.".

wthrowe · 2023-11-28T18:18:55Z

src/Evolution/Imex/Protocols/ImplicitSector.hpp

+///     * `tags_from_evolution` for tags to be made available from the
+///        evolution DataBox.  Volume quantities will be reduced to
+///        have one grid point, with the appropriate value for the
+///        point being solved for.


Done, and also added a check to the protocol verification below.

yoonso0-0

Looks good. Please squash

yoonso0-0 · 2023-11-29T18:51:15Z

src/Evolution/Imex/Protocols/ImplicitSector.hpp

+///   `std::vector<GuessResult>` indicating, for each point, whether
+///   the implicit equation has been solved analytically or whether
+///   the numerical solve should continue.  An empty return is
+///   equivalent to `imex::GuessResult::InitialGuess` for every point.


Maybe a short extension like "empty return is equivalent to returning GuessResult::InitialGuess everywhere and the numerical solve will be done for every point"?

nilsdeppe · 2023-11-29T19:13:24Z

I'm also happy! Please go ahead and squash

The examples will be added in a future commit

wthrowe added the new design technical with long discussions expected label Nov 7, 2023

yoonso0-0 self-requested a review November 8, 2023 06:57

nilsdeppe requested changes Nov 11, 2023

View reviewed changes

wthrowe force-pushed the imex/protocols branch from 6fe4445 to b156154 Compare November 15, 2023 19:07

yoonso0-0 reviewed Nov 15, 2023

View reviewed changes

wthrowe force-pushed the imex/protocols branch from b156154 to ba3abdc Compare November 28, 2023 19:36

wthrowe commented Nov 28, 2023

View reviewed changes

yoonso0-0 reviewed Nov 29, 2023

View reviewed changes

wthrowe added 2 commits November 30, 2023 16:36

Add IMEX GuessResult and convenience structs

853d24e

Add IMEX protocols

329ee0d

The examples will be added in a future commit

wthrowe force-pushed the imex/protocols branch from ba3abdc to 329ee0d Compare December 1, 2023 01:53

nilsdeppe approved these changes Dec 1, 2023

View reviewed changes

nilsdeppe merged commit 6ab1b8c into sxs-collaboration:develop Dec 1, 2023
21 of 22 checks passed

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Add IMEX protocols #5618

Add IMEX protocols #5618

wthrowe commented Nov 7, 2023 •

edited

Loading

nilsdeppe left a comment

nilsdeppe Nov 11, 2023

nilsdeppe Nov 11, 2023

nilsdeppe Nov 11, 2023

yoonso0-0 left a comment

yoonso0-0 Nov 15, 2023

yoonso0-0 Nov 15, 2023

yoonso0-0 Nov 15, 2023

wthrowe Nov 28, 2023

yoonso0-0 Nov 15, 2023

wthrowe Nov 28, 2023

yoonso0-0 Nov 15, 2023

wthrowe Nov 28, 2023

wthrowe left a comment

wthrowe Nov 28, 2023

wthrowe Nov 28, 2023

wthrowe Nov 28, 2023

yoonso0-0 left a comment

yoonso0-0 Nov 29, 2023

nilsdeppe commented Nov 29, 2023

Add IMEX protocols #5618

Add IMEX protocols #5618

Conversation

wthrowe commented Nov 7, 2023 • edited Loading

Proposed changes

Upgrade instructions

Code review checklist

Further comments

nilsdeppe left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

yoonso0-0 left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

wthrowe left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

yoonso0-0 left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

nilsdeppe commented Nov 29, 2023

wthrowe commented Nov 7, 2023 •

edited

Loading