FirstOrderInternalPenalty.hpp
1 // Distributed under the MIT License.
2 // See LICENSE.txt for details.
3 
4 #pragma once
5 
6 #include <cstddef>
7 #include <pup.h>
8 #include <string>
9 
12 #include "DataStructures/DataBox/Tag.hpp"
13 #include "DataStructures/Tensor/EagerMath/Magnitude.hpp"
17 #include "Domain/FaceNormal.hpp"
18 #include "Domain/Tags.hpp"
19 #include "NumericalAlgorithms/DiscontinuousGalerkin/NormalDotFlux.hpp"
20 #include "NumericalAlgorithms/DiscontinuousGalerkin/Protocols.hpp"
22 #include "Options/Options.hpp"
23 #include "Utilities/Gsl.hpp"
24 #include "Utilities/ProtocolHelpers.hpp"
25 #include "Utilities/TMPL.hpp"
26 
27 /// \cond
28 class DataVector;
29 namespace Tags {
30 template <typename>
31 struct Normalized;
32 } // namespace Tags
33 /// \endcond
34 
35 /// Numerical fluxes for elliptic systems
37 
38 /*!
39  * \brief The penalty factor in internal penalty fluxes
40  *
41  * The penalty factor is computed as
42  *
43  * \f{equation}
44  * \sigma = C \frac{N_\text{points}^2}{h}
45  * \f}
46  *
47  * where \f$N_\text{points} = 1 + N_p\f$ is the number of points (or one plus
48  * the polynomial degree) and \f$h\f$ is a measure of the element size. Both
49  * quantities are taken perpendicular to the face of the DG element that the
50  * penalty is being computed on. \f$C\f$ is the *penalty parameter*.
51  *
52  * \see `elliptic::dg::NumericalFluxes::FirstOrderInternalPenalty` for details
53  */
54 DataVector penalty(const DataVector& element_size, const size_t num_points,
55  const double penalty_parameter) noexcept {
56  return penalty_parameter * square(num_points) / element_size;
57 }
58 
59 /*!
60  * \ingroup DiscontinuousGalerkinGroup
61  * \ingroup NumericalFluxesGroup
62  * \brief The internal penalty flux for first-order elliptic equations.
63  *
64  * \details Computes the internal penalty numerical flux (see e.g.
65  * \cite HesthavenWarburton, section 7.2) dotted with the interface unit normal.
66  *
67  * We implement here a suggested generalization of the internal penalty flux
68  * for any set of elliptic PDEs up to second order in derivatives. It is
69  * designed for fluxes (i.e. principal parts of the PDEs) that may depend on the
70  * dynamic fields, but do so at most linearly. This is the case for the velocity
71  * potential equation of binary neutron stars, for example. The linearity is
72  * only necessary to impose inhomogeneous boundary conditions as contributions
73  * to the fixed source (see
74  * `elliptic::dg::Actions::ImposeInhomogeneousBoundaryConditionsOnSource`).
75  *
76  * We reduce the second-order elliptic PDEs to first-order form by introducing
77  * an _auxiliary_ variable \f$v\f$ for each _primal_ variable \f$u\f$ (i.e. for
78  * each variable whose second derivative appears in the equations). Then, the
79  * equations take the _flux-form_
80  *
81  * \f{align}
82  * -\partial_i F^i_A + S_A = f_A
83  * \f}
84  *
85  * where the fluxes \f$F^i_A(u,v)\f$ and sources \f$S_A(u,v)\f$ are indexed
86  * (with capital letters) by the variables and we have defined \f$f_A(x)\f$ as
87  * those sources that are independent of the variables. Note that the fluxes are
88  * functions of the primal and auxiliary variables; e.g. for a Poisson system
89  * \f$\{u,v_i\}\f$ on a spatial metric \f$\gamma_{ij}\f$ they are simply
90  * \f$F^i_u(v)=\sqrt{\gamma}\gamma^{ij} v_j\f$ and
91  * \f$F^i_{v_j}(u)=u\delta^i_j\f$ (see `Poisson::FirstOrderSystem`) or for an
92  * Elasticity system \f$\{\xi^i,S_{ij}\}\f$ with Young's tensor
93  * \f$Y^{ijkl}\f$ they are \f$F^i_{\xi^j}(S)=Y^{ijkl}S_{kl}\f$ and
94  * \f$F^i_{S_{jk}}(\xi)=\delta^i_{(j}\xi_{k)}\f$. Since each primal flux is
95  * commonly only a function of the corresponding auxiliary variable and
96  * vice-versa, we omit the other function arguments here to lighten the
97  * notational load.
98  *
99  * We now choose the internal penalty numerical fluxes \f$(n_i F^i_A)^*\f$ as
100  * follows for each primal variable \f$u\f$ and their corresponding auxiliary
101  * variable \f$v\f$:
102  *
103  * \f{align}
104  * (n_i F^i_u)^* &= \frac{1}{2} n_i \left( F^i_u(\partial_j
105  * F^j_v(u^\mathrm{int})) + F^i_u(\partial_j F^j_v(u^\mathrm{ext}))
106  * \right) - \sigma n_i \left(F^i_u(n_j F^j_v(u^\mathrm{int})) - F^i_u(
107  * n_j F^j_v(u^\mathrm{ext}))\right) \\
108  * (n_i F^i_v)^* &= \frac{1}{2} n_i \left(F^i_v(u^\mathrm{int}) +
109  * F^i_v(u^\mathrm{ext})\right)
110  * \f}
111  *
112  * Note that we have assumed \f$n^\mathrm{ext}_i=-n_i\f$ here, i.e. face normals
113  * don't depend on the dynamic variables (which may be discontinuous on element
114  * faces). This is the case for the problems we are expecting to solve, because
115  * those will be on fixed background metrics (e.g. a conformal metric for the
116  * XCTS system).
117  *
118  * Also note that the numerical fluxes intentionally don't depend on the
119  * auxiliary field values \f$v\f$. This property allows us to use the numerical
120  * fluxes also for the second-order (or _primal_) DG formulation, where we
121  * remove the need for an auxiliary variable. For the first-order system we
122  * could replace the divergence in \f$(n_i F^i_u)^*\f$ with \f$v\f$, which would
123  * result in a generalized "stabilized central flux" that is slightly less
124  * sparse than the internal penalty flux (see e.g. \cite HesthavenWarburton,
125  * section 7.2). We could also choose to ignore the fluxes in the penalty term,
126  * but preliminary tests suggest that this may hurt convergence.
127  *
128  * For a Poisson system (see above) this numeric flux reduces to the standard
129  * internal penalty flux (see e.g. \cite HesthavenWarburton, section 7.2, or
130  * \cite Arnold2002)
131  *
132  * \f{align}
133  * (n_i F^i_u)^* &= n_i v_i^* = \frac{1}{2} n_i \left(\partial_i u^\mathrm{int}
134  * + \partial_i u^\mathrm{ext}\right) - \sigma \left(u^\mathrm{int} -
135  * u^\mathrm{ext}\right) \\
136  * (n_i F^i_{v_j})^* &= n_j u^* = \frac{1}{2} n_j \left(u^\mathrm{int} +
137  * u^\mathrm{ext}\right)
138  * \f}
139  *
140  * where a sum over repeated indices is assumed, since the equation is
141  * formulated on a Euclidean geometry.
142  *
143  * This generalization of the internal penalty flux is based on unpublished work
144  * by Nils L. Fischer (nils.fischer@aei.mpg.de).
145  *
146  * The penalty factor \f$\sigma\f$ is responsible for removing zero eigenmodes
147  * and impacts the conditioning of the linear operator to be solved. It can be
148  * chosen as \f$\sigma=C\frac{N_\mathrm{points}^2}{h}\f$ where
149  * \f$N_\mathrm{points}\f$ is the number of collocation points (i.e. the
150  * polynomial degree plus 1), \f$h\f$ is a measure of the element size in
151  * inertial coordinates (both measured perpendicular to the interface under
152  * consideration) and \f$C\geq 1\f$ is a free parameter (see e.g.
153  * \cite HesthavenWarburton, section 7.2). For the element size we choose
154  * \f$h=\frac{J_\mathrm{volume}}{J_\mathrm{face}}\f$, i.e. the ratio of Jacobi
155  * determinants from logical to inertial coordinates in the element volume and
156  * on the element face, both evaluated on the face (see \cite Vincent2019qpd).
157  * Since both \f$N_\mathrm{points}\f$ and \f$h\f$ can be different on either
158  * side of the element boundary we take the maximum of \f$N_\mathrm{points}\f$
159  * and the pointwise minimum of \f$h\f$ across the element boundary as is done
160  * in \cite Vincent2019qpd. Note that we use the number of points
161  * \f$N_\mathrm{points}\f$ where \cite Vincent2019qpd uses the polynomial degree
162  * \f$N_\mathrm{points} - 1\f$ because we found unstable configurations on
163  * curved meshes when using the polynomial degree. Optimizing the penalty on
164  * curved meshes is subject to further investigation.
165  */
166 template <size_t Dim, typename FluxesComputerTag, typename FieldTagsList,
167  typename AuxiliaryFieldTagsList>
169 
170 template <size_t Dim, typename FluxesComputerTag, typename... FieldTags,
171  typename... AuxiliaryFieldTags>
172 struct FirstOrderInternalPenalty<Dim, FluxesComputerTag,
173  tmpl::list<FieldTags...>,
174  tmpl::list<AuxiliaryFieldTags...>>
175  : tt::ConformsTo<::dg::protocols::NumericalFlux> {
176  private:
177  using fluxes_computer_tag = FluxesComputerTag;
178  using FluxesComputer = typename fluxes_computer_tag::type;
179 
180  struct PerpendicularNumPoints : db::SimpleTag {
181  using type = size_t;
182  };
183 
184  struct ElementSize : db::SimpleTag {
185  using type = Scalar<DataVector>;
186  };
187 
188  template <typename Tag>
189  struct NormalDotDivFlux : db::PrefixTag, db::SimpleTag {
190  using tag = Tag;
192  };
193 
194  template <typename Tag>
195  struct NormalDotNormalDotFlux : db::PrefixTag, db::SimpleTag {
196  using tag = Tag;
198  };
199 
200  public:
201  struct PenaltyParameter {
202  using type = double;
203  static constexpr Options::String help = {
204  "The prefactor to the penalty term of the flux. Values closer to one "
205  "lead to better-conditioned problems, but on curved meshes the penalty "
206  "parameter may need to be increased to keep the problem well-defined."};
207  static type lower_bound() noexcept { return 1.; }
208  };
209  using options = tmpl::list<PenaltyParameter>;
210  static constexpr Options::String help = {
211  "The internal penalty flux for elliptic systems."};
212  static std::string name() noexcept { return "InternalPenalty"; }
213 
214  FirstOrderInternalPenalty() = default;
215  explicit FirstOrderInternalPenalty(const double penalty_parameter) noexcept
216  : penalty_parameter_(penalty_parameter) {}
217 
218  // clang-tidy: non-const reference
219  void pup(PUP::er& p) noexcept { p | penalty_parameter_; } // NOLINT
220 
221  using variables_tags = tmpl::list<FieldTags..., AuxiliaryFieldTags...>;
222 
223  using argument_tags = tmpl::push_front<
224  typename FluxesComputer::argument_tags, domain::Tags::Mesh<Dim>,
229  Frame::Inertial>>...,
231  fluxes_computer_tag>;
232  using volume_tags =
233  tmpl::push_front<get_volume_tags<FluxesComputer>, domain::Tags::Mesh<Dim>,
234  fluxes_computer_tag>;
235 
236  using package_field_tags =
237  tmpl::list<::Tags::NormalDotFlux<AuxiliaryFieldTags>...,
238  NormalDotDivFlux<AuxiliaryFieldTags>...,
239  NormalDotNormalDotFlux<AuxiliaryFieldTags>..., ElementSize>;
240  using package_extra_tags = tmpl::list<PerpendicularNumPoints>;
241 
242  template <typename... FluxesArgs>
243  void package_data(
244  const gsl::not_null<typename ::Tags::NormalDotFlux<
245  AuxiliaryFieldTags>::type*>... packaged_n_dot_aux_fluxes,
246  const gsl::not_null<typename NormalDotDivFlux<
247  AuxiliaryFieldTags>::type*>... packaged_n_dot_div_aux_fluxes,
248  const gsl::not_null<
249  typename NormalDotNormalDotFlux<AuxiliaryFieldTags>::
250  type*>... packaged_n_dot_principal_n_dot_aux_fluxes,
251  const gsl::not_null<Scalar<DataVector>*> element_size,
252  const gsl::not_null<size_t*> perpendicular_num_points,
253  const Mesh<Dim>& volume_mesh, const Direction<Dim>& direction,
254  const Scalar<DataVector>& face_normal_magnitude,
255  const typename ::Tags::NormalDotFlux<
256  AuxiliaryFieldTags>::type&... normal_dot_auxiliary_field_fluxes,
257  const typename ::Tags::div<
258  ::Tags::Flux<AuxiliaryFieldTags, tmpl::size_t<Dim>,
259  Frame::Inertial>>::type&... div_auxiliary_field_fluxes,
260  const tnsr::i<DataVector, Dim, Frame::Inertial>& interface_unit_normal,
261  const FluxesComputer& fluxes_computer,
262  const FluxesArgs&... fluxes_args) const noexcept {
263  *perpendicular_num_points = volume_mesh.extents(direction.dimension());
264  get(*element_size) = 2. / get(face_normal_magnitude);
265  auto principal_div_aux_field_fluxes = make_with_value<Variables<tmpl::list<
267  interface_unit_normal, std::numeric_limits<double>::signaling_NaN());
268  auto principal_ndot_aux_field_fluxes = make_with_value<Variables<tmpl::list<
270  interface_unit_normal, std::numeric_limits<double>::signaling_NaN());
271  fluxes_computer.apply(
273  &get<::Tags::Flux<FieldTags, tmpl::size_t<Dim>, Frame::Inertial>>(
274  principal_div_aux_field_fluxes))...,
275  fluxes_args..., div_auxiliary_field_fluxes...);
276  fluxes_computer.apply(
278  &get<::Tags::Flux<FieldTags, tmpl::size_t<Dim>, Frame::Inertial>>(
279  principal_ndot_aux_field_fluxes))...,
280  fluxes_args..., normal_dot_auxiliary_field_fluxes...);
281  const auto apply_normal_dot =
282  [&principal_div_aux_field_fluxes, &principal_ndot_aux_field_fluxes,
283  &interface_unit_normal](
284  auto field_tag_v, auto packaged_normal_dot_auxiliary_field_flux,
285  auto packaged_normal_dot_div_auxiliary_field_flux,
286  auto packaged_normal_dot_principal_ndot_aux_field_flux,
287  const auto& normal_dot_auxiliary_field_flux) noexcept {
288  using field_tag = decltype(field_tag_v);
289  // Compute n.F_v(u)
290  *packaged_normal_dot_auxiliary_field_flux =
291  normal_dot_auxiliary_field_flux;
292  // Compute n.F_u(div(F_v(u))) and n.F_u(n.F_v(u))
293  normal_dot_flux(
294  packaged_normal_dot_div_auxiliary_field_flux,
295  interface_unit_normal,
296  get<::Tags::Flux<field_tag, tmpl::size_t<Dim>, Frame::Inertial>>(
297  principal_div_aux_field_fluxes));
298  normal_dot_flux(
299  packaged_normal_dot_principal_ndot_aux_field_flux,
300  interface_unit_normal,
301  get<::Tags::Flux<field_tag, tmpl::size_t<Dim>, Frame::Inertial>>(
302  principal_ndot_aux_field_fluxes));
303  };
304  EXPAND_PACK_LEFT_TO_RIGHT(apply_normal_dot(
305  FieldTags{}, packaged_n_dot_aux_fluxes, packaged_n_dot_div_aux_fluxes,
306  packaged_n_dot_principal_n_dot_aux_fluxes,
307  normal_dot_auxiliary_field_fluxes));
308  }
309 
310  void operator()(
311  const gsl::not_null<typename ::Tags::NormalDotNumericalFlux<
312  FieldTags>::type*>... numerical_flux_for_fields,
313  const gsl::not_null<typename ::Tags::NormalDotNumericalFlux<
314  AuxiliaryFieldTags>::type*>... numerical_flux_for_auxiliary_fields,
315  const typename ::Tags::NormalDotFlux<
316  AuxiliaryFieldTags>::type&... normal_dot_auxiliary_flux_interiors,
317  const typename NormalDotDivFlux<
318  AuxiliaryFieldTags>::type&... normal_dot_div_auxiliary_flux_interiors,
319  const typename NormalDotNormalDotFlux<
320  AuxiliaryFieldTags>::type&... ndot_ndot_aux_flux_interiors,
321  const Scalar<DataVector>& element_size_interior,
322  const size_t perpendicular_num_points_interior,
323  const typename ::Tags::NormalDotFlux<AuxiliaryFieldTags>::
324  type&... minus_normal_dot_auxiliary_flux_exteriors,
325  const typename NormalDotDivFlux<
326  AuxiliaryFieldTags>::type&... minus_normal_dot_div_aux_flux_exteriors,
327  const typename NormalDotDivFlux<
328  AuxiliaryFieldTags>::type&... ndot_ndot_aux_flux_exteriors,
329  const Scalar<DataVector>& element_size_exterior,
330  const size_t perpendicular_num_points_exterior) const noexcept {
331  const auto penalty = ::elliptic::dg::NumericalFluxes::penalty(
332  min(get(element_size_interior), get(element_size_exterior)),
333  std::max(perpendicular_num_points_interior,
334  perpendicular_num_points_exterior),
335  penalty_parameter_);
336 
337  const auto assemble_numerical_fluxes = [&penalty](
338  const auto numerical_flux_for_field,
339  const auto numerical_flux_for_auxiliary_field,
340  const auto& normal_dot_auxiliary_flux_interior,
341  const auto& normal_dot_div_auxiliary_flux_interior,
342  const auto& ndot_ndot_aux_flux_interior,
343  const auto& minus_normal_dot_auxiliary_flux_exterior,
344  const auto& minus_normal_dot_div_aux_flux_exterior,
345  const auto& ndot_ndot_aux_flux_exterior) noexcept {
346  for (auto it = numerical_flux_for_auxiliary_field->begin();
347  it != numerical_flux_for_auxiliary_field->end(); it++) {
348  const auto index =
349  numerical_flux_for_auxiliary_field->get_tensor_index(it);
350  // We are working with the n.F_v(u) computed on either side of the
351  // interface, so (assuming the normal is independent of the dynamic
352  // variables) the data we get from the other element contains _minus_
353  // the normal from this element. So we cancel the minus sign when
354  // computing the average here.
355  *it = 0.5 * (normal_dot_auxiliary_flux_interior.get(index) -
356  minus_normal_dot_auxiliary_flux_exterior.get(index));
357  }
358  for (auto it = numerical_flux_for_field->begin();
359  it != numerical_flux_for_field->end(); it++) {
360  const auto index = numerical_flux_for_field->get_tensor_index(it);
361  *it = 0.5 * (normal_dot_div_auxiliary_flux_interior.get(index) -
362  minus_normal_dot_div_aux_flux_exterior.get(index)) -
363  penalty * (ndot_ndot_aux_flux_interior.get(index) -
364  ndot_ndot_aux_flux_exterior.get(index));
365  }
366  };
367  EXPAND_PACK_LEFT_TO_RIGHT(assemble_numerical_fluxes(
368  numerical_flux_for_fields, numerical_flux_for_auxiliary_fields,
369  normal_dot_auxiliary_flux_interiors,
370  normal_dot_div_auxiliary_flux_interiors, ndot_ndot_aux_flux_interiors,
371  minus_normal_dot_auxiliary_flux_exteriors,
372  minus_normal_dot_div_aux_flux_exteriors, ndot_ndot_aux_flux_exteriors));
373  }
374 
375  // This function computes the boundary contributions from Dirichlet boundary
376  // conditions. This data is what remains to be added to the boundaries when
377  // homogeneous (i.e. zero) boundary conditions are assumed in the calculation
378  // of the numerical fluxes, but we wish to impose inhomogeneous (i.e. nonzero)
379  // boundary conditions. Since this contribution does not depend on the
380  // numerical field values, but only on the Dirichlet boundary data, it may be
381  // added as contribution to the source of the elliptic systems. Then, it
382  // remains to solve the homogeneous problem with the modified source.
383  template <typename... FluxesArgs>
384  void compute_dirichlet_boundary(
385  const gsl::not_null<typename ::Tags::NormalDotNumericalFlux<
386  FieldTags>::type*>... numerical_flux_for_fields,
387  const gsl::not_null<typename ::Tags::NormalDotNumericalFlux<
388  AuxiliaryFieldTags>::type*>... numerical_flux_for_auxiliary_fields,
389  const typename FieldTags::type&... dirichlet_fields,
390  const Mesh<Dim>& volume_mesh, const Direction<Dim>& direction,
391  const tnsr::i<DataVector, Dim, Frame::Inertial>& interface_unit_normal,
392  const Scalar<DataVector>& face_normal_magnitude,
393  const FluxesComputer& fluxes_computer,
394  const FluxesArgs&... fluxes_args) const noexcept {
395  const auto penalty = ::elliptic::dg::NumericalFluxes::penalty(
396  2. / get(face_normal_magnitude),
397  volume_mesh.extents(direction.dimension()), penalty_parameter_);
398 
399  // Compute n.F_v(u)
400  auto dirichlet_auxiliary_field_fluxes =
401  make_with_value<Variables<tmpl::list<::Tags::Flux<
402  AuxiliaryFieldTags, tmpl::size_t<Dim>, Frame::Inertial>...>>>(
403  interface_unit_normal,
405  fluxes_computer.apply(
406  make_not_null(&get<::Tags::Flux<AuxiliaryFieldTags, tmpl::size_t<Dim>,
407  Frame::Inertial>>(
408  dirichlet_auxiliary_field_fluxes))...,
409  fluxes_args..., dirichlet_fields...);
410  const auto apply_normal_dot_aux =
411  [&interface_unit_normal, &dirichlet_auxiliary_field_fluxes ](
412  auto auxiliary_field_tag_v,
413  const auto numerical_flux_for_auxiliary_field) noexcept {
414  using auxiliary_field_tag = decltype(auxiliary_field_tag_v);
415  normal_dot_flux(
416  numerical_flux_for_auxiliary_field, interface_unit_normal,
417  get<::Tags::Flux<auxiliary_field_tag, tmpl::size_t<Dim>,
418  Frame::Inertial>>(dirichlet_auxiliary_field_fluxes));
419  };
420  EXPAND_PACK_LEFT_TO_RIGHT(apply_normal_dot_aux(
421  AuxiliaryFieldTags{}, numerical_flux_for_auxiliary_fields));
422 
423  // Compute 2 * sigma * n.F_u(n.F_v(u))
424  auto principal_dirichlet_auxiliary_field_fluxes =
425  make_with_value<Variables<tmpl::list<
427  interface_unit_normal,
429  fluxes_computer.apply(
431  &get<::Tags::Flux<FieldTags, tmpl::size_t<Dim>, Frame::Inertial>>(
432  principal_dirichlet_auxiliary_field_fluxes))...,
433  fluxes_args..., *numerical_flux_for_auxiliary_fields...);
434  const auto assemble_dirichlet_penalty = [
435  &interface_unit_normal, &penalty,
436  &principal_dirichlet_auxiliary_field_fluxes
437  ](auto field_tag_v, const auto numerical_flux_for_field) noexcept {
438  using field_tag = decltype(field_tag_v);
439  normal_dot_flux(
440  numerical_flux_for_field, interface_unit_normal,
441  get<::Tags::Flux<field_tag, tmpl::size_t<Dim>, Frame::Inertial>>(
442  principal_dirichlet_auxiliary_field_fluxes));
443  for (auto it = numerical_flux_for_field->begin();
444  it != numerical_flux_for_field->end(); it++) {
445  *it *= 2 * penalty;
446  }
447  };
449  assemble_dirichlet_penalty(FieldTags{}, numerical_flux_for_fields));
450  }
451 
452  private:
453  double penalty_parameter_ = std::numeric_limits<double>::signaling_NaN();
454 };
455 
456 } // namespace elliptic::dg::NumericalFluxes
domain::push_front
CoordinateMap< SourceFrame, TargetFrame, NewMap, Maps... > push_front(CoordinateMap< SourceFrame, TargetFrame, Maps... > old_map, NewMap new_map) noexcept
Creates a CoordinateMap by prepending the new map to the beginning of the old maps.
elliptic::dg::NumericalFluxes::penalty
DataVector penalty(const DataVector &element_size, const size_t num_points, const double penalty_parameter) noexcept
The penalty factor in internal penalty fluxes.
Definition: FirstOrderInternalPenalty.hpp:54
std::string
DataBoxTag.hpp
EXPAND_PACK_LEFT_TO_RIGHT
#define EXPAND_PACK_LEFT_TO_RIGHT(...)
Expand a parameter pack evaluating the terms from left to right.
Definition: TMPL.hpp:563
Frame::Inertial
Definition: IndexType.hpp:44
Divergence.hpp
get
constexpr Tag::type & get(Variables< TagList > &v) noexcept
Return Tag::type pointing into the contiguous array.
Definition: Variables.hpp:639
db::PrefixTag
Marks an item as being a prefix to another tag.
Definition: Tag.hpp:66
elliptic::dg::NumericalFluxes
Numerical fluxes for elliptic systems.
Definition: FirstOrderInternalPenalty.hpp:36
Options.hpp
FaceNormal.hpp
Tags.hpp
square
constexpr decltype(auto) square(const T &x)
Compute the square of x
Definition: ConstantExpressions.hpp:55
domain::Tags::Mesh
The computational grid of the Element in the DataBox.
Definition: Tags.hpp:107
db::SimpleTag
Tags for the DataBox inherit from this type.
Definition: Tag.hpp:23
Tags::NormalDotFlux
Prefix indicating a boundary unit normal vector dotted into the flux.
Definition: Prefixes.hpp:96
Tags::Flux
Prefix indicating a flux.
Definition: Prefixes.hpp:40
Direction
Definition: Direction.hpp:23
Tags::div
Prefix indicating the divergence.
Definition: Divergence.hpp:45
Tags::Normalized
Definition: Magnitude.hpp:138
make_with_value
std::remove_const_t< R > make_with_value(const T &input, const ValueType &value) noexcept
Given an object of type T, create an object of type R whose elements are initialized to value.
Definition: MakeWithValue.hpp:88
cstddef
TensorMetafunctions::remove_first_index
::Tensor< typename Tensor::type, tmpl::pop_front< typename Tensor::symmetry >, tmpl::pop_front< typename Tensor::index_list > > remove_first_index
remove the first index of a tensor
Definition: Metafunctions.hpp:113
DataVector
Stores a collection of function values.
Definition: DataVector.hpp:42
std::numeric_limits::signaling_NaN
T signaling_NaN(T... args)
Metafunctions.hpp
elliptic::dg::NumericalFluxes::FirstOrderInternalPenalty
The internal penalty flux for first-order elliptic equations.
Definition: FirstOrderInternalPenalty.hpp:168
std::min
T min(T... args)
Variables.hpp
Mesh
Holds the number of grid points, basis, and quadrature in each direction of the computational grid.
Definition: Mesh.hpp:47
Scalar
Tensor< T, Symmetry<>, index_list<> > Scalar
Definition: TypeAliases.hpp:21
Tags::Magnitude
Definition: Magnitude.hpp:99
Gsl.hpp
Options::String
const char *const String
The string used in option structs.
Definition: Options.hpp:32
domain::Tags::Direction
Definition: Tags.hpp:347
Tensor.hpp
make_not_null
gsl::not_null< T * > make_not_null(T *ptr) noexcept
Construct a not_null from a pointer. Often this will be done as an implicit conversion,...
Definition: Gsl.hpp:880
Prefixes.hpp
tt::ConformsTo
Indicate a class conforms to the Protocol.
Definition: ProtocolHelpers.hpp:22
TMPL.hpp
gsl::not_null
Require a pointer to not be a nullptr
Definition: Gsl.hpp:183
string