energymomentumtensor.h

Go to the documentation of this file.

/*
 *
 *    Copyright (c) 2015-2018
 *      SMASH Team
 *
 *    GNU General Public License (GPLv3 or later)
 *
 */

#ifndef SRC_INCLUDE_ENERGYMOMENTUMTENSOR_H_
#define SRC_INCLUDE_ENERGYMOMENTUMTENSOR_H_

#include <cmath>

#include <array>
#include <string>

#include "fourvector.h"
#include "particledata.h"

namespace smash {

class EnergyMomentumTensor {
 public:
  typedef std::array<double, 10> tmn_type;
  EnergyMomentumTensor() { Tmn_.fill(0.); }

  explicit EnergyMomentumTensor(const tmn_type &Tmn) {
    for (size_t i = 0; i < 10; i++) {
      Tmn_[i] = Tmn[i];
    }
  }

  double operator[](std::size_t i) const { return Tmn_[i]; }

  static std::int8_t tmn_index(std::int8_t mu, std::int8_t nu) {
    // clang-format off
    constexpr std::array<std::int8_t, 16> indices = {0, 1, 2, 3,
                                                     1, 4, 5, 6,
                                                     2, 5, 7, 8,
                                                     3, 6, 8, 9};
    // clang-format on
    if (mu < 4 && nu < 4 && mu >= 0 && nu >= 0) {
      return indices[mu + 4 * nu];
    } else {
      throw std::invalid_argument("Invalid indices: " + std::to_string(mu) +
                                  ", " + std::to_string(nu));
    }
  }

  EnergyMomentumTensor inline operator+=(const EnergyMomentumTensor &Tmn0);
  EnergyMomentumTensor inline operator-=(const EnergyMomentumTensor &Tmn0);
  EnergyMomentumTensor inline operator*=(double a);
  EnergyMomentumTensor inline operator/=(double a);

  using iterator = tmn_type::iterator;
  using const_iterator = tmn_type::const_iterator;

  FourVector landau_frame_4velocity() const;

  EnergyMomentumTensor boosted(const FourVector &u) const;

  void add_particle(const FourVector &mom);
  void add_particle(const ParticleData &p, double factor);
  void add_particle_for_derivatives(const ParticleData &, double,
                                    ThreeVector){};

  iterator begin() { return Tmn_.begin(); }

  iterator end() { return Tmn_.end(); }

  const_iterator begin() const { return Tmn_.begin(); }
  const_iterator end() const { return Tmn_.end(); }

  const_iterator cbegin() const { return Tmn_.cbegin(); }
  const_iterator cend() const { return Tmn_.cend(); }

 private:
  tmn_type Tmn_;
};

std::ostream &operator<<(std::ostream &, const EnergyMomentumTensor &);

EnergyMomentumTensor inline EnergyMomentumTensor::operator+=(
    const EnergyMomentumTensor &Tmn0) {
  for (size_t i = 0; i < 10; i++) {
    Tmn_[i] += Tmn0[i];
  }
  return *this;
}
EnergyMomentumTensor inline operator+(EnergyMomentumTensor a,
                                      const EnergyMomentumTensor &b) {
  a += b;
  return a;
}

EnergyMomentumTensor inline EnergyMomentumTensor::operator-=(
    const EnergyMomentumTensor &Tmn0) {
  for (size_t i = 0; i < 10; i++) {
    Tmn_[i] -= Tmn0[i];
  }
  return *this;
}
EnergyMomentumTensor inline operator-(EnergyMomentumTensor a,
                                      const EnergyMomentumTensor &b) {
  a -= b;
  return a;
}

EnergyMomentumTensor inline EnergyMomentumTensor::operator*=(const double a) {
  for (size_t i = 0; i < 10; i++) {
    Tmn_[i] *= a;
  }
  return *this;
}
inline EnergyMomentumTensor operator*(EnergyMomentumTensor a, const double b) {
  a *= b;
  return a;
}
inline EnergyMomentumTensor operator*(const double a, EnergyMomentumTensor b) {
  b *= a;
  return b;
}

EnergyMomentumTensor inline EnergyMomentumTensor::operator/=(const double a) {
  for (size_t i = 0; i < 10; i++) {
    Tmn_[i] /= a;
  }
  return *this;
}
EnergyMomentumTensor inline operator/(EnergyMomentumTensor a, const double b) {
  a /= b;
  return a;
}

}  // namespace smash

#endif  // SRC_INCLUDE_ENERGYMOMENTUMTENSOR_H_