Pass fluxes to elliptic source computers #2553

nilsvu · 2020-10-13T13:42:50Z

Proposed changes

More involved systems such as XCTS need all variables to compute source terms.

Upgrade instructions

The signature of functions in Elliptic/FirstOrderOperator.hpp have changed.

Types of changes:

Bugfix
New feature
Refactor

Component:

Code
Documentation
Build system
Continuous integration

Code review checklist

The PR passes all checks, including unit tests and clang-tidy.
For instructions on how to perform the CI checks locally refer to the Dev
guide on the Travis CI.
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.

Further comments

nilsdeppe · 2020-10-17T05:12:57Z

I don't understand the sentence "pass fluxes to sources". Things are either fluxes or sources, not both. Are the attempts at analogies between elliptic and conservative systems just bad? In the hyperbolic case you have

d_t u + d_i F^i(u) + B^i(u) d_i u= S(u)

where S(u) is not a function of F^i(u). You may happen to have some of the same terms showing up, but the sources definitely do not depend on the fluxes. Certain classes of hyperbolic-parabolic systems can be written in the form

d_t u + d_i F^i_A(u) + d_i F^i_D(u, d_j u) + B^i(u) d_i u= S(u)

whereA is for advective terms and D is for diffusion terms. Is that equation a better analogy to the elliptic systems?

nilsvu · 2020-10-17T12:20:05Z

@nilsdeppe Thanks for looking at this. I think the first-order conservative form still works for elliptic systems. Your variables u include both the "primal" fields and their derivatives, the "auxiliary" fields, that are evolved/solved along the primal fields. For example, a Poisson equation would use the gradient of the scalar field as auxiliary variable. For the XCTS momentum constraint I use the symmetric "strain" of the shift vector field as the auxiliary variable, i.e. the symmetrized gradient of the shift. Then the flux for the momentum constraint is the "longitudinal shift", which is some projection of the "strain" variable.
Now I agree the sources depend only on the variables u (primal and auxiliary). But as you say, terms show up that are identical to the fluxes. In particular, the sources of the momentum constraint include contractions with the longitudinal shift. I could always re-compute these terms from the variables, but it would be good to avoid that. So instead of passing primal and auxiliary variables to the sources (read "shift" and "shift strain") I pass the primal variables and their fluxes to the sources (read "shift and "longitudinal shift"). Do you think this is a reasonable thing to do? Any suggestions?

Note that I can't use argument_tags to retrieve pre-computed fluxes from the DataBox. Everything that depends on the variables must be passed into the functions. That's because the operator might be acting on some preconditioned or internal linear solver variables.

nilsdeppe · 2020-10-17T18:34:08Z

Okay, thanks for explaining that!

I was going to ask about argument_tags, but you already explained that can't work. My next thought, and this could very reasonably be a long-term change, would be if you can compute the sources and fluxes in one function. This would be similar to the new TimeDerivativeTerms structs in the evolution systems. If you always need the fluxes and sources at the same (iterative) time on the same grid, then you wouldn't need to do the song and dance that's done in the evolution where there are fluxes_impl and sources_impl functions that are wrapped by the TimeDerivativeTerms. The reason that's necessary in the evolution is because when we do FD and FV the sources and fluxes are computed on different grids (cell-center for sources, cell-faces for fluxes). This also in general means fewer allocations and compute tags, at least in the evolution code. Again, not anything that I think needs to be done immediately, but it would be good to consider whether something like that is feasible, what limitations there would be, and how it would interact with the current code.

nilsvu · 2020-10-19T13:05:59Z

@nilsdeppe Yes, I'm definitely planning to remove the compute tag for the sources in favor of a function that directly adds the sources to the -div(fluxes). At that point I might also change the interface for the fluxes and sources implementations to compute everything in one function. I'll see how this works out in the evolution code and then adapt it for the elliptic code :)

@kidder could you do the second review?

nilsdeppe

Sounds good! Just to be clear: only the volume fluxes are computed alongside the source terms so that user functions don't need to know whether we are using Gauss or GL points.

nilsdeppe · 2020-10-19T20:49:00Z

src/Elliptic/DiscontinuousGalerkin/InitializeFirstOrderOperator.hpp

@@ -71,7 +71,7 @@ struct InitializeFirstOrderOperator {
                                              vars_tag, PrimalVariables,


The commit message is incorrect. The systems don't need the fluxes to compute sources, but that having access to the fluxes is an optimization.

nilsdeppe · 2020-10-19T20:50:02Z

src/Elliptic/FirstOrderOperator.hpp

          using auxiliary_field_tag =
              tmpl::type_from<decltype(auxiliary_field_tag_v)>;
          get<::Tags::Source<auxiliary_field_tag>>(*sources) =
              get<auxiliary_field_tag>(vars);
        });
+    // Call into the sources computer to set primal field sources and possibly
+    // adjust auxiliary field sources


Can you add a comment here as to why the fluxes are being passed? As per our discussion, it is not at all clear that this code is correct without such a clarification.

nilsdeppe · 2020-10-19T20:50:46Z

src/Elliptic/FirstOrderOperator.hpp

@@ -120,28 +127,34 @@ auto first_order_fluxes(const Variables<VarsTags>& vars,
 *


add a discussion as to why the sources terms are getting the fluxes

nilsdeppe · 2020-10-19T20:51:25Z

src/Elliptic/Systems/Elasticity/Equations.hpp

@@ -94,7 +94,9 @@ struct Sources {
  using argument_tags = tmpl::list<>;
  static void apply(
      const gsl::not_null<tnsr::I<DataVector, Dim>*> source_for_displacement,
-      const tnsr::I<DataVector, Dim>& /*displacement*/) noexcept {
+      const gsl::not_null<tnsr::ii<DataVector, Dim>*> /*source_for_strain*/,
+      const tnsr::I<DataVector, Dim>& /*displacement*/,


aren't the new arguments fluxes? Maybe make that clear in the name?

The new source_for_strain argument is the source for the auxiliary equation (think the equation grad(u) = v that defines the auxiliary variable for a Poisson equation. The auxiliary variable for the elasticity system is the symmetric "strain", i.e. the symmetrized gradient of the displacement vector). The source computers gain this argument to get the opportunity to add non-principal contributions to the gradient (Christoffel-symbol terms) to the auxiliary equation.

The new stress argument is indeed the flux related to the displacement vector field. It's the strain contracted with the constitutive relation of the elastic material, and that's called "stress" in elasticity terms. I think in the context of the elasticity equations it makes sense to call it by its name.

Okay, I can't speak to what would be most intuitive for people used to the elasticity literature, so I'll have to take your word for it :)

nilsvu · 2020-10-21T10:02:10Z

@nilsdeppe Fixed the commit message and added docs. Please take a look :)

@kidder Could you do the second review?

nilsdeppe · 2020-10-21T16:02:24Z

@kidder already has a lot of PRs awaiting his review so please try to find someone else.

nilsvu · 2020-10-21T16:54:02Z

Sure! @wthrowe could you do me the honour? :)

wthrowe · 2020-10-21T19:58:11Z

Certainly. This looks good to me.

nilsvu · 2020-10-22T08:35:44Z

Thanks @wthrowe!

@nilsdeppe Did you look at the fixup already, or should I keep it around a bit longer?

nilsdeppe · 2020-10-22T15:08:11Z

src/Elliptic/FirstOrderOperator.hpp

+    // adjust auxiliary field sources
+    SourcesComputer::apply(
+        make_not_null(&get<::Tags::Source<PrimalFields>>(*sources))...,
+        make_not_null(&get<::Tags::Source<AuxiliaryFields>>(*sources))...,


You keep saying "instead of the auxiliary fields" but here you pass the aux fields. Please correct this

No this passes the sources for the auxiliary fields so terms can be added to them, such as Christoffel symbol terms (see this comment #2553 (comment)). Does that make sense?

Okay, that makes sense.

I still don't understand why you say "instead of" because the code didn't previously pass the aux vars. Why wouldn't you need the aux vars (I realize they weren't passed before)? Is that because technically they are just (roughly) the derivatives of the primal fields? If so, I'm a bit confused that stability and well-posedness aren't adversely affected by never using the aux variables, though explaining that is beyond the scope of a GitHub comment, I think ;)

Previously I didn't pass
(i) the auxiliary sources, because they are "usually" just the auxiliary vars, i.e. when no Christoffel-terms need to be added (e.g. consider the Poisson auxiliary equation -grad(u) + v = 0 where the sources are just v. Remember the auxiliary equation essentially defines the aux vars, which is why the aux vars are always its sources). The aux sources are set to the aux vars just a few lines above.
(ii) the auxiliary variables, because none of the systems that are currently implemented needs to compute non-trivial sources: Both the Poisson and the (linear) elasticity equations are only sourced by a fixed function f(x).

In the comments I say "instead of" because the sources are a function of the variables (primal and auxiliary), but instead of the aux vars I pass the primal fluxes as an optimization, as we discussed before.

Does this clear things up?

Yes, but I think this lengthy discussion suggests that you should either remove the instead of or expand on it in both the docs and the git message

Adjusted the commit message and the code comment in the fixup commit, please take a look.

nilsdeppe · 2020-10-22T15:09:27Z

src/Elliptic/Systems/Elasticity/Equations.hpp

@@ -94,7 +94,9 @@ struct Sources {
  using argument_tags = tmpl::list<>;
  static void apply(
      const gsl::not_null<tnsr::I<DataVector, Dim>*> source_for_displacement,
-      const tnsr::I<DataVector, Dim>& /*displacement*/) noexcept {
+      const gsl::not_null<tnsr::ii<DataVector, Dim>*> /*source_for_strain*/,
+      const tnsr::I<DataVector, Dim>& /*displacement*/,


Okay, I can't speak to what would be most intuitive for people used to the elasticity literature, so I'll have to take your word for it :)

nilsvu · 2020-10-26T08:46:26Z

@nilsdeppe could you check the fixup? I adjusted the commit message and the code comment

nilsdeppe

Looks good. Please go ahead and squash!

More involved systems such as XCTS need all variables to compute source terms. We pass the computed volume fluxes to the source computers as an optimization so they don't have to be re-computed.

nilsvu · 2020-10-29T22:27:19Z

Rebased and squashed. Thanks for your reviews @nilsdeppe and @wthrowe.

nilsvu · 2020-11-03T10:10:40Z

@wthrowe I think this is ready to go when you get the chance

nilsvu force-pushed the pass_all_vars_to_elliptic_computers branch 4 times, most recently from c87c621 to 2439033 Compare October 14, 2020 07:40

nilsvu added the in progress Don't review, used for sharing code and getting feedback label Oct 14, 2020

nilsvu force-pushed the pass_all_vars_to_elliptic_computers branch from 2439033 to 6cea811 Compare October 15, 2020 08:18

nilsvu changed the title ~~Pass all variables to elliptic fluxes and sources~~ Pass fluxes to elliptic source computers Oct 15, 2020

nilsvu force-pushed the pass_all_vars_to_elliptic_computers branch from 6cea811 to f570b71 Compare October 15, 2020 14:01

nilsvu removed the in progress Don't review, used for sharing code and getting feedback label Oct 16, 2020

nilsvu requested a review from kidder October 19, 2020 13:06

nilsdeppe reviewed Oct 19, 2020

View reviewed changes

nilsvu force-pushed the pass_all_vars_to_elliptic_computers branch from f570b71 to 4c46932 Compare October 21, 2020 10:01

nilsvu removed the request for review from kidder October 22, 2020 08:34

nilsdeppe reviewed Oct 22, 2020

View reviewed changes

nilsvu force-pushed the pass_all_vars_to_elliptic_computers branch from 4c46932 to 7fe91fc Compare October 22, 2020 19:43

nilsdeppe reviewed Oct 27, 2020

View reviewed changes

nilsvu added 2 commits October 29, 2020 16:18

Delete unnecessary Frame::Inertial params from Poisson system

4ba34a8

Pass fluxes to elliptic source computers

8ef75ee

More involved systems such as XCTS need all variables to compute source terms. We pass the computed volume fluxes to the source computers as an optimization so they don't have to be re-computed.

nilsvu force-pushed the pass_all_vars_to_elliptic_computers branch from 7fe91fc to 8ef75ee Compare October 29, 2020 15:19

nilsdeppe approved these changes Oct 31, 2020

View reviewed changes

wthrowe approved these changes Nov 3, 2020

View reviewed changes

wthrowe merged commit 1b4247a into sxs-collaboration:develop Nov 3, 2020

nilsvu deleted the pass_all_vars_to_elliptic_computers branch December 22, 2020 23:56

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Pass fluxes to elliptic source computers #2553

Pass fluxes to elliptic source computers #2553

nilsvu commented Oct 13, 2020 •

edited

nilsdeppe commented Oct 17, 2020

nilsvu commented Oct 17, 2020

nilsdeppe commented Oct 17, 2020

nilsvu commented Oct 19, 2020

nilsdeppe left a comment

nilsdeppe Oct 19, 2020

nilsdeppe Oct 19, 2020

nilsdeppe Oct 19, 2020

nilsdeppe Oct 19, 2020

nilsvu Oct 21, 2020

nilsdeppe Oct 22, 2020

nilsvu commented Oct 21, 2020

nilsdeppe commented Oct 21, 2020

nilsvu commented Oct 21, 2020

wthrowe commented Oct 21, 2020

nilsvu commented Oct 22, 2020

nilsdeppe Oct 22, 2020

nilsvu Oct 22, 2020

nilsdeppe Oct 22, 2020

nilsvu Oct 22, 2020

nilsdeppe Oct 22, 2020

nilsvu Oct 22, 2020

nilsdeppe Oct 22, 2020

nilsvu commented Oct 26, 2020

nilsdeppe left a comment

nilsvu commented Oct 29, 2020

nilsvu commented Nov 3, 2020

		@@ -71,7 +71,7 @@ struct InitializeFirstOrderOperator {
		vars_tag, PrimalVariables,

		@@ -120,28 +127,34 @@ auto first_order_fluxes(const Variables<VarsTags>& vars,
		*

Pass fluxes to elliptic source computers #2553

Pass fluxes to elliptic source computers #2553

Conversation

nilsvu commented Oct 13, 2020 • edited

Proposed changes

Upgrade instructions

Types of changes:

Component:

Code review checklist

Further comments

nilsdeppe commented Oct 17, 2020

nilsvu commented Oct 17, 2020

nilsdeppe commented Oct 17, 2020

nilsvu commented Oct 19, 2020

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

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

nilsvu commented Oct 21, 2020

nilsdeppe commented Oct 21, 2020

nilsvu commented Oct 21, 2020

wthrowe commented Oct 21, 2020

nilsvu commented Oct 22, 2020

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

nilsvu commented Oct 26, 2020

nilsdeppe left a comment

Choose a reason for hiding this comment

nilsvu commented Oct 29, 2020

nilsvu commented Nov 3, 2020

nilsvu commented Oct 13, 2020 •

edited