LORENE
time_slice.C
1 /*
2  * Methods of class Time_slice
3  *
4  * (see file time_slice.h for documentation).
5  *
6  */
7 
8 /*
9  * Copyright (c) 2004 Eric Gourgoulhon, Jose Luis Jaramillo & Jerome Novak
10  *
11  * This file is part of LORENE.
12  *
13  * LORENE is free software; you can redistribute it and/or modify
14  * it under the terms of the GNU General Public License version 2
15  * as published by the Free Software Foundation.
16  *
17  * LORENE is distributed in the hope that it will be useful,
18  * but WITHOUT ANY WARRANTY; without even the implied warranty of
19  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
20  * GNU General Public License for more details.
21  *
22  * You should have received a copy of the GNU General Public License
23  * along with LORENE; if not, write to the Free Software
24  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
25  *
26  */
27 
28 
29 
30 /*
31  * $Id: time_slice.C,v 1.17 2016/12/05 16:18:19 j_novak Exp $
32  * $Log: time_slice.C,v $
33  * Revision 1.17 2016/12/05 16:18:19 j_novak
34  * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
35  *
36  * Revision 1.16 2014/10/13 08:53:47 j_novak
37  * Lorene classes and functions now belong to the namespace Lorene.
38  *
39  * Revision 1.15 2014/10/06 15:13:21 j_novak
40  * Modified #include directives to use c++ syntax.
41  *
42  * Revision 1.14 2008/12/02 15:02:21 j_novak
43  * Implementation of the new constrained formalism, following Cordero et al. 2009
44  * paper. The evolution eqs. are solved as a first-order system. Not tested yet!
45  *
46  * Revision 1.13 2004/05/31 09:08:18 e_gourgoulhon
47  * Method sauve and constructor from binary file are now operational.
48  *
49  * Revision 1.12 2004/05/27 15:25:04 e_gourgoulhon
50  * Added constructors from binary file, as well as corresponding
51  * functions sauve and save.
52  *
53  * Revision 1.11 2004/05/12 15:24:20 e_gourgoulhon
54  * Reorganized the #include 's, taking into account that
55  * time_slice.h contains now an #include "metric.h".
56  *
57  * Revision 1.10 2004/05/10 09:08:34 e_gourgoulhon
58  * Added "adm_mass_evol.downdate(jtime)" in method del_deriv.
59  * Added printing of ADM mass in operator>>(ostream&).
60  *
61  * Revision 1.9 2004/05/09 20:57:34 e_gourgoulhon
62  * Added data member adm_mass_evol.
63  *
64  * Revision 1.8 2004/05/05 14:26:25 e_gourgoulhon
65  * Minor modif. in operator>>(ostream& ).
66  *
67  * Revision 1.7 2004/04/07 07:58:21 e_gourgoulhon
68  * Constructor as Minkowski slice: added call to std_spectral_base()
69  * after setting the lapse to 1.
70  *
71  * Revision 1.6 2004/04/01 16:09:02 j_novak
72  * Trace of K_ij is now member of Time_slice (it was member of Time_slice_conf).
73  * Added new methods for checking 3+1 Einstein equations (preliminary).
74  *
75  * Revision 1.5 2004/03/29 11:59:23 e_gourgoulhon
76  * Added operator>>.
77  *
78  * Revision 1.4 2004/03/28 21:29:45 e_gourgoulhon
79  * Evolution_std's renamed with suffix "_evol"
80  * Method gam() modified
81  * Added special constructor for derived classes.
82  *
83  * Revision 1.3 2004/03/26 13:33:02 j_novak
84  * New methods for accessing/updating members (nn(), beta(), gam_uu(), k_uu(), ...)
85  *
86  * Revision 1.2 2004/03/26 08:22:56 e_gourgoulhon
87  * Modifications to take into account the new setting of class
88  * Evolution.
89  *
90  * Revision 1.1 2004/03/24 14:57:17 e_gourgoulhon
91  * First version
92  *
93  *
94  * $Header: /cvsroot/Lorene/C++/Source/Time_slice/time_slice.C,v 1.17 2016/12/05 16:18:19 j_novak Exp $
95  *
96  */
97 
98 // C headers
99 #include <cassert>
100 
101 // Lorene headers
102 #include "time_slice.h"
103 #include "utilitaires.h"
104 
105 
106 
107  //--------------//
108  // Constructors //
109  //--------------//
110 
111 
112 namespace Lorene {
113 // Standard constructor (Hamiltonian-like)
114 // ---------------------------------------
115 Time_slice::Time_slice(const Scalar& lapse_in, const Vector& shift_in,
116  const Sym_tensor& gamma_in, const Sym_tensor& kk_in,
117  int depth_in)
118  : depth(depth_in),
119  scheme_order(depth_in-1),
120  jtime(0),
121  the_time(0., depth_in),
122  gam_dd_evol(depth_in),
123  gam_uu_evol(depth_in),
124  k_dd_evol(depth_in),
125  k_uu_evol(depth_in),
126  n_evol(lapse_in, depth_in),
127  beta_evol(shift_in, depth_in),
128  trk_evol(depth_in),
129  adm_mass_evol() {
130 
131  set_der_0x0() ;
132 
133  double time_init = the_time[jtime] ;
134 
135  p_gamma = new Metric(gamma_in) ;
136 
137  if (gamma_in.get_index_type(0) == COV) {
138  gam_dd_evol.update(gamma_in, jtime, time_init) ;
139  }
140  else {
141  gam_uu_evol.update(gamma_in, jtime, time_init) ;
142  }
143 
144  if (kk_in.get_index_type(0) == COV) {
145  k_dd_evol.update(kk_in, jtime, time_init) ;
146  }
147  else {
148  k_uu_evol.update(kk_in, jtime, time_init) ;
149  }
150 
151  trk_evol.update( kk_in.trace(*p_gamma), jtime, the_time[jtime] ) ;
152 
153 }
154 
155 
156 // Standard constructor (Lagrangian-like)
157 // ---------------------------------------
158 Time_slice::Time_slice(const Scalar& lapse_in, const Vector& shift_in,
159  const Evolution_std<Sym_tensor>& gamma_in)
160  : depth(gamma_in.get_size()),
161  scheme_order(gamma_in.get_size()-1),
162  jtime(0),
163  the_time(0., gamma_in.get_size()),
164  gam_dd_evol( gamma_in.get_size() ),
165  gam_uu_evol( gamma_in.get_size() ),
166  k_dd_evol( gamma_in.get_size() ),
167  k_uu_evol( gamma_in.get_size() ),
168  n_evol(lapse_in, gamma_in.get_size() ),
169  beta_evol(shift_in, gamma_in.get_size() ),
170  trk_evol(gamma_in.get_size() ),
171  adm_mass_evol() {
172 
173  cerr <<
174  "Time_slice constuctor from evolution of gamma not implemented yet !\n" ;
175  abort() ;
176 
177  set_der_0x0() ;
178 
179 }
180 
181 // Constructor as standard time slice of flat spacetime (Minkowski)
182 //-----------------------------------------------------------------
183 Time_slice::Time_slice(const Map& mp, const Base_vect& triad, int depth_in)
184  : depth(depth_in),
185  scheme_order(depth_in-1),
186  jtime(0),
187  the_time(0., depth_in),
188  gam_dd_evol(depth_in),
189  gam_uu_evol(depth_in),
190  k_dd_evol(depth_in),
191  k_uu_evol(depth_in),
192  n_evol(depth_in),
193  beta_evol(depth_in),
194  trk_evol(depth_in),
195  adm_mass_evol() {
196 
197  double time_init = the_time[jtime] ;
198 
199  const Base_vect_spher* ptriad_s =
200  dynamic_cast<const Base_vect_spher*>(&triad) ;
201  bool spher = (ptriad_s != 0x0) ;
202 
203  if (spher) {
204  gam_dd_evol.update( mp.flat_met_spher().cov(), jtime, time_init) ;
205  }
206  else {
207  assert( dynamic_cast<const Base_vect_cart*>(&triad) != 0x0) ;
208  gam_dd_evol.update( mp.flat_met_cart().cov(), jtime, time_init) ;
209  }
210 
211 
212  // K_ij identically zero:
213  Sym_tensor ktmp(mp, COV, triad) ;
214  ktmp.set_etat_zero() ;
215  k_dd_evol.update(ktmp, jtime, time_init) ;
216 
217  // Lapse identically one:
218  Scalar tmp(mp) ;
219  tmp.set_etat_one() ;
220  tmp.std_spectral_base() ;
221  n_evol.update(tmp, jtime, time_init) ;
222 
223  // shift identically zero:
224  Vector btmp(mp, CON, triad) ;
225  btmp.set_etat_zero() ;
226  beta_evol.update(btmp, jtime, time_init) ;
227 
228  // trace(K) identically zero:
229  tmp.set_etat_zero() ;
230  trk_evol.update(tmp, jtime, time_init) ;
231 
232  set_der_0x0() ;
233 }
234 
235 // Constructor from binary file
236 // ----------------------------
237 
238 Time_slice::Time_slice(const Map& mp, const Base_vect& triad, FILE* fich,
239  bool partial_read, int depth_in)
240  : depth(depth_in),
241  the_time(depth_in),
242  gam_dd_evol(depth_in),
243  gam_uu_evol(depth_in),
244  k_dd_evol(depth_in),
245  k_uu_evol(depth_in),
246  n_evol(depth_in),
247  beta_evol(depth_in),
248  trk_evol(depth_in),
249  adm_mass_evol() {
250 
251  // Reading various integer parameters
252  // ----------------------------------
253 
254  int depth_file ;
255  fread_be(&depth_file, sizeof(int), 1, fich) ;
256  if (depth_file != depth_in) {
257  cout <<
258  "Time_slice constructor from file: the depth read in file \n"
259  << " is different from that given in the argument list : \n"
260  << " depth_file = " << depth_file
261  << " <-> depth_in " << depth_in << " !" << endl ;
262  abort() ;
263  }
264  fread_be(&scheme_order, sizeof(int), 1, fich) ;
265  fread_be(&jtime, sizeof(int), 1, fich) ;
266 
267  // Reading the_time
268  // ----------------
269  int jmin = jtime - depth + 1 ;
270  int indicator ;
271  for (int j=jmin; j<=jtime; j++) {
272  fread_be(&indicator, sizeof(int), 1, fich) ;
273  if (indicator == 1) {
274  double xx ;
275  fread_be(&xx, sizeof(double), 1, fich) ;
276  the_time.update(xx, j, xx) ;
277  }
278  }
279 
280  // Reading of various fields
281  // -------------------------
282 
283  // N
284  for (int j=jmin; j<=jtime; j++) {
285  fread_be(&indicator, sizeof(int), 1, fich) ;
286  if (indicator == 1) {
287  Scalar nn_file(mp, *(mp.get_mg()), fich) ;
288  n_evol.update(nn_file, j, the_time[j]) ;
289  }
290  }
291 
292  // beta
293  for (int j=jmin; j<=jtime; j++) {
294  fread_be(&indicator, sizeof(int), 1, fich) ;
295  if (indicator == 1) {
296  Vector beta_file(mp, triad, fich) ;
297  beta_evol.update(beta_file, j, the_time[j]) ;
298  }
299  }
300 
301  // Case of a full reading
302  // -----------------------
303  if (!partial_read) {
304  cout <<
305  "Time_slice constructor from file: the case of full reading\n"
306  << " is not ready yet !" << endl ;
307  abort() ;
308  }
309 
310  set_der_0x0() ;
311 
312 }
313 
314 
315 // Copy constructor
316 // ----------------
318  : depth(tin.depth),
319  scheme_order(tin.scheme_order),
320  jtime(tin.jtime),
321  the_time(tin.the_time),
322  gam_dd_evol(tin.gam_dd_evol),
323  gam_uu_evol(tin.gam_uu_evol),
324  k_dd_evol(tin.k_dd_evol),
325  k_uu_evol(tin.k_uu_evol),
326  n_evol(tin.n_evol),
327  beta_evol(tin.beta_evol),
328  trk_evol(tin.trk_evol),
329  adm_mass_evol(tin.adm_mass_evol) {
330 
331  set_der_0x0() ;
332 }
333 
334 // Special constructor for derived classes
335 //----------------------------------------
337  : depth(depth_in),
338  scheme_order(depth_in-1),
339  jtime(0),
340  the_time(0., depth_in),
341  gam_dd_evol(depth_in),
342  gam_uu_evol(depth_in),
343  k_dd_evol(depth_in),
344  k_uu_evol(depth_in),
345  n_evol(depth_in),
346  beta_evol(depth_in),
347  trk_evol(depth_in),
348  adm_mass_evol() {
349 
350  set_der_0x0() ;
351 }
352 
353 
354 
355 
356  //--------------//
357  // Destructor //
358  //--------------//
359 
361 
363 
364 }
365 
366  //---------------------//
367  // Memory management //
368  //---------------------//
369 
370 void Time_slice::del_deriv() const {
371 
372  if (p_gamma != 0x0) delete p_gamma ;
373 
374  set_der_0x0() ;
375 
376  adm_mass_evol.downdate(jtime) ;
377 }
378 
379 
381 
382  p_gamma = 0x0 ;
383 
384 }
385 
386 
387  //-----------------------//
388  // Mutators / assignment //
389  //-----------------------//
390 
392 
393  depth = tin.depth;
394  scheme_order = tin.scheme_order ;
395  jtime = tin.jtime;
396  the_time = tin.the_time;
397  gam_dd_evol = tin.gam_dd_evol;
398  gam_uu_evol = tin.gam_uu_evol;
399  k_dd_evol = tin.k_dd_evol;
400  k_uu_evol = tin.k_uu_evol;
401  n_evol = tin.n_evol;
402  beta_evol = tin.beta_evol;
403  trk_evol = tin.trk_evol ;
404 
405  del_deriv() ;
406 
407 }
408 
409 
410  //------------------//
411  // output //
412  //------------------//
413 
414 ostream& Time_slice::operator>>(ostream& flux) const {
415 
416  flux << "\n------------------------------------------------------------\n"
417  << "Lorene class : " << typeid(*this).name() << '\n' ;
418  flux << "Number of stored slices : " << depth
419  << " order of time scheme : " << scheme_order << '\n'
420  << "Time label t = " << the_time[jtime]
421  << " index of time step j = " << jtime << '\n' << '\n' ;
422  if (adm_mass_evol.is_known(jtime)) {
423  flux << "ADM mass : " << adm_mass() << endl ;
424  }
425 
426  flux << "Max. of absolute values of the various fields in each domain: \n" ;
427  if (gam_dd_evol.is_known(jtime)) {
428  maxabs( gam_dd_evol[jtime], "gam_{ij}", flux) ;
429  }
430  if (gam_uu_evol.is_known(jtime)) {
431  maxabs( gam_uu_evol[jtime], "gam^{ij}", flux) ;
432  }
433  if (k_dd_evol.is_known(jtime)) {
434  maxabs( k_dd_evol[jtime], "K_{ij}", flux) ;
435  }
436  if (k_uu_evol.is_known(jtime)) {
437  maxabs( k_uu_evol[jtime], "K^{ij}", flux) ;
438  }
439  if (n_evol.is_known(jtime)) {
440  maxabs( n_evol[jtime], "N", flux) ;
441  }
442  if (beta_evol.is_known(jtime)) {
443  maxabs( beta_evol[jtime], "beta^i", flux) ;
444  }
445  if (trk_evol.is_known(jtime)) {
446  maxabs( trk_evol[jtime], "tr K", flux) ;
447  }
448 
449  if (p_gamma != 0x0) flux << "Metric gamma is up to date" << endl ;
450 
451  return flux ;
452 
453 }
454 
455 
456 ostream& operator<<(ostream& flux, const Time_slice& sigma) {
457 
458  sigma >> flux ;
459  return flux ;
460 
461 }
462 
463 
464 void Time_slice::save(const char* rootname) const {
465 
466  // Opening of file
467  // ---------------
468  char* filename = new char[ strlen(rootname)+10 ] ;
469  strcpy(filename, rootname) ;
470  char nomj[7] ;
471  sprintf(nomj, "%06d", jtime) ;
472  strcat(filename, nomj) ;
473  strcat(filename, ".d") ;
474 
475  FILE* fich = fopen(filename, "w") ;
476  if (fich == 0x0) {
477  cout << "Problem in opening file " << filename << " ! " << endl ;
478  perror(" reason") ;
479  abort() ;
480  }
481 
482  // Write grid, mapping, triad and depth
483  // ------------------------------------
484  const Map& map = nn().get_mp() ;
485  const Mg3d& mgrid = *(map.get_mg()) ;
486  const Base_vect& triad = *(beta().get_triad()) ;
487 
488  mgrid.sauve(fich) ;
489  map.sauve(fich) ;
490  triad.sauve(fich) ;
491 
492  fwrite_be(&depth, sizeof(int), 1, fich) ;
493 
494  // Write all binary data by means of virtual function sauve
495  // --------------------------------------------------------
496  bool partial_save = false ;
497  sauve(fich, partial_save) ;
498 
499  // Close the file
500  // --------------
501 
502  fclose(fich) ;
503 
504  delete [] filename ;
505 
506 }
507 
508 
509 
510 void Time_slice::sauve(FILE* fich, bool partial_save) const {
511 
512  // Writing various integer parameters
513  // ----------------------------------
514 
515  fwrite_be(&depth, sizeof(int), 1, fich) ;
516  fwrite_be(&scheme_order, sizeof(int), 1, fich) ;
517  fwrite_be(&jtime, sizeof(int), 1, fich) ;
518 
519  // Writing the_time
520  // ----------------
521  int jmin = jtime - depth + 1 ;
522  for (int j=jmin; j<=jtime; j++) {
523  int indicator = (the_time.is_known(j)) ? 1 : 0 ;
524  fwrite_be(&indicator, sizeof(int), 1, fich) ;
525  if (indicator == 1) {
526  double xx = the_time[j] ;
527  fwrite_be(&xx, sizeof(double), 1, fich) ;
528  }
529  }
530 
531  // Writing of various fields
532  // -------------------------
533 
534  // N
535  nn() ; // forces the update at the current time step
536  for (int j=jmin; j<=jtime; j++) {
537  int indicator = (n_evol.is_known(j)) ? 1 : 0 ;
538  fwrite_be(&indicator, sizeof(int), 1, fich) ;
539  if (indicator == 1) n_evol[j].sauve(fich) ;
540  }
541 
542  // beta
543  beta() ; // forces the update at the current time step
544  for (int j=jmin; j<=jtime; j++) {
545  int indicator = (beta_evol.is_known(j)) ? 1 : 0 ;
546  fwrite_be(&indicator, sizeof(int), 1, fich) ;
547  if (indicator == 1) beta_evol[j].sauve(fich) ;
548  }
549 
550  // Case of a complete save
551  // -----------------------
552  if (!partial_save) {
553 
554  cout << "Time_slice::sauve: the full writing is not ready yet !"
555  << endl ;
556  abort() ;
557  }
558 
559 
560 }
561 
562 
563 
564 
565 
566 
567 
568 
569 
570 
571 
572 
573 
574 
575 }
void operator=(const Time_slice &)
Assignment to another Time_slice.
Definition: time_slice.C:391
virtual const Vector & beta() const
shift vector at the current time step (jtime )
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.
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:115
virtual void sauve(FILE *) const
Save in a file.
Definition: map.C:227
virtual void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition: scalar.C:330
Lorene prototypes.
Definition: app_hor.h:67
Evolution_full< Tbl > adm_mass_evol
ADM mass at each time step, since the creation of the slice.
Definition: time_slice.h:236
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:783
void save(const char *rootname) const
Saves in a binary file.
Definition: time_slice.C:464
int scheme_order
Order of the finite-differences scheme for the computation of time derivatives.
Definition: time_slice.h:190
Tensor field of valence 0 (or component of a tensorial field).
Definition: scalar.h:393
virtual void sauve(FILE *) const
Save in a file.
Definition: base_vect.C:153
Base class for coordinate mappings.
Definition: map.h:688
virtual void update(const TyT &new_value, int j, double time_j)
Sets a new value at a given time step.
int jtime
Time step index of the latest slice.
Definition: time_slice.h:193
virtual void std_spectral_base()
Sets the spectral bases of the Valeur va to the standard ones for a scalar field. ...
Definition: scalar.C:790
Evolution_std< Sym_tensor > gam_dd_evol
Values at successive time steps of the covariant components of the induced metric ...
Definition: time_slice.h:201
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:380
Vectorial bases (triads) with respect to which the tensorial components are defined.
Definition: base_vect.h:105
const Metric_flat & flat_met_cart() const
Returns the flat metric associated with the Cartesian coordinates and with components expressed in th...
Definition: map.C:334
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:216
Spacelike time slice of a 3+1 spacetime.
Definition: time_slice.h:176
const Base_vect * get_triad() const
Returns the vectorial basis (triad) on which the components are defined.
Definition: tensor.h:879
void sauve(FILE *fd, bool save_base=false) const
Saves into a file.
Definition: mg3d.C:500
int get_index_type(int i) const
Gives the type (covariant or contravariant) of the index number i .
Definition: tensor.h:899
Evolution_std< Sym_tensor > gam_uu_evol
Values at successive time steps of the contravariant components of the induced metric ...
Definition: time_slice.h:206
bool is_known(int j) const
Checks whether the value a given time step has been set.
Definition: evolution.C:335
void set_etat_one()
Sets the logical state to ETATUN (one).
Definition: scalar.C:340
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:211
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.
Definition: fwrite_be.C:73
Metric * p_gamma
Pointer on the induced metric at the current time step (jtime)
Definition: time_slice.h:242
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator<<).
Definition: time_slice.C:414
Time evolution with partial storage (*** under development ***).
Definition: evolution.h:371
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:196
Spherical orthonormal vectorial bases (triads).
Definition: base_vect.h:308
int fread_be(int *aa, int size, int nb, FILE *fich)
Reads integer(s) from a binary file according to the big endian convention.
Definition: fread_be.C:72
Multi-domain grid.
Definition: grilles.h:279
virtual const Sym_tensor & cov() const
Read-only access to the covariant representation.
Definition: metric_flat.C:137
Evolution_std< Scalar > trk_evol
Values at successive time steps of the trace K of the extrinsic curvature.
Definition: time_slice.h:227
Evolution_std< Scalar > n_evol
Values at successive time steps of the lapse function N.
Definition: time_slice.h:219
virtual void set_etat_zero()
Sets the logical state of all components to ETATZERO (zero state).
Definition: tensor.C:506
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: time_slice.C:370
Evolution_std< Vector > beta_evol
Values at successive time steps of the shift vector .
Definition: time_slice.h:222
const Map & get_mp() const
Returns the mapping.
Definition: tensor.h:874
int depth
Number of stored time slices.
Definition: time_slice.h:182
virtual ~Time_slice()
Destructor.
Definition: time_slice.C:360
const Metric_flat & flat_met_spher() const
Returns the flat metric associated with the spherical coordinates and with components expressed in th...
Definition: map.C:324
Class intended to describe valence-2 symmetric tensors.
Definition: sym_tensor.h:226
Tensor trace(int ind1, int ind2) const
Trace on two different type indices.
virtual void sauve(FILE *fich, bool partial_save) const
Total or partial saves in a binary file.
Definition: time_slice.C:510
virtual const Scalar & nn() const
Lapse function N at the current time step (jtime )