LORENE
|
Spacelike time slice of a 3+1 spacetime. More...
#include <time_slice.h>
Public Member Functions | |
Time_slice (const Scalar &lapse_in, const Vector &shift_in, const Sym_tensor &gamma_in, const Sym_tensor &kk_in, int depth_in=3) | |
General constructor (Hamiltonian-like). | |
Time_slice (const Scalar &lapse_in, const Vector &shift_in, const Evolution_std< Sym_tensor > &gamma_in) | |
General constructor (Lagrangian-like). | |
Time_slice (const Map &mp, const Base_vect &triad, int depth_in=3) | |
Constructor as standard time slice of flat spacetime (Minkowski). | |
Time_slice (const Map &mp, const Base_vect &triad, FILE *fich, bool partial_read, int depth_in=3) | |
Constructor from binary file. | |
Time_slice (const Time_slice &) | |
Copy constructor. | |
virtual | ~Time_slice () |
Destructor. | |
void | operator= (const Time_slice &) |
Assignment to another Time_slice . | |
void | set_scheme_order (int ord) |
Sets the order of the finite-differences scheme. | |
int | get_scheme_order () const |
Gets the order of the finite-differences scheme. | |
int | get_latest_j () const |
Gets the latest value of time step index. | |
const Evolution_std< double > & | get_time () const |
Gets the time coordinate t at successive time steps. | |
virtual const Scalar & | nn () const |
Lapse function N at the current time step (jtime ) | |
virtual const Vector & | beta () const |
shift vector ![]() jtime ) | |
const Metric & | gam () const |
Induced metric ![]() jtime ) | |
virtual const Sym_tensor & | gam_dd () const |
Induced metric (covariant components ![]() jtime ) | |
virtual const Sym_tensor & | gam_uu () const |
Induced metric (contravariant components ![]() jtime ) | |
virtual const Sym_tensor & | k_dd () const |
Extrinsic curvature tensor (covariant components ![]() jtime ) | |
virtual const Sym_tensor & | k_uu () const |
Extrinsic curvature tensor (contravariant components ![]() jtime ) | |
virtual const Scalar & | trk () const |
Trace K of the extrinsic curvature at the current time step (jtime ) | |
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. | |
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. | |
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. | |
virtual double | adm_mass () const |
Returns the ADM mass (geometrical units) at the current step. | |
void | save (const char *rootname) const |
Saves in a binary file. | |
Protected Member Functions | |
Time_slice (int depth_in) | |
Special constructor for derived classes. | |
virtual void | del_deriv () const |
Deletes all the derived quantities. | |
void | set_der_0x0 () const |
Sets to 0x0 all the pointers on derived quantities. | |
virtual ostream & | operator>> (ostream &) const |
Operator >> (virtual function called by the operator<<). | |
virtual void | sauve (FILE *fich, bool partial_save) const |
Total or partial saves in a binary file. | |
Protected Attributes | |
int | depth |
Number of stored time slices. | |
int | scheme_order |
Order of the finite-differences scheme for the computation of time derivatives. | |
int | jtime |
Time step index of the latest slice. | |
Evolution_std< double > | the_time |
Time label of each slice. | |
Evolution_std< Sym_tensor > | gam_dd_evol |
Values at successive time steps of the covariant components of the induced metric ![]() | |
Evolution_std< Sym_tensor > | gam_uu_evol |
Values at successive time steps of the contravariant components of the induced metric ![]() | |
Evolution_std< Sym_tensor > | k_dd_evol |
Values at successive time steps of the covariant components of the extrinsic curvature tensor ![]() | |
Evolution_std< Sym_tensor > | k_uu_evol |
Values at successive time steps of the contravariant components of the extrinsic curvature tensor ![]() | |
Evolution_std< Scalar > | n_evol |
Values at successive time steps of the lapse function N. | |
Evolution_std< Vector > | beta_evol |
Values at successive time steps of the shift vector ![]() | |
Evolution_std< Scalar > | trk_evol |
Values at successive time steps of the trace K of the extrinsic curvature. | |
Evolution_full< Tbl > | adm_mass_evol |
ADM mass at each time step, since the creation of the slice. | |
Metric * | p_gamma |
Pointer on the induced metric at the current time step (jtime ) | |
Friends | |
ostream & | operator<< (ostream &, const Time_slice &) |
Display. | |
Lorene::Time_slice::Time_slice | ( | const Scalar & | lapse_in, |
const Vector & | shift_in, | ||
const Sym_tensor & | gamma_in, | ||
const Sym_tensor & | kk_in, | ||
int | depth_in = 3 |
||
) |
General constructor (Hamiltonian-like).
lapse_in | lapse function N |
shift_in | shift vector |
gamma_in | induced metric (covariant or contravariant components) |
kk_in | extrinsic curvature (covariant or contravariant components) |
depth_in | number of stored time slices; this parameter is used to set the scheme_order member with scheme_order = depth_in - 1. scheme_order can be changed afterwards by the method set_scheme_order(int) . |
Definition at line 112 of file time_slice.C.
References gam_dd_evol, gam_uu_evol, jtime, k_dd_evol, k_uu_evol, p_gamma, set_der_0x0(), the_time, trk_evol, and Lorene::Evolution_std< TyT >::update().
Lorene::Time_slice::Time_slice | ( | const Scalar & | lapse_in, |
const Vector & | shift_in, | ||
const Evolution_std< Sym_tensor > & | gamma_in | ||
) |
General constructor (Lagrangian-like).
lapse_in | lapse function N |
shift_in | shift vector |
gamma_in | induced metric (covariant or contravariant components) at various time steps; note that the scheme_order member is set to gamma_in.size - 1. scheme_order can be changed afterwards by the method set_scheme_order(int) . |
Definition at line 155 of file time_slice.C.
References set_der_0x0().
Constructor as standard time slice of flat spacetime (Minkowski).
mp | Mapping on which the various Lorene fields will be constructed |
triad | vector basis with respect to which the various tensor components will be defined |
depth_in | number of stored time slices; this parameter is used to set the scheme_order member with scheme_order = depth_in - 1. scheme_order can be changed afterwards by the method set_scheme_order(int) . |
Definition at line 180 of file time_slice.C.
References beta_evol, Lorene::Map::flat_met_cart(), Lorene::Map::flat_met_spher(), gam_dd_evol, jtime, k_dd_evol, n_evol, set_der_0x0(), the_time, trk_evol, and Lorene::Evolution_std< TyT >::update().
Lorene::Time_slice::Time_slice | ( | const Map & | mp, |
const Base_vect & | triad, | ||
FILE * | fich, | ||
bool | partial_read, | ||
int | depth_in = 3 |
||
) |
Constructor from binary file.
The binary file must have been created by method save
.
mp | Mapping on which the various Lorene fields will be constructed |
triad | vector basis with respect to which the various tensor components will be defined |
fich | file containing the saved Time_slice |
partial_read | indicates whether the full object must be read in file or whether the final construction is devoted to a constructor of a derived class |
depth_in | number of stored time slices; the given value must coincide with that stored in the file. |
Definition at line 235 of file time_slice.C.
References beta_evol, depth, Lorene::fread_be(), Lorene::Map::get_mg(), jtime, n_evol, scheme_order, set_der_0x0(), the_time, and Lorene::Evolution_std< TyT >::update().
Lorene::Time_slice::Time_slice | ( | const Time_slice & | tin | ) |
|
explicitprotected |
Special constructor for derived classes.
Definition at line 333 of file time_slice.C.
References set_der_0x0().
|
virtual |
|
virtual |
Returns the ADM mass (geometrical units) at the current step.
Moreover this method updates adm_mass_evol
if necessary.
Reimplemented in Lorene::Tslice_dirac_max.
Definition at line 74 of file tslice_adm_mass.C.
References adm_mass_evol, Lorene::Vector::flux(), gam_dd(), Lorene::Map::get_mg(), jtime, the_time, Lorene::Evolution_std< TyT >::update(), and Lorene::Map::val_r().
shift vector jtime
)
Definition at line 87 of file time_slice_access.C.
References beta_evol, Lorene::Evolution< TyT >::is_known(), and jtime.
Tbl Lorene::Time_slice::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.
strain_tensor | : a pointer on the strain_tensor ![]() ![]() ![]() |
energy_density | : a pointer on the energy density E (see check_hamiltonian_constraint ) |
ost | : output stream for a formatted output of the result |
Definition at line 139 of file tslice_check_einstein.C.
References beta(), Lorene::contract(), gam(), jtime, k_dd(), k_dd_evol, k_uu(), Lorene::maxabs(), nn(), Lorene::Evolution< TyT >::position(), scheme_order, Lorene::Evolution< TyT >::time_derive(), and trk().
Tbl Lorene::Time_slice::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.
energy_density | : a pointer on the energy density E measured by the Eulerian observer of 4-velocity ![]() |
ost | : output stream for a formatted output of the result |
Definition at line 79 of file tslice_check_einstein.C.
References Lorene::contract(), gam(), k_dd(), k_uu(), Lorene::maxabs(), and trk().
Tbl Lorene::Time_slice::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.
momentum_density | : a pointer on the momentum density ![]() ![]() ![]() |
ost | : output stream for a formatted output of the result |
Definition at line 109 of file tslice_check_einstein.C.
References gam(), k_uu(), Lorene::maxabs(), and trk().
|
protectedvirtual |
Deletes all the derived quantities.
Reimplemented in Lorene::Time_slice_conf.
Definition at line 367 of file time_slice.C.
References adm_mass_evol, Lorene::Evolution< TyT >::downdate(), jtime, p_gamma, and set_der_0x0().
Induced metric jtime
)
Definition at line 95 of file time_slice_access.C.
|
virtual |
Induced metric (covariant components jtime
)
Reimplemented in Lorene::Time_slice_conf.
Definition at line 107 of file time_slice_access.C.
References Lorene::Metric::cov(), gam_dd_evol, gam_uu_evol, Lorene::Evolution< TyT >::is_known(), jtime, p_gamma, the_time, and Lorene::Evolution_std< TyT >::update().
|
virtual |
Induced metric (contravariant components jtime
)
Reimplemented in Lorene::Time_slice_conf.
Definition at line 122 of file time_slice_access.C.
References gam(), gam_dd_evol, gam_uu_evol, Lorene::Evolution< TyT >::is_known(), jtime, the_time, and Lorene::Evolution_std< TyT >::update().
|
inline |
Gets the latest value of time step index.
Definition at line 343 of file time_slice.h.
References jtime.
|
inline |
Gets the order of the finite-differences scheme.
Definition at line 340 of file time_slice.h.
References scheme_order.
|
inline |
Gets the time coordinate t at successive time steps.
Definition at line 346 of file time_slice.h.
References the_time.
|
virtual |
Extrinsic curvature tensor (covariant components jtime
)
Reimplemented in Lorene::Time_slice_conf.
Definition at line 135 of file time_slice_access.C.
References beta(), gam(), gam_dd(), gam_dd_evol, Lorene::Evolution< TyT >::is_known(), jtime, k_dd_evol, nn(), scheme_order, the_time, Lorene::Evolution< TyT >::time_derive(), and Lorene::Evolution_std< TyT >::update().
|
virtual |
Extrinsic curvature tensor (contravariant components jtime
)
Reimplemented in Lorene::Time_slice_conf.
Definition at line 157 of file time_slice_access.C.
References beta(), gam(), gam_uu(), gam_uu_evol, Lorene::Evolution< TyT >::is_known(), jtime, k_uu_evol, nn(), scheme_order, the_time, Lorene::Evolution< TyT >::time_derive(), and Lorene::Evolution_std< TyT >::update().
Lapse function N at the current time step (jtime
)
Reimplemented in Lorene::Time_slice_conf.
Definition at line 79 of file time_slice_access.C.
References Lorene::Evolution< TyT >::is_known(), jtime, and n_evol.
void Lorene::Time_slice::operator= | ( | const Time_slice & | tin | ) |
Assignment to another Time_slice
.
Definition at line 388 of file time_slice.C.
References beta_evol, del_deriv(), depth, gam_dd_evol, gam_uu_evol, jtime, k_dd_evol, k_uu_evol, n_evol, scheme_order, Lorene::Evolution< TyT >::the_time, the_time, and trk_evol.
Operator >> (virtual function called by the operator<<).
Reimplemented in Lorene::Isol_hor.
Definition at line 411 of file time_slice.C.
References adm_mass(), adm_mass_evol, beta_evol, depth, gam_dd_evol, gam_uu_evol, Lorene::Evolution< TyT >::is_known(), jtime, k_dd_evol, k_uu_evol, Lorene::maxabs(), n_evol, p_gamma, scheme_order, the_time, and trk_evol.
Total or partial saves in a binary file.
This protected method is to be called either from public method save
or from method sauve
of a derived class.
fich | binary file |
partial_save | indicates whether the whole object must be saved. |
Reimplemented in Lorene::Isol_hor, and Lorene::Tslice_dirac_max.
Definition at line 507 of file time_slice.C.
References beta(), beta_evol, depth, Lorene::fwrite_be(), Lorene::Evolution< TyT >::is_known(), jtime, n_evol, nn(), scheme_order, and the_time.
Saves in a binary file.
The saved data is sufficient to restore the whole time slice via the constructor from file.
rootname | root for the file name; the current time step index will be appended to it. |
Definition at line 461 of file time_slice.C.
References beta(), depth, Lorene::fwrite_be(), Lorene::Map::get_mg(), jtime, nn(), Lorene::Base_vect::sauve(), Lorene::Map::sauve(), and sauve().
|
protected |
Sets to 0x0
all the pointers on derived quantities.
Definition at line 377 of file time_slice.C.
References p_gamma.
Sets the order of the finite-differences scheme.
Definition at line 331 of file time_slice.h.
References scheme_order.
Trace K of the extrinsic curvature at the current time step (jtime
)
Reimplemented in Lorene::Time_slice_conf, and Lorene::Tslice_dirac_max.
Definition at line 177 of file time_slice_access.C.
References gam(), Lorene::Evolution< TyT >::is_known(), jtime, k_dd(), k_uu(), k_uu_evol, the_time, trk_evol, and Lorene::Evolution_std< TyT >::update().
|
friend |
Display.
Definition at line 453 of file time_slice.C.
|
mutableprotected |
ADM mass at each time step, since the creation of the slice.
At a given time step j
, adm_mass_evol
[j] is a 1-D Tbl
of size the number nz
of domains, containing the "ADM mass" evaluated at the outer boundary of each domain. The true ADM mass is thus the last value, i.e. adm_mass_evol
[j](nz-1).
Definition at line 233 of file time_slice.h.
|
mutableprotected |
Values at successive time steps of the shift vector
Definition at line 219 of file time_slice.h.
|
protected |
Number of stored time slices.
Definition at line 179 of file time_slice.h.
|
mutableprotected |
Values at successive time steps of the covariant components of the induced metric
Definition at line 198 of file time_slice.h.
|
mutableprotected |
Values at successive time steps of the contravariant components of the induced metric
Definition at line 203 of file time_slice.h.
|
protected |
Time step index of the latest slice.
Definition at line 190 of file time_slice.h.
|
mutableprotected |
Values at successive time steps of the covariant components of the extrinsic curvature tensor
Definition at line 208 of file time_slice.h.
|
mutableprotected |
Values at successive time steps of the contravariant components of the extrinsic curvature tensor
Definition at line 213 of file time_slice.h.
|
mutableprotected |
Values at successive time steps of the lapse function N.
Definition at line 216 of file time_slice.h.
|
mutableprotected |
Pointer on the induced metric at the current time step (jtime
)
Definition at line 239 of file time_slice.h.
|
protected |
Order of the finite-differences scheme for the computation of time derivatives.
This order is not constant and can be adjusted via set_scheme_order()
.
Definition at line 187 of file time_slice.h.
|
protected |
Time label of each slice.
Definition at line 193 of file time_slice.h.
|
mutableprotected |
Values at successive time steps of the trace K of the extrinsic curvature.
Definition at line 224 of file time_slice.h.