28char time_slice_conf_C[] =
"$Header: /cvsroot/Lorene/C++/Source/Time_slice/time_slice_conf.C,v 1.18 2014/10/13 08:53:47 j_novak Exp $" ;
104#include "time_slice.h"
105#include "utilitaires.h"
149 "Deviation of det tgam^{ij} from 1/f")) ;
152 "Time_slice_conf::Time_slice_conf : the input h^{ij} does not"
153 <<
" ensure \n" <<
" det tgam_{ij} = f ! \n"
189 0.3333333333333333) ) ;
193 tmp.std_spectral_base() ;
224 tmp.std_spectral_base() ;
232 stmp.set_etat_zero() ;
238 tmp.set_etat_zero() ;
307 "Time_slice constructor from file: the case of full reading\n"
308 <<
" is not ready yet !" <<
endl ;
323 psi_evol(
tin.psi_evol),
324 npsi_evol(
tin.npsi_evol),
325 hh_evol(
tin.hh_evol),
326 hata_evol(
tin.hata_evol),
327 A_hata_evol(
tin.A_hata_evol),
328 B_hata_evol(
tin.B_hata_evol){
396 "Time_slice_conf::operator=(const Time_slice& ) : not implemented yet !"
503 "Deviation of det tgam^{ij} from 1/f")) ;
506 "Time_slice_conf::set_hh : the input h^{ij} does not"
507 <<
" ensure \n" <<
" det tgam_{ij} = f ! \n"
549 if (
tmp.get_dzpuis() == 3)
555 if (
tmp.get_dzpuis() == 3)
648 + 0.3333333333333333*
trk()*
gam().con(),
670 if (
tmp.get_dzpuis() == 3)
684 if (
tmp.get_dzpuis() == 3)
782 - 0.6666666666666666 *
beta().divergence(
ff) *
hh()
783 +
beta().ope_killing_conf(
ff) ;
843 for (
int i=1;
i<=3;
i++)
844 for (
int j=
i;
j<=3;
j++)
874 for (
int i=1;
i<=3;
i++)
904void Time_slice_conf::check_psi_dot(Tbl& tlnpsi_dot, Tbl& tdiff, Tbl& tdiff_rel)
const {
908 Scalar lnpsi_dot(
psi().get_mp()) ;
917 tlnpsi_dot =
max(
abs(lnpsi_dot)) ;
922 + 0.1666666666666666 * (
beta().divergence(
ff) -
nn() *
trk() ) ;
926 Tbl tref =
max(
abs(diff)) + tlnpsi_dot ;
932 tdiff_rel = tdiff / tref ;
945 flux <<
"Triad on which the components of the flat metric are defined:\n"
946 << *(
ff.get_triad()) <<
'\n' ;
962 "Conformal metric tilde gamma is up to date" << endl ;
988 for (
int j=jmin; j<=
jtime; j++) {
990 fwrite_be(&indicator,
sizeof(
int), 1, fich) ;
991 if (indicator == 1)
psi_evol[j].sauve(fich) ;
996 for (
int j=jmin; j<=
jtime; j++) {
998 fwrite_be(&indicator,
sizeof(
int), 1, fich) ;
999 if (indicator == 1)
hata_evol[j].sauve(fich) ;
1004 for (
int j=jmin; j<=
jtime; j++) {
1006 fwrite_be(&indicator,
sizeof(
int), 1, fich) ;
1011 for (
int j=jmin; j<=
jtime; j++) {
1013 fwrite_be(&indicator,
sizeof(
int), 1, fich) ;
1019 if (!partial_save) {
1020 cout <<
"Time_slice_conf::sauve: the full writing is not ready yet !"
Vectorial bases (triads) with respect to which the tensorial components are defined.
Time evolution with partial storage (*** under development ***).
virtual void update(const TyT &new_value, int j, double time_j)
Sets a new value at a given time step.
void downdate(int j)
Suppresses a stored value.
TyT time_derive(int j, int n=2) const
Computes the time derivative at time step j by means of a n-th order scheme, from the values at steps...
bool is_known(int j) const
Checks whether the value a given time step has been set.
Base class for coordinate mappings.
const Base_vect_spher & get_bvect_spher() const
Returns the orthonormal vectorial basis associated with the coordinates of the mapping.
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Flat metric for tensor calculation.
Metric for tensor calculation.
virtual const Sym_tensor & con() const
Read-only access to the contravariant representation.
virtual const Scalar & determinant() const
Returns the determinant.
Tensor field of valence 0 (or component of a tensorial field).
virtual void std_spectral_base()
Sets the spectral bases of the Valeur va to the standard ones for a scalar field.
Transverse and traceless symmetric tensors of rank 2.
Class intended to describe valence-2 symmetric tensors.
void dec_dzpuis()
dzpuis -= 1 ;
Symmetric tensors (with respect to two of their arguments).
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 ).
Time_slice_conf(const Scalar &lapse_in, const Vector &shift_in, const Metric_flat &ff_in, const Scalar &psi_in, const Sym_tensor &hh_in, const Sym_tensor &hata_in, const Scalar &trk_in, int depth_in=3)
Constructor from conformal decomposition.
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.
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.
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 .
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 void del_deriv() const
Deletes all the derived quantities.
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.
Evolution_std< Scalar > n_evol
Values at successive time steps of the lapse function N.
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 .
Evolution_std< Vector > beta_evol
Values at successive time steps of the shift vector .
Tensor field of valence 1.
Sym_tensor ope_killing_conf(const Metric &gam) const
Computes the conformal Killing operator associated with a given metric.
Tbl max(const Cmp &)
Maximum values of a Cmp in each domain.
Cmp pow(const Cmp &, int)
Power .
Cmp abs(const Cmp &)
Absolute value.
Cmp log(const Cmp &)
Neperian logarithm.
int fread_be(int *aa, int size, int nb, FILE *fich)
Reads integer(s) from a binary file according to the big endian convention.
int fwrite_be(const int *aa, int size, int nb, FILE *fich)
Writes integer(s) into a binary file according to the big endian convention.
const Map & get_mp() const
Returns the mapping.
Tenseur contract(const Tenseur &, int id1, int id2)
Self contraction of two indices of a Tenseur .
Tbl maxabs(const Tensor &aa, const char *comment=0x0, ostream &ost=cout, bool verb=true)
Maxima in each domain of the absolute values of the tensor components.