Skip to content

Commit

Permalink
Update initial data class
Browse files Browse the repository at this point in the history
  • Loading branch information
dinatraykova committed Mar 28, 2024
1 parent 2171d67 commit ee5bbb4
Show file tree
Hide file tree
Showing 3 changed files with 78 additions and 56 deletions.
58 changes: 8 additions & 50 deletions Examples/Fluid_Kerr/InitialFluidData.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -17,7 +17,6 @@
#include "simd.hpp"

//! Class which sets the initial fluid matter config
// aka Kevin data
class InitialFluidData
{
public:
Expand All @@ -33,61 +32,20 @@ class InitialFluidData
};

//! The constructor
InitialFluidData(params_t a_params, double a_dx)
: m_dx(a_dx), m_params(a_params)
{
}
InitialFluidData(params_t a_params, double a_dx);
// : m_dx(a_dx), m_params(a_params)
//{
//}

//! Function to compute the value of all the initial vars on the grid
template <class data_t> void compute(Cell<data_t> current_cell) const
{
// load vars
FluidCCZ4RHS<PerfectFluid<>>::Vars<data_t> vars;
VarsTools::assign(vars, 0.);
// where am i?
Coordinates<data_t> coords(current_cell, m_dx, m_params.center);
data_t rr = coords.get_radius();
data_t rr2 = rr * rr;

data_t x = coords.x;
double y = coords.y;
double z = coords.z;

// Tensor<1, data_t> vi, Sj;
data_t chi_regularised = simd_max(vars.chi, 1e-6);

// vi[0] = 0.;
// vi[1] = 0.;
// vi[2] = 0.;

// calculate the field value
vars.rho = m_params.rho0 * (exp(-pow(rr / m_params.awidth, 2.0))) +
m_params.delta;
data_t v2 = 0.;
FOR(i, j)
v2 += vars.h[i][j] * vars.vi[i] * vars.vi[j] / chi_regularised;
vars.eps = 0.;
data_t P = vars.rho * (1. + vars.eps) / 3.;
data_t WW = 1. / (1. - v2);
data_t hh = 1. + vars.eps + P / vars.rho;

vars.D = vars.rho * sqrt(WW);
vars.tau = vars.rho * hh * WW - P - vars.D;
FOR(i)
{
vars.Sj[i] = 0.;
FOR(j)
vars.Sj[i] += vars.rho * hh * WW * vars.h[i][j] * vars.vi[j] /
chi_regularised;
}

// store the vars
current_cell.store_vars(vars);
}
//! The constructor
template <class data_t> void compute(Cell<data_t> current_cell) const;

protected:
double m_dx;
const params_t m_params;
};

#include "InitialFluidData.impl.hpp"

#endif /* INITIALFLUIDDATA_HPP_ */
67 changes: 67 additions & 0 deletions Examples/Fluid_Kerr/InitialFluidData.impl.hpp
Original file line number Diff line number Diff line change
@@ -0,0 +1,67 @@
/* GRChombo
* Copyright 2012 The GRChombo collaboration.
* Please refer to LICENSE in GRChombo's root directory.
*/

#if !defined(INITIALFLUIDDATA_HPP_)
#error "This file should only be included through PrimitiveRecovery.hpp"
#endif

#ifndef INITIALFLUIDDATA_IMPL_HPP_
#define INITIALFLUIDDATA_IMPL_HPP_

inline InitialFluidData::InitialFluidData(params_t a_params, double a_dx)
: m_dx(a_dx), m_params(a_params)
{
}

template <class data_t>
void InitialFluidData::compute(Cell<data_t> current_cell) const
{
// Set up EoS
// data_t P_over_rho = 0.;
// my_eos.compute_eos(P_over_rho, vars);
// load vars
FluidCCZ4RHS<PerfectFluid<>>::Vars<data_t> vars;
VarsTools::assign(vars, 0.);
// where am i?
Coordinates<data_t> coords(current_cell, m_dx, m_params.center);
data_t rr = coords.get_radius();
data_t rr2 = rr * rr;

data_t x = coords.x;
double y = coords.y;
double z = coords.z;

// Tensor<1, data_t> vi, Sj;
data_t chi_regularised = simd_max(vars.chi, 1e-6);

// vi[0] = 0.;
// vi[1] = 0.;
// vi[2] = 0.;

// calculate the field value
vars.rho =
m_params.rho0 * (exp(-pow(rr / m_params.awidth, 2.0))) + m_params.delta;
data_t v2 = 0.;
FOR(i, j) v2 += vars.h[i][j] * vars.vi[i] * vars.vi[j] / chi_regularised;
vars.eps = 0.;
data_t P_over_rho = vars.rho * (1. + vars.eps) / 3.;
data_t WW = 1. / (1. - v2);
data_t hh = 1. + vars.eps + P_over_rho;

vars.D = vars.rho * sqrt(WW);
vars.tau = vars.rho * (hh * WW - P_over_rho) - vars.D;
FOR(i)
{
vars.Sj[i] = 0.;
FOR(j)
vars.Sj[i] +=
vars.rho * hh * WW * vars.h[i][j] * vars.vi[j] / chi_regularised;
}

// store the vars
current_cell.store_vars(vars);
}

#endif /* INITIALFLUIDDATA_IMPL_HPP_ */
9 changes: 3 additions & 6 deletions Examples/Fluid_Kerr/PerfectFluid.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -25,13 +25,10 @@
This class is an example of a matter_t object which calculates the
matter type specific elements for the RHS update and the evaluation
of the constraints. This includes the Energy Momentum Tensor, and
the matter evolution terms. In this case, a scalar field,
the matter elements are phi and (minus) its conjugate momentum, Pi.
It is templated over a potential function potential_t which the
the matter evolution terms.
It is templated over an equation of state function eos_t, which the
user must specify in a class, although a default is provided which
sets dVdphi and V_of_phi to zero.
It assumes minimal coupling of the field to gravity.
\sa MatterCCZ4(), ConstraintsMatter()
sets P_over_rho to (1 + eps) / 3.
*/
template <class eos_t = DefaultEoS> class PerfectFluid
{
Expand Down

0 comments on commit ee5bbb4

Please sign in to comment.