26#ifndef __TIME_SLICE_H_
27#define __TIME_SLICE_H_
160#include "star_rot_dirac.h"
161#include "evolution.h"
299 bool partial_read,
int depth_in = 3) ;
332 assert ((0<= ord)&&(ord < 4)) ;
402 ostream& ost = cout,
bool verb=
true)
const ;
420 ostream& ost = cout,
bool verb=
true)
const ;
445 const Scalar* energy_density = 0x0,
446 ostream& ost = cout,
bool verb=
true)
const ;
458 virtual ostream&
operator>>(ostream& )
const ;
470 void save(
const char* rootname)
const ;
481 virtual void sauve(FILE* fich,
bool partial_save)
const ;
485ostream& operator<<(ostream& ,
const Time_slice& ) ;
610 const Scalar& trk_in,
int depth_in = 3) ;
664 bool partial_read,
int depth_in = 3) ;
863 virtual const Vector&
vec_X(
int method_poisson=6)
const ;
873 int iter_max = 200,
double precis = 1.e-12,
874 double relax = 0.8,
int methode_poisson = 6) ;
881 virtual void set_AB_hata(
const Scalar& A_in,
const Scalar& B_in) ;
917 const Scalar& trk_point,
double pdt,
double precis = 1.e-12,
918 int method_poisson_vect = 6,
const char* graph_device = 0x0,
919 const Scalar* ener_dens = 0x0,
const Vector* mom_dens = 0x0,
920 const Scalar* trace_stress = 0x0 ) ;
939 void check_psi_dot(
Tbl& tlnpsi_dot,
Tbl& tdiff,
Tbl& tdiff_rel)
const ;
945 virtual ostream&
operator>>(ostream& )
const ;
955 virtual void sauve(FILE* fich,
bool partial_save)
const ;
1078 bool partial_read,
int depth_in = 3) ;
1141 const Scalar& trk_point,
double pdt,
double precis = 1.e-12,
1142 int method_poisson_vect = 6,
const char* graph_device = 0x0,
1143 const Scalar* ener_dens = 0x0,
const Vector* mom_dens = 0x0,
1144 const Scalar* trace_stress = 0x0 ) ;
1202 const Scalar* trace_stress=0x0)
const ;
1237 void evolve(
double pdt,
int nb_time_steps,
int niter_elliptic,
1238 double relax_elliptic,
int check_mod,
int save_mod,
1239 int method_poisson_vect = 6,
int nopause = 1,
1240 const char* graph_device = 0x0,
bool verbose=
true,
1241 const Scalar* ener_euler = 0x0,
1242 const Vector* mom_euler = 0x0,
const Scalar* s_euler = 0x0,
1263 void initialize_sources_copy()
const ;
1333 virtual ostream&
operator>>(ostream& )
const ;
1343 virtual void sauve(FILE* fich,
bool partial_save)
const ;
Vectorial bases (triads) with respect to which the tensorial components are defined.
Time evolution with full storage (*** under development ***).
Time evolution with partial storage (*** under development ***).
Base class for coordinate mappings.
Flat metric for tensor calculation.
Metric for tensor calculation.
Tensor field of valence 0 (or component of a tensorial field).
Class for relativistic rotating stars in Dirac gauge and maximal slicing.
Transverse symmetric tensors of rank 2.
Transverse and traceless symmetric tensors of rank 2.
Class intended to describe valence-2 symmetric tensors.
Spacelike time slice of a 3+1 spacetime with conformal decomposition.
virtual const Scalar & npsi() const
Factor at the current time step (jtime ).
Evolution_std< Sym_tensor > hh_evol
Values at successive time steps of the components .
virtual ~Time_slice_conf()
Destructor.
void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
virtual const Sym_tensor & gam_dd() const
Induced metric (covariant components ) at the current time step (jtime )
virtual const Sym_tensor & hh(Param *=0x0, Param *=0x0) const
Deviation of the conformal metric from the flat metric : .
virtual void set_psi_del_n(const Scalar &psi_in)
Sets the conformal factor relating the physical metric to the conformal one: .
virtual void set_hata_TT(const Sym_tensor_tt &hata_tt)
Sets the TT part of (see member hata_evol ).
virtual const Sym_tensor & k_dd() const
Extrinsic curvature tensor (covariant components ) at the current time step (jtime )
virtual void set_npsi_del_psi(const Scalar &npsi_in)
Sets the factor at the current time step (jtime ) and deletes the value of .
virtual const Scalar & nn() const
Lapse function N at the current time step (jtime )
virtual void set_psi_del_npsi(const Scalar &psi_in)
Sets the conformal factor relating the physical metric to the conformal one: .
virtual const Sym_tensor & hata() const
Conformal representation of the traceless part of the extrinsic curvature: .
virtual const Sym_tensor & k_uu() const
Extrinsic curvature tensor (contravariant components ) at the current time step (jtime )
const Scalar & ln_psi() const
Logarithm of at the current time step (jtime ).
Evolution_std< Sym_tensor > hata_evol
Values at successive time steps of the components .
virtual const Vector & hdirac() const
Vector which vanishes in Dirac gauge.
virtual const Scalar & A_hata() const
Returns the potential A of .
virtual const Vector & vec_X(int method_poisson=6) const
Vector representing the longitudinal part of .
void operator=(const Time_slice_conf &)
Assignment to another Time_slice_conf.
virtual const Scalar & psi() const
Conformal factor relating the physical metric to the conformal one: .
Metric * p_tgamma
Pointer on the conformal metric at the current time step (jtime)
virtual const Scalar & B_hata() const
Returns the potential of .
virtual void del_deriv() const
Deletes all the derived quantities.
Evolution_std< Scalar > npsi_evol
Values at successive time steps of the factor .
Vector * p_vec_X
Pointer on the vector representing the longitudinal part of .
virtual Sym_tensor aa() const
Conformal representation of the traceless part of the extrinsic curvature: .
virtual void set_hata(const Sym_tensor &hata_in)
Sets the conformal representation of the traceless part of the extrinsic curvature: .
Vector * p_hdirac
Pointer on the vector (which vanishes in Dirac gauge), at the current time step (jtime).
const Scalar & psi4() const
Factor at the current time step (jtime ).
virtual void set_npsi_del_n(const Scalar &npsi_in)
Sets the factor at the current time step (jtime ) and deletes the value of N.
Evolution_std< Scalar > psi_evol
Values at successive time steps of the conformal factor relating the physical metric to the conform...
virtual void set_hata_from_XAB(Param *par_bc=0x0, Param *par_mat=0x0)
Sets the conformal representation of the traceless part of the extrinsic curvature from its potentia...
virtual const Metric & tgam() const
Conformal metric Returns the value at the current time step (jtime ).
virtual const Sym_tensor & gam_uu() const
Induced metric (contravariant components ) at the current time step (jtime )
Evolution_std< Scalar > A_hata_evol
Potential A associated with the symmetric tensor .
const Metric_flat & ff
Pointer on the flat metric with respect to which the conformal decomposition is performed.
void compute_X_from_momentum_constraint(const Vector &hat_S, const Sym_tensor_tt &hata_tt, int iter_max=200, double precis=1.e-12, double relax=0.8, int methode_poisson=6)
Computes the vector from the conformally-rescaled momentum , using the momentum constraint.
Scalar * p_ln_psi
Pointer on the logarithm of at the current time step (jtime)
Evolution_std< Scalar > B_hata_evol
Potential associated with the symmetric tensor .
virtual void set_hh(const Sym_tensor &hh_in)
Sets the deviation of the conformal metric from the flat metric : .
virtual const Scalar & trk() const
Trace K of the extrinsic curvature at the current time step (jtime )
Scalar * p_psi4
Pointer on the factor at the current time step (jtime)
Spacelike time slice of a 3+1 spacetime.
virtual const Sym_tensor & k_dd() const
Extrinsic curvature tensor (covariant components ) at the current time step (jtime )
int jtime
Time step index of the latest slice.
Evolution_std< Scalar > trk_evol
Values at successive time steps of the trace K of the extrinsic curvature.
Metric * p_gamma
Pointer on the induced metric at the current time step (jtime)
Evolution_std< Sym_tensor > k_dd_evol
Values at successive time steps of the covariant components of the extrinsic curvature tensor .
virtual void sauve(FILE *fich, bool partial_save) const
Total or partial saves in a binary file.
void operator=(const Time_slice &)
Assignment to another Time_slice.
Evolution_full< Tbl > adm_mass_evol
ADM mass at each time step, since the creation of the slice.
virtual const Scalar & trk() const
Trace K of the extrinsic curvature at the current time step (jtime )
virtual const Scalar & nn() const
Lapse function N at the current time step (jtime )
void save(const char *rootname) const
Saves in a binary file.
virtual const Sym_tensor & gam_dd() const
Induced metric (covariant components ) at the current time step (jtime )
void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Evolution_std< Sym_tensor > gam_dd_evol
Values at successive time steps of the covariant components of the induced metric .
friend ostream & operator<<(ostream &, const Time_slice &)
Display.
Tbl check_hamiltonian_constraint(const Scalar *energy_density=0x0, ostream &ost=cout, bool verb=true) const
Checks the level at which the hamiltonian constraint is verified.
Evolution_std< Sym_tensor > gam_uu_evol
Values at successive time steps of the contravariant components of the induced metric .
int scheme_order
Order of the finite-differences scheme for the computation of time derivatives.
virtual const Vector & beta() const
shift vector at the current time step (jtime )
virtual ~Time_slice()
Destructor.
virtual void del_deriv() const
Deletes all the derived quantities.
int get_latest_j() const
Gets the latest value of time step index.
const Metric & gam() const
Induced metric at the current time step (jtime )
virtual const Sym_tensor & k_uu() const
Extrinsic curvature tensor (contravariant components ) at the current time step (jtime )
int depth
Number of stored time slices.
const Evolution_std< double > & get_time() const
Gets the time coordinate t at successive time steps.
Evolution_std< Scalar > n_evol
Values at successive time steps of the lapse function N.
virtual const Sym_tensor & gam_uu() const
Induced metric (contravariant components ) at the current time step (jtime )
virtual double adm_mass() const
Returns the ADM mass (geometrical units) at the current step.
int get_scheme_order() const
Gets the order of the finite-differences scheme.
Evolution_std< double > the_time
Time label of each slice.
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator<<).
Evolution_std< Sym_tensor > k_uu_evol
Values at successive time steps of the contravariant components of the extrinsic curvature tensor .
Tbl check_dynamical_equations(const Sym_tensor *strain_tensor=0x0, const Scalar *energy_density=0x0, ostream &ost=cout, bool verb=true) const
Checks the level at which the dynamical equations are verified.
void set_scheme_order(int ord)
Sets the order of the finite-differences scheme.
Evolution_std< Vector > beta_evol
Values at successive time steps of the shift vector .
Tbl check_momentum_constraint(const Vector *momentum_density=0x0, ostream &ost=cout, bool verb=true) const
Checks the level at which the momentum constraints are verified.
Spacelike time slice of a 3+1 spacetime with conformal decomposition in the maximal slicing and Dirac...
virtual const Scalar & trh() const
Computes the trace h, with respect to the flat metric ff , of .
virtual double adm_mass() const
Returns the ADM mass at (geometrical units) the current step.
virtual const Sym_tensor & hh(Param *par_bc=0x0, Param *par_mat=0x0) const
Deviation of the conformal metric from the flat metric : .
virtual Scalar solve_psi(const Scalar *ener_dens=0x0) const
Solves the elliptic equation for the conformal factor $\Psi$ (Hamiltonian constraint).
Evolution_std< Scalar > source_B_hh_evol
The potential of the source of equation for .
Evolution_std< Scalar > source_A_hh_evol
The A potential of the source of equation for .
void evolve(double pdt, int nb_time_steps, int niter_elliptic, double relax_elliptic, int check_mod, int save_mod, int method_poisson_vect=6, int nopause=1, const char *graph_device=0x0, bool verbose=true, const Scalar *ener_euler=0x0, const Vector *mom_euler=0x0, const Scalar *s_euler=0x0, const Sym_tensor *strain_euler=0x0)
Time evolution by resolution of Einstein equations.
virtual const Scalar & B_hh() const
Returns the potential of .
virtual void initial_data_cts(const Sym_tensor &uu, const Scalar &trk_in, const Scalar &trk_point, double pdt, double precis=1.e-12, int method_poisson_vect=6, const char *graph_device=0x0, const Scalar *ener_dens=0x0, const Vector *mom_dens=0x0, const Scalar *trace_stress=0x0)
Computes valid initial data by solving the constraint equations in the conformal thin-sandwich approa...
Evolution_std< Scalar > A_hh_evol
The A potential of .
virtual ~Tslice_dirac_max()
Destructor.
virtual const Scalar & trk() const
Trace K of the extrinsic curvature at the current time step (jtime ).
void hh_det_one(int j, Param *par_bc=0x0, Param *par_mat=0x0) const
Computes from the values of A and and using the condition , which fixes the trace of .
Evolution_std< Scalar > source_B_hata_evol
The potential of the source of equation for .
virtual void set_hh(const Sym_tensor &hh_in)
Sets the deviation of the conformal metric from the flat metric : .
void compute_sources(const Sym_tensor *strain_tensor=0x0) const
Computes the sources source_A_XXX_evol and source_B_XXX_evol , for the solution of the evolution equa...
virtual Scalar solve_npsi(const Scalar *ener_dens=0x0, const Scalar *trace_stress=0x0) const
Solves the elliptic equation for (maximal slicing condition + Hamiltonian constraint)
virtual void set_khi_mu(const Scalar &khi_in, const Scalar &mu_in)
Sets the potentials and of the TT part of (see the documentation of Sym_tensor_tt for details).
virtual void sauve(FILE *fich, bool partial_save) const
Total or partial saves in a binary file.
virtual void set_trh(const Scalar &trh_in)
Sets the trace, with respect to the flat metric ff , of .
virtual const Scalar & A_hh() const
Returns the potential A of .
virtual const Vector & hdirac() const
Vector which vanishes in Dirac gauge.
virtual Vector solve_beta(int method=6) const
Solves the elliptic equation for the shift vector from (Eq.
virtual void set_AB_hh(const Scalar &A_in, const Scalar &B_in)
Sets the potentials A and of the TT part of (see the documentation of Sym_tensor for details).
Evolution_std< Scalar > B_hh_evol
The potential of .
Evolution_std< Scalar > source_A_hata_evol
The potential A of the source of equation for .
Tensor field of valence 1.