LORENE
time_slice.h
1 /*
2  * Definition of Lorene class Time_slice
3  *
4  */
5 
6 /*
7  * Copyright (c) 2004 Eric Gourgoulhon, Jose Luis Jaramillo & Jerome Novak
8  *
9  * This file is part of LORENE.
10  *
11  * LORENE is free software; you can redistribute it and/or modify
12  * it under the terms of the GNU General Public License version 2
13  * as published by the Free Software Foundation.
14  *
15  * LORENE is distributed in the hope that it will be useful,
16  * but WITHOUT ANY WARRANTY; without even the implied warranty of
17  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18  * GNU General Public License for more details.
19  *
20  * You should have received a copy of the GNU General Public License
21  * along with LORENE; if not, write to the Free Software
22  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
23  *
24  */
25 
26 #ifndef __TIME_SLICE_H_
27 #define __TIME_SLICE_H_
28 
29 /*
30  * $Id: time_slice.h,v 1.32 2014/10/13 08:52:37 j_novak Exp $
31  * $Log: time_slice.h,v $
32  * Revision 1.32 2014/10/13 08:52:37 j_novak
33  * Lorene classes and functions now belong to the namespace Lorene.
34  *
35  * Revision 1.31 2012/02/06 12:59:07 j_novak
36  * Correction of some errors.
37  *
38  * Revision 1.30 2010/10/20 07:58:09 j_novak
39  * Better implementation of the explicit time-integration. Not fully-tested yet.
40  *
41  * Revision 1.29 2008/12/04 18:22:49 j_novak
42  * Enhancement of the dzpuis treatment + various bug fixes.
43  *
44  * Revision 1.28 2008/12/02 15:02:21 j_novak
45  * Implementation of the new constrained formalism, following Cordero et al. 2009
46  * paper. The evolution eqs. are solved as a first-order system. Not tested yet!
47  *
48  * Revision 1.27 2008/08/19 06:41:59 j_novak
49  * Minor modifications to avoid warnings with gcc 4.3. Most of them concern
50  * cast-type operations, and constant strings that must be defined as const char*
51  *
52  * Revision 1.26 2007/11/06 14:47:06 j_novak
53  * New constructor from a rotating star in Dirac gauge (class Star_rot_Dirac).
54  * Evolution can take into account matter terms.
55  *
56  * Revision 1.25 2007/04/25 15:20:59 j_novak
57  * Corrected an error in the initialization of tildeB in
58  * Tslice_dirac_max::initial_dat_cts. + New method for solve_hij_AB.
59  *
60  * Revision 1.24 2007/03/21 14:51:48 j_novak
61  * Introduction of potentials A and tilde(B) of h^{ij} into Tslice_dirac_max.
62  *
63  * Revision 1.23 2005/03/28 19:44:00 f_limousin
64  * Function tgam() is now virtual.
65  *
66  * Revision 1.22 2004/06/24 20:36:07 e_gourgoulhon
67  * Class Time_slice_conf: added method check_psi_dot.
68  *
69  * Revision 1.21 2004/06/14 20:46:35 e_gourgoulhon
70  * Added argument method_poisson to Tslice_dirac_max::solve_hij.
71  *
72  * Revision 1.20 2004/05/31 20:28:20 e_gourgoulhon
73  * -- Class Time_slice : added inline functions get_latest_j() and
74  * get_time()
75  * -- Class Tslice_dirac_max: method hh_det_one takes now a time step
76  * as argument, to compute h^{ij} from khi and mu at an arbitrary
77  * time step and not only the latest one.
78  *
79  * Revision 1.19 2004/05/27 15:22:28 e_gourgoulhon
80  * Added functions save and sauve, as well as constructors from file.
81  *
82  * Revision 1.18 2004/05/20 20:30:37 e_gourgoulhon
83  * Added arguments check_mod and save_mod to method Tsclice_dirac_max::evolve.
84  *
85  * Revision 1.17 2004/05/17 19:52:16 e_gourgoulhon
86  * -- Method initial_data_cts: added arguments graph_device and
87  * method_poisson_vect.
88  * -- Method Tslice_dirac_max::solve_beta : added argument method
89  * -- Method Tslice_dirac_max::solve_hij : added argument graph_device
90  * -- Method Tslice_dirac_max::evolve : added arguments
91  * method_poisson_vect, nopause and graph_device.
92  *
93  * Revision 1.16 2004/05/12 15:16:25 e_gourgoulhon
94  * Added #include "metric.h" before #include "evolution.h".
95  *
96  * Revision 1.15 2004/05/09 20:56:09 e_gourgoulhon
97  * Added member adm_mass_evol and corresponding virtual method adm_mass().
98  *
99  * Revision 1.14 2004/05/06 15:23:10 e_gourgoulhon
100  * initial_data_cts is know a virtual function of class Time_slice_conf
101  * and is implemented also for class Tslice_dirac_max.
102  *
103  * Revision 1.13 2004/05/03 14:46:11 e_gourgoulhon
104  * Class Tslice_dirac_max: -- changed prototype of method solve_hij
105  * -- added new method evolve
106  *
107  * Revision 1.12 2004/04/30 14:36:15 j_novak
108  * Added the method Tslice_dirac_max::solve_hij(...)
109  * NOT READY YET!!!
110  *
111  * Revision 1.11 2004/04/30 10:51:38 e_gourgoulhon
112  * Class Tslice_dirac_max: added methods solve_n, solve_q and solve_beta
113  * for resolution of the elliptic part of Einstein equations.
114  *
115  * Revision 1.10 2004/04/29 17:07:27 e_gourgoulhon
116  * Added argument pdt to Time_slice_conf::initial_data_cts.
117  *
118  * Revision 1.9 2004/04/08 16:42:11 e_gourgoulhon
119  * Many changes:
120  * -- class Time_slice_conf: added methods set_*, changed argument list of
121  * method initial_data_cts.
122  * -- class Tslice_dirac_max: added methods set_* and hh_det_one().
123  *
124  * Revision 1.8 2004/04/05 21:21:51 e_gourgoulhon
125  * class Time_slice_conf: added method initial_data_cts (computation of
126  * initial data from conformally thin sandwich method).
127  * classes Time_slice_conf and Tslice_dirac_max: added constructor as
128  * standard time slice of Minkowski spacetime.
129  *
130  * Revision 1.7 2004/04/01 16:09:01 j_novak
131  * Trace of K_ij is now member of Time_slice (it was member of Time_slice_conf).
132  * Added new methods for checking 3+1 Einstein equations (preliminary).
133  *
134  * Revision 1.6 2004/03/30 14:00:30 j_novak
135  * New class Tslide_dirac_max (first version).
136  *
137  * Revision 1.5 2004/03/29 11:58:53 e_gourgoulhon
138  * Many modif. to class Time_slice_conf.
139  * Minor modif. to class Time_slice.
140  *
141  * Revision 1.4 2004/03/28 21:33:14 e_gourgoulhon
142  * Constructor Time_slice::Time_slice(int depth_in) declared "explicit".
143  *
144  * Revision 1.3 2004/03/28 21:27:57 e_gourgoulhon
145  * Class Time_slice: - renamed the Evolution_std with suffix "_evol".
146  * - added protected constructor for derived classes
147  * Added class Time_slice_conf.
148  *
149  * Revision 1.2 2004/03/26 13:33:02 j_novak
150  * New methods for accessing/updating members (nn(), beta(), gam_uu(), k_uu(), ...)
151  *
152  * Revision 1.1 2004/03/24 14:56:18 e_gourgoulhon
153  * First version
154  *
155  *
156  * $Header: /cvsroot/Lorene/C++/Include/time_slice.h,v 1.32 2014/10/13 08:52:37 j_novak Exp $
157  *
158  */
159 
160 #include "star_rot_dirac.h"
161 #include "evolution.h"
162 
163  //---------------------------//
164  // class Time_slice //
165  //---------------------------//
166 
167 namespace Lorene {
173 class Time_slice {
174 
175  // Data :
176  // -----
177  protected:
179  int depth ;
180 
188 
190  int jtime ;
191 
194 
199 
204 
209 
214 
217 
220 
225 
234 
235  // Derived data :
236  // ------------
237  protected:
239  mutable Metric* p_gamma ;
240 
241  // Constructors - Destructor
242  // -------------------------
243  public:
244 
256  Time_slice(const Scalar& lapse_in, const Vector& shift_in,
257  const Sym_tensor& gamma_in, const Sym_tensor& kk_in,
258  int depth_in = 3) ;
259 
269  Time_slice(const Scalar& lapse_in, const Vector& shift_in,
270  const Evolution_std<Sym_tensor>& gamma_in) ;
271 
282  Time_slice(const Map& mp, const Base_vect& triad, int depth_in = 3) ;
283 
298  Time_slice(const Map& mp, const Base_vect& triad, FILE* fich,
299  bool partial_read, int depth_in = 3) ;
300 
301  Time_slice(const Time_slice& ) ;
302 
303  protected:
307  explicit Time_slice(int depth_in) ;
308 
309  public:
310  virtual ~Time_slice() ;
311 
312 
313  // Memory management
314  // -----------------
315  protected:
316 
318  virtual void del_deriv() const ;
319 
321  void set_der_0x0() const ;
322 
323 
324  // Mutators / assignment
325  // ---------------------
326  public:
328  void operator=(const Time_slice&) ;
329 
331  void set_scheme_order(int ord) {
332  assert ((0<= ord)&&(ord < 4)) ;
333  scheme_order = ord ; } ;
334 
335  // Accessors
336  // ---------
337  public:
338 
340  int get_scheme_order() const { return scheme_order ; } ;
341 
343  int get_latest_j() const {return jtime; } ;
344 
346  const Evolution_std<double>& get_time() const {return the_time; } ;
347 
349  virtual const Scalar& nn() const ;
350 
352  virtual const Vector& beta() const ;
353 
355  const Metric& gam() const ;
356 
360  virtual const Sym_tensor& gam_dd() const ;
361 
365  virtual const Sym_tensor& gam_uu() const ;
366 
370  virtual const Sym_tensor& k_dd() const ;
371 
375  virtual const Sym_tensor& k_uu() const ;
376 
380  virtual const Scalar& trk() const ;
381 
382 
383  // Computational functions
384  // -----------------------
385  public:
401  Tbl check_hamiltonian_constraint(const Scalar* energy_density = 0x0,
402  ostream& ost = cout, bool verb=true) const ;
403 
419  Tbl check_momentum_constraint(const Vector* momentum_density = 0x0,
420  ostream& ost = cout, bool verb=true) const ;
421 
444  Tbl check_dynamical_equations(const Sym_tensor* strain_tensor = 0x0,
445  const Scalar* energy_density = 0x0,
446  ostream& ost = cout, bool verb=true) const ;
447 
452  virtual double adm_mass() const ;
453 
454  // Outputs
455  // -------
456  protected:
458  virtual ostream& operator>>(ostream& ) const ;
459 
461  friend ostream& operator<<(ostream& , const Time_slice& ) ;
462 
463  public:
470  void save(const char* rootname) const ;
471 
472  protected:
481  virtual void sauve(FILE* fich, bool partial_save) const ;
482 
483 };
484 
485 ostream& operator<<(ostream& , const Time_slice& ) ;
486 
487 
488 
489  //---------------------------//
490  // class Time_slice_conf //
491  //---------------------------//
492 
498 class Time_slice_conf : public Time_slice {
499 
500  // Data :
501  // -----
502  protected:
503 
507  const Metric_flat& ff ;
508 
518 
523 
524 
531 
543 
548 
553 
554  // Derived data :
555  // ------------
556  protected:
560  mutable Metric* p_tgamma ;
561 
563  mutable Scalar* p_psi4 ;
564 
566  mutable Scalar* p_ln_psi ;
567 
571  mutable Vector* p_hdirac ;
572 
577  mutable Vector* p_vec_X ;
578 
579  // Constructors - Destructor
580  // -------------------------
581  public:
582 
607  Time_slice_conf(const Scalar& lapse_in, const Vector& shift_in,
608  const Metric_flat& ff_in, const Scalar& psi_in,
609  const Sym_tensor& hh_in, const Sym_tensor& hata_in,
610  const Scalar& trk_in, int depth_in = 3) ;
611 
612 
627  Time_slice_conf(const Scalar& lapse_in, const Vector& shift_in,
628  const Sym_tensor& gamma_in, const Sym_tensor& kk_in,
629  const Metric_flat& ff_in, int depth_in = 3) ;
630 
643  Time_slice_conf(const Map& mp, const Base_vect& triad,
644  const Metric_flat& ff_in, int depth_in = 3) ;
645 
662  Time_slice_conf(const Map& mp, const Base_vect& triad,
663  const Metric_flat& ff_in, FILE* fich,
664  bool partial_read, int depth_in = 3) ;
665 
666  Time_slice_conf(const Time_slice_conf& ) ;
667 
668  virtual ~Time_slice_conf() ;
669 
670 
671  // Memory management
672  // -----------------
673  protected:
674 
676  virtual void del_deriv() const ;
677 
679  void set_der_0x0() const ;
680 
681 
682  // Mutators / assignment
683  // ---------------------
684  public:
686  void operator=(const Time_slice_conf&) ;
687 
689  void operator=(const Time_slice&) ;
690 
701  virtual void set_psi_del_npsi(const Scalar& psi_in) ;
702 
713  virtual void set_psi_del_n(const Scalar& psi_in) ;
714 
719  virtual void set_npsi_del_psi(const Scalar& npsi_in) ;
720 
725  virtual void set_npsi_del_n(const Scalar& npsi_in) ;
726 
735  virtual void set_hh(const Sym_tensor& hh_in) ;
736 
743  virtual void set_hata(const Sym_tensor& hata_in) ;
744 
748  virtual void set_hata_TT(const Sym_tensor_tt& hata_tt) ;
749 
756  virtual void set_hata_from_XAB(Param* par_bc=0x0, Param* par_mat=0x0) ;
757 
758  // Accessors
759  // ---------
760  public:
761 
762  // Virtual functions from base class Time_slice:
763  // ---------------------------------------------
764 
766  virtual const Scalar& nn() const ;
767 
771  virtual const Sym_tensor& gam_dd() const ;
772 
776  virtual const Sym_tensor& gam_uu() const ;
777 
781  virtual const Sym_tensor& k_dd() const ;
782 
786  virtual const Sym_tensor& k_uu() const ;
787 
788  // Virtual functions from this class:
789  // ----------------------------------
790 
795  virtual const Scalar& A_hata() const ;
796 
801  virtual const Scalar& B_hata() const ;
802 
810  virtual const Scalar& psi() const ;
811 
813  const Scalar& psi4() const ;
814 
816  const Scalar& ln_psi() const ;
817 
820  virtual const Scalar& npsi() const ;
821 
826  virtual const Metric& tgam() const ;
827 
834  virtual const Sym_tensor& hh(Param* = 0x0, Param* = 0x0) const ;
835 
841  virtual const Sym_tensor& hata() const ;
842 
848  virtual Sym_tensor aa() const ;
849 
853  virtual const Scalar& trk() const ;
854 
858  virtual const Vector& hdirac() const ;
859 
863  virtual const Vector& vec_X(int method_poisson=6) const ;
864 
865  // Computational methods
866  // ---------------------
867  public:
871  void compute_X_from_momentum_constraint
872  (const Vector& hat_S, const Sym_tensor_tt& hata_tt,
873  int iter_max = 200, double precis = 1.e-12,
874  double relax = 0.8, int methode_poisson = 6) ;
875 
881  virtual void set_AB_hata(const Scalar& A_in, const Scalar& B_in) ;
882 
916  virtual void initial_data_cts(const Sym_tensor& uu, const Scalar& trk_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 ) ;
921 
926  virtual double adm_mass() const ;
927 
939  void check_psi_dot(Tbl& tlnpsi_dot, Tbl& tdiff, Tbl& tdiff_rel) const ;
940 
941  // Outputs
942  // -------
943  protected:
945  virtual ostream& operator>>(ostream& ) const ;
946 
955  virtual void sauve(FILE* fich, bool partial_save) const ;
956 
957 } ;
958  //----------------------------//
959  // class Tslice_dirac_max //
960  //----------------------------//
961 
969 
970  // Data :
971  // -----
972  protected:
978 
984 
990 
996 
1002 
1008 
1010  mutable Evolution_std<Scalar> trh_evol ;
1011 
1012 
1013  // Constructors - Destructor
1014  // -------------------------
1015  public:
1040  Tslice_dirac_max(const Scalar& lapse_in, const Vector& shift_in,
1041  const Metric_flat& ff_in, const Scalar& psi_in,
1042  const Sym_tensor_trans& hh_in, const Sym_tensor& hata_in,
1043  int depth_in = 3) ;
1044 
1057  Tslice_dirac_max(const Map& mp, const Base_vect& triad,
1058  const Metric_flat& ff_in, int depth_in = 3) ;
1059 
1076  Tslice_dirac_max(const Map& mp, const Base_vect& triad,
1077  const Metric_flat& ff_in, FILE* fich,
1078  bool partial_read, int depth_in = 3) ;
1079 
1081  Tslice_dirac_max(const Star_rot_Dirac& star, double pdt, int depth_in = 3) ;
1082 
1084 
1085  virtual ~Tslice_dirac_max() ;
1086 
1087 
1088  // Mutators / assignment
1089  // ---------------------
1090  public:
1092  void operator=(const Tslice_dirac_max&) ;
1093 
1094  // Virtual functions from base class Time_slice_conf:
1095  // -------------------------------------------------
1096 
1105  virtual void set_hh(const Sym_tensor& hh_in) ;
1106 
1140  virtual void initial_data_cts(const Sym_tensor& uu, const Scalar& trk_in,
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 ) ;
1145 
1146  // Virtual functions from this class:
1147  // ----------------------------------
1148 
1156  virtual void set_khi_mu(const Scalar& khi_in, const Scalar& mu_in) ;
1157 
1164  virtual void set_AB_hh(const Scalar& A_in, const Scalar& B_in) ;
1165 
1172  virtual void set_trh(const Scalar& trh_in) ;
1173 
1184  virtual Scalar solve_psi(const Scalar* ener_dens=0x0) const ;
1185 
1201  virtual Scalar solve_npsi(const Scalar* ener_dens=0x0,
1202  const Scalar* trace_stress=0x0) const ;
1203 
1214  virtual Vector solve_beta(int method = 6) const ;
1215 
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,
1243  const Sym_tensor* strain_euler = 0x0) ;
1244 
1249  virtual double adm_mass() const ;
1250 
1251  protected:
1260  void compute_sources(const Sym_tensor* strain_tensor = 0x0) const ;
1261 
1263  void initialize_sources_copy() const ;
1264 
1272  void hh_det_one(int j, Param* par_bc = 0x0, Param* par_mat = 0x0) const ;
1273 
1279  void hh_det_one(const Sym_tensor_tt& hijtt, Param* par_mat = 0x0) const ;
1280 
1281  // Accessors
1282  // ---------
1283  public:
1284  // Virtual functions from base class Time_slice_conf:
1285  // -------------------------------------------------
1286 
1293  virtual const Sym_tensor& hh(Param* par_bc = 0x0, Param* par_mat = 0x0) const ;
1294 
1299  virtual const Scalar& trk() const ;
1300 
1305  virtual const Vector& hdirac() const ;
1306 
1307  // Virtual functions from this class:
1308  // ----------------------------------
1309 
1314  virtual const Scalar& A_hh() const ;
1315 
1320  virtual const Scalar& B_hh() const ;
1321 
1326  virtual const Scalar& trh() const ;
1327 
1328 
1329  // Outputs
1330  // -------
1331  protected:
1333  virtual ostream& operator>>(ostream& ) const ;
1334 
1343  virtual void sauve(FILE* fich, bool partial_save) const ;
1344 
1345 };
1346 }
1347 #endif
void operator=(const Time_slice &)
Assignment to another Time_slice.
Definition: time_slice.C:388
virtual const Vector & beta() const
shift vector at the current time step (jtime )
Metric for tensor calculation.
Definition: metric.h:90
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).
Definition: time_slice.C:112
Evolution_std< Scalar > npsi_evol
Values at successive time steps of the factor .
Definition: time_slice.h:522
Evolution_std< Scalar > psi_evol
Values at successive time steps of the conformal factor relating the physical metric to the conform...
Definition: time_slice.h:517
Scalar * p_ln_psi
Pointer on the logarithm of at the current time step (jtime)
Definition: time_slice.h:566
Lorene prototypes.
Definition: app_hor.h:64
Evolution_full< Tbl > adm_mass_evol
ADM mass at each time step, since the creation of the slice.
Definition: time_slice.h:233
void save(const char *rootname) const
Saves in a binary file.
Definition: time_slice.C:461
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.
Flat metric for tensor calculation.
Definition: metric.h:261
int scheme_order
Order of the finite-differences scheme for the computation of time derivatives.
Definition: time_slice.h:187
Tensor field of valence 0 (or component of a tensorial field).
Definition: scalar.h:387
Evolution_std< Scalar > A_hata_evol
Potential A associated with the symmetric tensor .
Definition: time_slice.h:547
Base class for coordinate mappings.
Definition: map.h:670
int jtime
Time step index of the latest slice.
Definition: time_slice.h:190
Evolution_std< Sym_tensor > gam_dd_evol
Values at successive time steps of the covariant components of the induced metric ...
Definition: time_slice.h:198
Tensor field of valence 1.
Definition: vector.h:188
void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: time_slice.C:377
Vectorial bases (triads) with respect to which the tensorial components are defined.
Definition: base_vect.h:105
virtual const Sym_tensor & gam_dd() const
Induced metric (covariant components ) at the current time step (jtime )
Vector * p_hdirac
Pointer on the vector (which vanishes in Dirac gauge), at the current time step (jtime).
Definition: time_slice.h:571
Evolution_std< Scalar > B_hata_evol
Potential associated with the symmetric tensor .
Definition: time_slice.h:552
Evolution_std< Sym_tensor > k_uu_evol
Values at successive time steps of the contravariant components of the extrinsic curvature tensor ...
Definition: time_slice.h:213
virtual const Scalar & trk() const
Trace K of the extrinsic curvature at the current time step (jtime )
virtual const Sym_tensor & k_uu() const
Extrinsic curvature tensor (contravariant components ) at the current time step (jtime ) ...
friend ostream & operator<<(ostream &, const Time_slice &)
Display.
Definition: time_slice.C:453
Spacelike time slice of a 3+1 spacetime.
Definition: time_slice.h:173
Evolution_std< Scalar > source_A_hata_evol
The potential A of the source of equation for .
Definition: time_slice.h:1001
Evolution_std< Scalar > source_A_hh_evol
The A potential of the source of equation for .
Definition: time_slice.h:989
virtual const Sym_tensor & gam_uu() const
Induced metric (contravariant components ) at the current time step (jtime )
Evolution_std< Sym_tensor > gam_uu_evol
Values at successive time steps of the contravariant components of the induced metric ...
Definition: time_slice.h:203
int get_scheme_order() const
Gets the order of the finite-differences scheme.
Definition: time_slice.h:340
Parameter storage.
Definition: param.h:125
Evolution_std< Sym_tensor > k_dd_evol
Values at successive time steps of the covariant components of the extrinsic curvature tensor ...
Definition: time_slice.h:208
const Metric & gam() const
Induced metric at the current time step (jtime )
Spacelike time slice of a 3+1 spacetime with conformal decomposition.
Definition: time_slice.h:498
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.
Scalar * p_psi4
Pointer on the factor at the current time step (jtime)
Definition: time_slice.h:563
Evolution_std< Scalar > source_B_hh_evol
The potential of the source of equation for .
Definition: time_slice.h:995
Metric * p_gamma
Pointer on the induced metric at the current time step (jtime)
Definition: time_slice.h:239
Class for relativistic rotating stars in Dirac gauge and maximal slicing.
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator<<).
Definition: time_slice.C:411
Transverse symmetric tensors of rank 2.
Definition: sym_tensor.h:608
virtual const Sym_tensor & k_dd() const
Extrinsic curvature tensor (covariant components ) at the current time step (jtime ) ...
Metric * p_tgamma
Pointer on the conformal metric at the current time step (jtime)
Definition: time_slice.h:560
virtual double adm_mass() const
Returns the ADM mass (geometrical units) at the current step.
Evolution_std< double > the_time
Time label of each slice.
Definition: time_slice.h:193
int get_latest_j() const
Gets the latest value of time step index.
Definition: time_slice.h:343
const Metric_flat & ff
Pointer on the flat metric with respect to which the conformal decomposition is performed.
Definition: time_slice.h:507
Vector * p_vec_X
Pointer on the vector representing the longitudinal part of .
Definition: time_slice.h:577
Evolution_std< Scalar > B_hh_evol
The potential of .
Definition: time_slice.h:983
Evolution_std< Scalar > trk_evol
Values at successive time steps of the trace K of the extrinsic curvature.
Definition: time_slice.h:224
Evolution_std< Scalar > n_evol
Values at successive time steps of the lapse function N.
Definition: time_slice.h:216
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: time_slice.C:367
Basic array class.
Definition: tbl.h:161
void set_scheme_order(int ord)
Sets the order of the finite-differences scheme.
Definition: time_slice.h:331
Evolution_std< Sym_tensor > hata_evol
Values at successive time steps of the components .
Definition: time_slice.h:542
const Evolution_std< double > & get_time() const
Gets the time coordinate t at successive time steps.
Definition: time_slice.h:346
Evolution_std< Vector > beta_evol
Values at successive time steps of the shift vector .
Definition: time_slice.h:219
int depth
Number of stored time slices.
Definition: time_slice.h:179
Spacelike time slice of a 3+1 spacetime with conformal decomposition in the maximal slicing and Dirac...
Definition: time_slice.h:968
virtual ~Time_slice()
Destructor.
Definition: time_slice.C:357
Evolution_std< Sym_tensor > hh_evol
Values at successive time steps of the components .
Definition: time_slice.h:530
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.
Class intended to describe valence-2 symmetric tensors.
Definition: sym_tensor.h:223
Evolution_std< Scalar > source_B_hata_evol
The potential of the source of equation for .
Definition: time_slice.h:1007
Transverse and traceless symmetric tensors of rank 2.
Definition: sym_tensor.h:938
virtual void sauve(FILE *fich, bool partial_save) const
Total or partial saves in a binary file.
Definition: time_slice.C:507
virtual const Scalar & nn() const
Lapse function N at the current time step (jtime )
Time evolution with full storage (*** under development ***).
Definition: evolution.h:270
Evolution_std< Scalar > A_hh_evol
The A potential of .
Definition: time_slice.h:977