Line data    Source code

       1           0 : // Distributed under the MIT License.
       2             : // See LICENSE.txt for details.
       3             : 
       4             : #pragma once
       5             : 
       6             : #include "DataStructures/DataVector.hpp"
       7             : #include "DataStructures/Tensor/Tensor.hpp"
       8             : #include "Utilities/Gsl.hpp"
       9             : 
      10             : namespace grmhd::GhValenciaDivClean {
      11             : 
      12             : /*!
      13             :  * \brief Add in the trace-reversed stress-energy source term to the \f$\Pi\f$
      14             :  * evolved variable of the Generalized Harmonic system.
      15             :  *
      16             :  * \details The only stress energy source term in the Generalized Harmonic
      17             :  * evolution equations is in the equation for \f$\Pi_{a b}\f$:
      18             :  *
      19             :  * \f[
      20             :  * \partial_t \Pi_{ab} + \text{(spatial derivative terms)} =
      21             :  * \text{(GH source terms)}
      22             :  * - 16 \pi \alpha (T_{ab} - \frac{1}{2} g_{a b} T^c{}_c)
      23             :  * \f]
      24             :  *
      25             :  * (note that this function takes as argument the trace-reversed stress energy
      26             :  * tensor)
      27             :  *
      28             :  * This function adds that contribution to the existing value of `dt_pi`. The
      29             :  * spacetime terms in the GH equation should be computed before passing the
      30             :  * `dt_pi` to this function for updating.
      31             :  *
      32             :  * \see `gh::TimeDerivative` for details about the spacetime
      33             :  * part of the time derivative calculation.
      34             :  */
      35           1 : void add_stress_energy_term_to_dt_pi(
      36             :     gsl::not_null<tnsr::aa<DataVector, 3>*> dt_pi,
      37             :     const tnsr::aa<DataVector, 3>& trace_reversed_stress_energy,
      38             :     const Scalar<DataVector>& lapse);
      39             : 
      40             : /*!
      41             :  * \brief Calculate the trace-reversed stress-energy tensor \f$(T_{\mu \nu} -
      42             :  * 1/2 g_{\mu \nu} g^{\lambda \sigma} T_{\lambda \sigma}) \f$ associated with
      43             :  * the matter part of the GRMHD system.
      44             :  *
      45             :  * \details The stress energy tensor is needed to compute the backreaction of
      46             :  * the matter to the spacetime degrees of freedom. The stress energy is
      47             :  * calculated as described in \cite Moesta2013dna :
      48             :  *
      49             :  * \f[
      50             :  * T_{\mu \nu} = \rho h^* u_\mu u_\nu + p^* g_{\mu \nu} - b_\mu b_\nu,
      51             :  * \f]
      52             :  *
      53             :  * where \f$u_\mu\f$ is the four-velocity, \f$\rho\f$ is the rest mass density,
      54             :  * \f$h^*\f$ is the magnetically modified enthalpy, \f$p^*\f$ is the
      55             :  * magnetically modified pressure, and \f$b_{\mu}\f$ is the comoving magnetic
      56             :  * field (note that we deviate from the notation of \cite Moesta2013dna by
      57             :  * denoting the pressure with a lower-case \f$p\f$ instead of an upper-case
      58             :  * \f$P\f$).
      59             :  *
      60             :  * The spatial components of the four velocity \f$u_\mu\f$ are
      61             :  *
      62             :  * \f[
      63             :  * u_i = W v_i,
      64             :  * \f]
      65             :  *
      66             :  * and the time component is
      67             :  *
      68             :  * \f[
      69             :  * u_0 = - \alpha W + \beta^i u_i.
      70             :  * \f].
      71             :  *
      72             :  * The magnetically modified enthalpy is
      73             :  *
      74             :  * \f[
      75             :  * h^* = 1 + \epsilon + (p + b^2) / \rho.
      76             :  * \f]
      77             :  *
      78             :  * The magnetically modified pressure is
      79             :  *
      80             :  * \f[
      81             :  * p^* = p + b^2 / 2.
      82             :  * \f]
      83             :  *
      84             :  *
      85             :  * The comoving magnetic field is computed via
      86             :  *
      87             :  * \f{align}{
      88             :  * b_i &= B_i / W + v_i W v^k B_k\\
      89             :  * b_0 &= - \alpha W v^i B_i + \beta^i b_i
      90             :  * \f}
      91             :  *
      92             :  * Therefore, the trace-reversed stress energy simplifies to
      93             :  *
      94             :  * \f[
      95             :  * (T_{\mu \nu} - \frac{1}{2} g_{\mu \nu} g^{\lambda \sigma} T_{\lambda \sigma})
      96             :  * = \rho h^* u_\mu u_\nu  + \left(\frac{1}{2} \rho h^* - p\right) g_{\mu \nu}
      97             :  *  - b_\mu b_\nu
      98             :  * \f]
      99             :  */
     100           1 : void trace_reversed_stress_energy(
     101             :     gsl::not_null<tnsr::aa<DataVector, 3>*> stress_energy,
     102             :     gsl::not_null<tnsr::a<DataVector, 3>*> four_velocity_one_form_buffer,
     103             :     gsl::not_null<tnsr::a<DataVector, 3>*>
     104             :         comoving_magnetic_field_one_form_buffer,
     105             :     const Scalar<DataVector>& rest_mass_density,
     106             :     const tnsr::i<DataVector, 3, Frame::Inertial>& spatial_velocity_one_form,
     107             :     const tnsr::i<DataVector, 3, Frame::Inertial>& magnetic_field_one_form,
     108             :     const Scalar<DataVector>& magnetic_field_squared,
     109             :     const Scalar<DataVector>& magnetic_field_dot_spatial_velocity,
     110             :     const Scalar<DataVector>& lorentz_factor,
     111             :     const Scalar<DataVector>& one_over_w_squared,
     112             :     const Scalar<DataVector>& pressure,
     113             :     const Scalar<DataVector>& specific_internal_energy,
     114             :     const tnsr::aa<DataVector, 3, Frame::Inertial>& spacetime_metric,
     115             :     const tnsr::I<DataVector, 3, Frame::Inertial>& shift,
     116             :     const Scalar<DataVector>& lapse);
     117             : }  // namespace grmhd::GhValenciaDivClean