athena.hpp

Go to the documentation of this file.

#ifndef ATHENA_HPP_
#define ATHENA_HPP_
//========================================================================================
// Athena++ astrophysical MHD code
// Copyright(C) 2014 James M. Stone <jmstone@princeton.edu> and other code contributors
// Licensed under the 3-clause BSD License, see LICENSE file for details
//========================================================================================
//  \brief contains Athena++ general purpose types, structures, enums, etc.
 
// C headers
 
// C++ headers
#include <cmath>
#include <cstdint>  // std::int64_t
 
// Athena++ headers
#include "athena_arrays.hpp"
#include "defs.hpp"
 
// primitive type alias that allows code to run with either floats or doubles
#if SINGLE_PRECISION_ENABLED
using Real = float;
#ifdef MPI_PARALLEL
#define MPI_ATHENA_REAL MPI_FLOAT
#define MPI_ATHENA_REAL_INT MPI_FLOAT_INT
#endif
#else
using Real = double;
#ifdef MPI_PARALLEL
#define MPI_ATHENA_REAL MPI_DOUBLE
#define MPI_ATHENA_REAL_INT MPI_DOUBLE_INT
#endif
#endif
 
// for OpenMP 4.0 SIMD vectorization, control width of SIMD lanes
#if defined(__AVX512F__)
#define SIMD_WIDTH 8
#elif defined(__AVX__)
#define SIMD_WIDTH 4
#elif defined(__SSE2__)
#define SIMD_WIDTH 2
#else
#define SIMD_WIDTH 4
#endif
 
#define CACHELINE_BYTES 64
 
// forward declarations needed for function pointer type aliases
class MeshBlock;
class Coordinates;
class ParameterInput;
class HydroDiffusion;
class FieldDiffusion;
 
//--------------------------------------------------------------------------------------
//  \brief stores logical location and level of MeshBlock
 
struct LogicalLocation { // aggregate and POD type
  // These values can exceed the range of std::int32_t even if the root grid has only a
  // single MeshBlock if >30 levels of AMR are used, since the corresponding max index =
  // 1*2^31 > INT_MAX = 2^31 -1 for most 32-bit signed integer type impelementations
  std::int64_t lx1, lx2, lx3;
  int level;
 
  // operators useful for sorting
  bool operator==(LogicalLocation &ll) {
    return ((ll.level == level) && (ll.lx1 == lx1) && (ll.lx2 == lx2) && (ll.lx3 == lx3));
  }
  static bool Lesser(const LogicalLocation &left, const LogicalLocation &right) {
    return left.level < right.level;
  }
  static bool Greater(const LogicalLocation & left, const LogicalLocation &right) {
    return left.level > right.level;
  }
};
 
//----------------------------------------------------------------------------------------
//  \brief physical size and number of cells in a Mesh or a MeshBlock
 
struct RegionSize {  // aggregate and POD type; do NOT reorder member declarations:
  Real x1min, x2min, x3min;
  Real x1max, x2max, x3max;
  Real x1rat, x2rat, x3rat; // ratio of dxf(i)/dxf(i-1)
  // the size of the root grid or a MeshBlock should not exceed std::int32_t limits
  int nx1, nx2, nx3;        // number of active cells (not including ghost zones)
};
 
//---------------------------------------------------------------------------------------
//  \brief container for face-centered fields
 
struct FaceField {
  AthenaArray<Real> x1f, x2f, x3f;
  FaceField() = default;
  FaceField(int ncells3, int ncells2, int ncells1,
            AthenaArray<Real>::DataStatus init=AthenaArray<Real>::DataStatus::allocated) :
      x1f(ncells3, ncells2, ncells1+1, init), x2f(ncells3, ncells2+1, ncells1, init),
      x3f(ncells3+1, ncells2, ncells1, init) {}
};
 
//----------------------------------------------------------------------------------------
//  \brief container for edge-centered fields
 
struct EdgeField {
  AthenaArray<Real> x1e, x2e, x3e;
  EdgeField() = default;
  EdgeField(int ncells3, int ncells2, int ncells1,
            AthenaArray<Real>::DataStatus init=AthenaArray<Real>::DataStatus::allocated) :
      x1e(ncells3+1, ncells2+1, ncells1, init), x2e(ncells3+1, ncells2, ncells1+1, init),
      x3e(ncells3, ncells2+1, ncells1+1, init) {}
};
 
//----------------------------------------------------------------------------------------
// enums used everywhere
// (not specifying underlying integral type (C++11) for portability & performance)
 
// TODO(felker): C++ Core Guidelines Enum.5: Don’t use ALL_CAPS for enumerators
// (avoid clashes with preprocessor macros). Enumerated type definitions in this file and:
// athena_fft.hpp, io_wrapper.hpp, bvals.hpp, hydro_diffusion.hpp, field_diffusion.hpp,
// task_list.hpp, ???
 
//------------------
// named, weakly typed / unscoped enums:
//------------------
 
// enumerators only used for indexing AthenaArray and regular arrays; typename and
// explicitly specified enumerator values aare unnecessary, but provided for clarity:
 
// array indices for conserved: density, momemtum, total energy, face-centered field
//enum {NMASS = 1+NPHASE*NVAPOR};
enum ConsIndex {IDN=0, IM1=1+NVAPOR, IM2=2+NVAPOR, IM3=3+NVAPOR, IEN=4+NVAPOR};
enum MagneticIndex {IB1=0, IB2=1, IB3=2};
 
// array indices for 1D primitives: velocity, transverse components of field
enum PrimIndex {IVX=1+NVAPOR, IVY=2+NVAPOR, IVZ=3+NVAPOR, IPR=4+NVAPOR, 
  IBY=(NHYDRO), IBZ=((NHYDRO)+1), IV1=IVX, IV2=IVY, IV3=IVZ};
 
// array indices for face-centered electric fields returned by Riemann solver
enum ElectricIndex {X1E2=0, X1E3=1, X2E3=0, X2E1=1, X3E1=0, X3E2=1};
 
// array indices for metric and triangular matrices in GR
enum MetricIndex {I00=0, I01=1, I02=2, I03=3, I11=4, I12=5, I13=6, I22=7, I23=8, I33=9,
                  NMETRIC=10};
enum TriangleIndex {T00=0, T10=1, T11=2, T20=3, T21=4, T22=5, T30=6, T31=7, T32=8, T33=9,
                    NTRIANGULAR=10};
 
// enumerator types that are used for variables and function parameters:
 
// needed for arrays dimensioned over grid directions
// enumerator type only used in Mesh::EnrollUserMeshGenerator()
enum CoordinateDirection {X1DIR=0, X2DIR=1, X3DIR=2};
 
//------------------
// strongly typed / scoped enums (C++11):
//------------------
// KGF: Except for the 2x MG* enums, these may be unnessary w/ the new class inheritance
// Now, only passed to BoundaryVariable::InitBoundaryData(); could replace w/ bool switch
enum class BoundaryQuantity {cc, fc, cc_flcor, fc_flcor, mggrav, mggrav_f};
enum class HydroBoundaryQuantity {cons, prim};
enum class BoundaryCommSubset {mesh_init, gr_amr, all};
// TODO(felker): consider generalizing/renaming to QuantityFormulation
enum class FluidFormulation {evolve, background, disabled}; // rename background -> fixed?
enum class UserHistoryOperation {sum, max, min};
enum class VariableType {prim, cons, chem};
 
//----------------------------------------------------------------------------------------
// function pointer prototypes for user-defined modules set at runtime
 
using BValFunc = void (*)(
    MeshBlock *pmb, Coordinates *pco, AthenaArray<Real> &prim, FaceField &b,
    Real time, Real dt,
    int is, int ie, int js, int je, int ks, int ke, int ngh);
using AMRFlagFunc = int (*)(MeshBlock *pmb);
using MeshGenFunc = Real (*)(Real x, RegionSize rs);
using SrcTermFunc = void (*)(
    MeshBlock *pmb, const Real time, const Real dt,
    const AthenaArray<Real> &prim, const AthenaArray<Real> &bcc, AthenaArray<Real> &cons);
using TimeStepFunc = Real (*)(MeshBlock *pmb);
using HistoryOutputFunc = Real (*)(MeshBlock *pmb, int iout);
using MetricFunc = void (*)(
    Real x1, Real x2, Real x3, ParameterInput *pin,
    AthenaArray<Real> &g, AthenaArray<Real> &g_inv,
    AthenaArray<Real> &dg_dx1, AthenaArray<Real> &dg_dx2, AthenaArray<Real> &dg_dx3);
using MGBoundaryFunc = void (*)(
    AthenaArray<Real> &dst,Real time, int nvar,
    int is, int ie, int js, int je, int ks, int ke, int ngh,
    Real x0, Real y0, Real z0, Real dx, Real dy, Real dz);
using ViscosityCoeffFunc = void (*)(
    HydroDiffusion *phdif, MeshBlock *pmb,
    const  AthenaArray<Real> &w, const AthenaArray<Real> &bc,
    int is, int ie, int js, int je, int ks, int ke);
using ConductionCoeffFunc = void (*)(
    HydroDiffusion *phdif, MeshBlock *pmb,
    const AthenaArray<Real> &w, const AthenaArray<Real> &bc,
    int is, int ie, int js, int je, int ks, int ke);
using FieldDiffusionCoeffFunc = void (*)(
    FieldDiffusion *pfdif, MeshBlock *pmb,
    const AthenaArray<Real> &w,
    const AthenaArray<Real> &bmag,
    int is, int ie, int js, int je, int ks, int ke);
 
#endif // ATHENA_HPP_