LORENE
tslice_dirac_max.C
1  /*
2  * Methods of class Tslice_dirac_max
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: tslice_dirac_max.C,v 1.27 2016/12/05 16:18:19 j_novak Exp $
32  * $Log: tslice_dirac_max.C,v $
33  * Revision 1.27 2016/12/05 16:18:19 j_novak
34  * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
35  *
36  * Revision 1.26 2014/10/13 08:53:48 j_novak
37  * Lorene classes and functions now belong to the namespace Lorene.
38  *
39  * Revision 1.25 2014/10/06 15:13:22 j_novak
40  * Modified #include directives to use c++ syntax.
41  *
42  * Revision 1.24 2008/12/04 18:22:49 j_novak
43  * Enhancement of the dzpuis treatment + various bug fixes.
44  *
45  * Revision 1.23 2008/12/02 15:02:22 j_novak
46  * Implementation of the new constrained formalism, following Cordero et al. 2009
47  * paper. The evolution eqs. are solved as a first-order system. Not tested yet!
48  *
49  * Revision 1.22 2007/11/06 14:47:07 j_novak
50  * New constructor from a rotating star in Dirac gauge (class Star_rot_Dirac).
51  * Evolution can take into account matter terms.
52  *
53  * Revision 1.21 2007/09/25 16:54:11 j_novak
54  * *** empty log message ***
55  *
56  * Revision 1.20 2007/09/25 16:52:15 j_novak
57  * *** empty log message ***
58  *
59  * Revision 1.19 2007/06/05 07:38:37 j_novak
60  * Better treatment of dzpuis for A and tilde(B) potentials. Some errors in the bases manipulation have been also corrected.
61  *
62  * Revision 1.18 2007/04/25 15:21:01 j_novak
63  * Corrected an error in the initialization of tildeB in
64  * Tslice_dirac_max::initial_dat_cts. + New method for solve_hij_AB.
65  *
66  * Revision 1.17 2007/03/21 14:51:50 j_novak
67  * Introduction of potentials A and tilde(B) of h^{ij} into Tslice_dirac_max.
68  *
69  * Revision 1.16 2004/12/28 14:21:48 j_novak
70  * Added the method Sym_tensor_trans::trace_from_det_one
71  *
72  * Revision 1.15 2004/07/08 12:29:01 j_novak
73  * use of new method Tensor::annule_extern_cn
74  *
75  * Revision 1.14 2004/06/30 08:02:40 j_novak
76  * Added filtering in l of khi_new and mu_new. ki_source is forced to go to
77  * zero at least as r^2.
78  *
79  * Revision 1.13 2004/06/17 06:59:41 e_gourgoulhon
80  * -- Method initial_data_cts: re-organized treatment of vanishing uu.
81  * -- Method hh_det_one: replaced the attenuation with tempo by a call
82  * to the new method Tensor::annule_extern_c2.
83  *
84  * Revision 1.12 2004/06/08 14:05:06 j_novak
85  * Added the attenuation of khi and mu in the last domain in ::det_one(). They are set to zero in the CED.
86  *
87  * Revision 1.11 2004/05/31 20:31:31 e_gourgoulhon
88  * -- Method hh_det_one takes now a time step as argument, to compute
89  * h^{ij} from khi and mu at some arbitrary time step and not only at
90  * the latest one.
91  * -- h^{ij} is no longer saved in binary files (method sauve);
92  * accordingly, the constructor from file calls the new version of
93  * hh_det_one to restore h^{ij}.
94  *
95  * Revision 1.10 2004/05/31 09:08:18 e_gourgoulhon
96  * Method sauve and constructor from binary file are now operational.
97  *
98  * Revision 1.9 2004/05/27 15:25:04 e_gourgoulhon
99  * Added constructors from binary file, as well as corresponding
100  * functions sauve and save.
101  *
102  * Revision 1.8 2004/05/17 19:54:10 e_gourgoulhon
103  * Method initial_data_cts: added arguments graph_device and method_poisson_vect.
104  *
105  * Revision 1.7 2004/05/12 15:24:20 e_gourgoulhon
106  * Reorganized the #include 's, taking into account that
107  * time_slice.h contains now an #include "metric.h".
108  *
109  * Revision 1.6 2004/05/10 09:16:32 e_gourgoulhon
110  * -- Method initial_data_cts: added a call to del_deriv() at the end.
111  * -- Methods set_trh and hh_det_one: added "adm_mass_evol.downdate(jtime)".
112  * -- Method trh() : the update is now performed via a call to hh_det_one().
113  *
114  * Revision 1.5 2004/05/06 15:23:55 e_gourgoulhon
115  * Added method initial_data_cts.
116  *
117  * Revision 1.4 2004/05/03 08:15:48 e_gourgoulhon
118  * Method hh_det_one(): added check at the end (deviation from det = 1).
119  *
120  * Revision 1.3 2004/04/08 16:44:19 e_gourgoulhon
121  * Added methods set_* and hh_det_one().
122  *
123  * Revision 1.2 2004/04/05 21:22:49 e_gourgoulhon
124  * Added constructor as standard time slice of Minkowski spacetime.
125  *
126  * Revision 1.1 2004/03/30 14:00:31 j_novak
127  * New class Tslide_dirac_max (first version).
128  *
129  *
130  * $Header: /cvsroot/Lorene/C++/Source/Time_slice/tslice_dirac_max.C,v 1.27 2016/12/05 16:18:19 j_novak Exp $
131  *
132  */
133 
134 // C headers
135 #include <cassert>
136 
137 // Lorene headers
138 #include "time_slice.h"
139 #include "utilitaires.h"
140 
141 
142 
143  //--------------//
144  // Constructors //
145  //--------------//
146 
147 
148 // Constructor from conformal decomposition
149 // ----------------------------------------
150 
151 namespace Lorene {
152 Tslice_dirac_max::Tslice_dirac_max(const Scalar& lapse_in, const Vector& shift_in,
153  const Metric_flat& ff_in, const Scalar& psi_in,
154  const Sym_tensor_trans& hh_in, const Sym_tensor& hata_in,
155  int depth_in)
156  : Time_slice_conf( lapse_in, shift_in, ff_in, psi_in, hh_in, hata_in,
157  0*lapse_in, depth_in),
158  A_hh_evol(depth_in), B_hh_evol(depth_in), source_A_hh_evol(depth_in),
159  source_B_hh_evol(depth_in), source_A_hata_evol(depth_in) ,
160  source_B_hata_evol(depth_in), trh_evol(hh_in.the_trace(), depth_in)
161 { }
162 
163 // Constructor as standard time slice of flat spacetime (Minkowski)
164 // ----------------------------------------------------------------
165 
167  const Metric_flat& ff_in, int depth_in)
168  : Time_slice_conf(mp, triad, ff_in, depth_in),
169  A_hh_evol(depth_in), B_hh_evol(depth_in),
170  source_A_hh_evol(depth_in), source_B_hh_evol(depth_in),
171  source_A_hata_evol(depth_in), source_B_hata_evol(depth_in),
172  trh_evol(depth_in) {
173 
174  double time_init = the_time[jtime] ;
175 
176  // All potentials identically zero:
177  Scalar tmp(mp) ;
178  tmp.set_etat_zero() ;
179  A_hh_evol.update(tmp, jtime, time_init) ;
180  B_hh_evol.update(tmp, jtime, time_init) ;
181 
182  source_A_hh_evol.update(tmp, jtime, time_init) ;
183  source_B_hh_evol.update(tmp, jtime, time_init) ;
184 
185  source_A_hata_evol.update(tmp, jtime, time_init) ;
186  source_B_hata_evol.update(tmp, jtime, time_init) ;
187 
188  // tr h identically zero:
189  trh_evol.update(tmp, jtime, time_init) ;
190 }
191 
192 
193 // Constructor from binary file
194 // ----------------------------
195 
197  const Metric_flat& ff_in, FILE* fich,
198  bool partial_read, int depth_in)
199  : Time_slice_conf(mp, triad, ff_in, fich, true, depth_in),
200  A_hh_evol(depth_in), B_hh_evol(depth_in),
201  source_A_hh_evol(depth_in), source_B_hh_evol(depth_in),
202  source_A_hata_evol(depth_in), source_B_hata_evol(depth_in),
203  trh_evol(depth_in) {
204 
205  if (partial_read) {
206  cout <<
207  "Constructor of Tslice_dirac_max from file: the case of partial reading\n"
208  << " is not ready yet !"
209  << endl ;
210  abort() ;
211  }
212 
213  // Reading of various fields
214  // -------------------------
215 
216  int jmin = jtime - depth + 1 ;
217  int indicator ;
218 
219  // h^{ij}
220  for (int j=jmin; j<=jtime; j++) {
221  fread_be(&indicator, sizeof(int), 1, fich) ;
222  if (indicator == 1) {
223  Sym_tensor hh_file(mp, triad, fich) ;
224  hh_evol.update(hh_file, j, the_time[j]) ;
225  }
226  }
227 
228  // A - hh
229  for (int j=jmin; j<=jtime; j++) {
230  fread_be(&indicator, sizeof(int), 1, fich) ;
231  if (indicator == 1) {
232  Scalar A_hh_file(mp, *(mp.get_mg()), fich) ;
233  A_hh_evol.update(A_hh_file, j, the_time[j]) ;
234  }
235  }
236 
237  // B - hh
238  for (int j=jmin; j<=jtime; j++) {
239  fread_be(&indicator, sizeof(int), 1, fich) ;
240  if (indicator == 1) {
241  Scalar B_hh_file(mp, *(mp.get_mg()), fich) ;
242  B_hh_evol.update(B_hh_file, j, the_time[j]) ;
243  }
244  }
245 }
246 
247 // Constructor from a rotating star
248 // --------------------------------
249 
250 Tslice_dirac_max::Tslice_dirac_max(const Star_rot_Dirac& star, double pdt, int depth_in)
251  : Time_slice_conf(star.get_nn(), star.get_beta(), star.get_mp().flat_met_spher(),
252  0.*star.get_nn(), star.get_hh(), 0.*star.get_aa(),
253  0.*star.get_nn(), depth_in),
254  A_hh_evol(depth_in), B_hh_evol(depth_in),
255  source_A_hh_evol(depth_in), source_B_hh_evol(depth_in),
256  source_A_hata_evol(depth_in), source_B_hata_evol(depth_in),
257  trh_evol(depth_in) {
258 
259  Scalar tmp = exp(star.get_ln_psi()) ;
260  tmp.std_spectral_base() ;
261  psi_evol.downdate(jtime) ;
262  psi_evol.update(tmp, jtime, the_time[jtime]) ;
263 
264  Sym_tensor tmp2 = psi4()*psi()*psi()*star.get_aa() ;
265  hata_evol.downdate(jtime) ;
266  hata_evol.update(tmp2, jtime, the_time[jtime]) ;
267 
268  A_hh() ;
269  B_hh() ;
270  A_hata() ;
271  B_hata() ;
272  compute_sources() ;
273 
274  // Update of various fields
275  // -------------------------
276  double ttime1 = the_time[jtime] ;
277  int jtime1 = jtime ;
278  for (int j=1; j < depth; j++) {
279  jtime1++ ;
280  ttime1 += pdt ;
281  psi_evol.update(psi_evol[jtime], jtime1, ttime1) ;
282  n_evol.update(n_evol[jtime], jtime1, ttime1) ;
283  beta_evol.update(beta_evol[jtime], jtime1, ttime1) ;
284  hh_evol.update(hh_evol[jtime], jtime1, ttime1) ;
285  trk_evol.update(trk_evol[jtime], jtime1, ttime1) ;
286  A_hh_evol.update(A_hh_evol[jtime], jtime1, ttime1) ;
287  B_hh_evol.update(B_hh_evol[jtime], jtime1, ttime1) ;
288  A_hata_evol.update(A_hata_evol[jtime], jtime1, ttime1) ;
289  B_hata_evol.update(B_hata_evol[jtime], jtime1, ttime1) ;
290  trh_evol.update(trh_evol[jtime], jtime1, ttime1) ;
291  k_dd_evol.update(k_dd_evol[jtime], jtime1, ttime1) ;
292  the_time.update(ttime1, jtime1, ttime1) ;
293  }
294  jtime += depth - 1 ;
295 
297 }
298 
299 // Copy constructor
300 // ----------------
301 
303  : Time_slice_conf(tin),
304  A_hh_evol(tin.A_hh_evol),
305  B_hh_evol(tin.B_hh_evol),
306  source_A_hh_evol(tin.source_A_hh_evol),
307  source_B_hh_evol(tin.source_B_hh_evol),
308  source_A_hata_evol(tin.source_A_hata_evol),
309  source_B_hata_evol(tin.source_B_hata_evol),
310  trh_evol(tin.trh_evol) { }
311 
312 
313  //--------------//
314  // Destructor //
315  //--------------//
316 
318 
319 
320  //-----------------------//
321  // Mutators / assignment //
322  //-----------------------//
323 
325 
327 
328  A_hh_evol = tin.A_hh_evol ;
329  B_hh_evol = tin.B_hh_evol ;
334  trh_evol = tin.trh_evol ;
335 
336 }
337 
338 
340 
341  Time_slice_conf::set_hh(hh_in) ;
342 
343  // Reset of quantities depending on h^{ij}:
344  A_hh_evol.downdate(jtime) ;
345  B_hh_evol.downdate(jtime) ;
346  source_A_hh_evol.downdate(jtime) ;
347  source_B_hh_evol.downdate(jtime) ;
348  source_A_hata_evol.downdate(jtime) ;
349  source_B_hata_evol.downdate(jtime) ;
350  trh_evol.downdate(jtime) ;
351 
352 }
353 
354 
356  const Scalar& trk_in, const Scalar& trk_point,
357  double pdt, double precis, int method_poisson_vect,
358  const char* graph_device, const Scalar* p_ener_dens,
359  const Vector* p_mom_dens, const Scalar* p_trace_stress) {
360 
361 
362  Time_slice_conf::initial_data_cts(uu, trk_in, trk_point, pdt, precis,
363  method_poisson_vect, graph_device,
364  p_ener_dens, p_mom_dens, p_trace_stress) ;
365 
366  int nz = trk_in.get_mp().get_mg()->get_nzone() ;
367  // Setting khi and mu for j <= jtime, taking into account u^{ij} = dh^{ij}/dt
368  //--------------------------------------------------------------------------
369  for (int j = jtime-depth+1 ; j <= jtime; j++) {
370 
371  // A and tildeB are computed from the value of hh
372  Scalar tmp = hh_evol[j].compute_A(true) ;
373  assert (tmp.get_etat() != ETATNONDEF) ;
374  if (tmp.get_etat() != ETATZERO) {
375  assert(tmp.get_mp().get_mg()->get_type_r(nz-1) == UNSURR) ;
376  tmp.annule_domain(nz-1) ;
377  }
378  tmp.set_dzpuis(0) ;
379  A_hh_evol.update(tmp, j, the_time[j]) ;
380 
381  tmp = hh_evol[jtime].compute_tilde_B_tt(true) ;
382  assert (tmp.get_etat() != ETATNONDEF) ;
383  if (tmp.get_etat() != ETATZERO) {
384  assert(tmp.get_mp().get_mg()->get_type_r(nz-1) == UNSURR) ;
385  tmp.annule_domain(nz-1) ;
386  }
387  tmp.set_dzpuis(0) ;
388  B_hh_evol.update(tmp, j, the_time[j]) ;
389  }
390 
391  cout << endl <<
392  "Tslice_dirac_max::initial_data_cts : variation of A and tilde(B) for J = "
393  << jtime << " :\n" ;
394  maxabs(A_hh_evol[jtime] - A_hh_evol[jtime-1], "A(h)^J - A(h)^{J-1}") ;
395 
396  maxabs(B_hh_evol[jtime] - B_hh_evol[jtime-1], "B(h)^J - B(h)^{J-1}") ;
397 
398  maxabs(A_hata_evol[jtime] - A_hata_evol[jtime-1], "A(hat{A})^J - A(hat{A})^{J-1}") ;
399 
400  maxabs(B_hata_evol[jtime] - B_hata_evol[jtime-1], "B(hat{A})^J - B(hat{A})^{J-1}") ;
401 
402  // Reset of derived quantities (at the new time step jtime)
403  // ---------------------------
404  del_deriv() ;
405 
406  compute_sources() ;
407  initialize_sources_copy() ; //## should be a call to the Runge-Kutta integrator
408 }
409 
410 
411 void Tslice_dirac_max::set_khi_mu(const Scalar& khi_in, const Scalar& mu_in) {
412 
413  const Map& mp = khi_in.get_mp() ;
414 
415  Sym_tensor_tt hh_tt(mp, mp.get_bvect_spher(), mp.flat_met_spher());
416  hh_tt.set_khi_mu(khi_in, mu_in, 2) ;
417 
418  Sym_tensor_trans hh_tmp(mp, mp.get_bvect_spher(), mp.flat_met_spher());
419  hh_tmp.trace_from_det_one(hh_tt) ;
420 
421  // Result set to trh_evol and hh_evol
422  // ----------------------------------
423  Scalar tmp = hh_tmp.the_trace() ;
424  tmp.dec_dzpuis(4) ;
425  trh_evol.update(tmp, jtime, the_time[jtime]) ;
426 
427  // The longitudinal part of h^{ij}, which is zero by virtue of Dirac gauge :
428  Vector wzero(mp, CON, *(ff.get_triad())) ;
429  wzero.set_etat_zero() ;
430 
431  // Temporary Sym_tensor with longitudinal part set to zero :
432  Sym_tensor hh_new(mp, CON, *(ff.get_triad())) ;
433 
434  hh_new.set_longit_trans(wzero, hh_tmp) ;
435 
436  hh_evol.update(hh_new, jtime, the_time[jtime]) ;
437 
438 }
439 
440 void Tslice_dirac_max::set_trh(const Scalar& trh_in) {
441 
442  trh_evol.update(trh_in, jtime, the_time[jtime]) ;
443  cout << "Tslice_dirac_max::set_trh : #### WARNING : \n"
444  << " this method does not check whether det(tilde gamma) = 1"
445  << endl ;
446 
447  // Reset of quantities depending on the trace:
448  hh_evol.downdate(jtime) ;
449  if (p_tgamma != 0x0) {
450  delete p_tgamma ;
451  p_tgamma = 0x0 ;
452  }
453  if (p_hdirac != 0x0) {
454  delete p_hdirac ;
455  p_hdirac = 0x0 ;
456  }
457  if (p_gamma != 0x0) {
458  delete p_gamma ;
459  p_gamma = 0x0 ;
460  }
461  source_A_hh_evol.downdate(jtime) ;
462  source_B_hh_evol.downdate(jtime) ;
463  source_A_hata_evol.downdate(jtime) ;
464  source_B_hata_evol.downdate(jtime) ;
465  gam_dd_evol.downdate(jtime) ;
466  gam_uu_evol.downdate(jtime) ;
467  adm_mass_evol.downdate(jtime) ;
468 
469 }
470 
471 
472  //----------------------------------------------------//
473  // Update of fields from base class Time_slice_conf //
474  //----------------------------------------------------//
475 
476 
477 const Sym_tensor& Tslice_dirac_max::hh(Param* par_bc, Param* par_mat) const {
478 
479  if (!( hh_evol.is_known(jtime) ) ) {
480 
481  assert (A_hh_evol.is_known(jtime)) ;
482  assert (B_hh_evol.is_known(jtime)) ;
483 
484  // Computation of h^{ij} to ensure det tgam_{ij} = det f_{ij} :
485 
486  hh_det_one(jtime, par_bc, par_mat) ;
487  }
488 
489  return hh_evol[jtime] ;
490 
491 }
492 
493 
495 
496  if( !(trk_evol.is_known(jtime)) ) {
497 
498  Scalar resu(ff.get_mp()) ;
499  resu.set_etat_zero() ;
500 
501  trk_evol.update(resu, jtime, the_time[jtime]) ;
502 
503  }
504 
505  return trk_evol[jtime] ;
506 
507 }
508 
509 
511 
512  if (p_hdirac == 0x0) {
513  p_hdirac = new Vector(ff.get_mp(), CON, ff.get_triad() ) ;
515  }
516 
517  return *p_hdirac ;
518 
519 }
520 
521 
522 
523 
524  //-----------------------------------//
525  // Update of fields from this class //
526  //-----------------------------------//
527 
529 
530  if (!( A_hh_evol.is_known(jtime) ) ) {
531  assert( hh_evol.is_known(jtime) ) ;
532 
533  A_hh_evol.update( hh_evol[jtime].compute_A(true), jtime, the_time[jtime] ) ;
534 
535  }
536  return A_hh_evol[jtime] ;
537 }
538 
540 
541  if (!( B_hh_evol.is_known(jtime) ) ) {
542  assert( hh_evol.is_known(jtime) ) ;
543 
544  B_hh_evol.update( hh_evol[jtime].compute_tilde_B_tt(true), jtime, the_time[jtime] ) ;
545 
546  }
547  return B_hh_evol[jtime] ;
548 }
549 
551 
552  if( !(trh_evol.is_known(jtime)) ) {
553 
554  // Computation of tr(h) to ensure det tgam_{ij} = det f_{ij} :
555  hh_det_one(jtime) ;
556 
557  }
558 
559  return trh_evol[jtime] ;
560 
561 }
562 
563 
564 
565  //------------------//
566  // output //
567  //------------------//
568 
569 ostream& Tslice_dirac_max::operator>>(ostream& flux) const {
570 
572 
573  flux << "Dirac gauge and maximal slicing" << '\n' ;
574 
575  if (A_hh_evol.is_known(jtime)) {
576  maxabs( A_hh_evol[jtime], "A_hh", flux) ;
577  }
578  if (B_hh_evol.is_known(jtime)) {
579  maxabs( B_hh_evol[jtime], "B_hh", flux) ;
580  }
581  if (trh_evol.is_known(jtime)) {
582  maxabs( trh_evol[jtime], "tr h", flux) ;
583  }
584 
585  return flux ;
586 
587 }
588 
589 
590 void Tslice_dirac_max::sauve(FILE* fich, bool partial_save) const {
591 
592  if (partial_save) {
593  cout <<
594  "Tslice_dirac_max::sauve : the partial_save case is not ready yet !"
595  << endl ;
596  abort() ;
597  }
598 
599  // Writing of quantities common to all derived classes of Time_slice_conf
600  // ----------------------------------------------------------------------
601 
602  Time_slice_conf::sauve(fich, true) ;
603 
604  // Writing of the other fields
605  // ---------------------------
606 
607  int jmin = jtime - depth + 1 ;
608 
609  // h^{ij}
610  assert( hh_evol.is_known(jtime) ) ;
611  for (int j=jmin; j<=jtime; j++) {
612  int indicator = (hh_evol.is_known(j)) ? 1 : 0 ;
613  fwrite_be(&indicator, sizeof(int), 1, fich) ;
614  if (indicator == 1) hh_evol[j].sauve(fich) ;
615  }
616 
617  // A_hh
618  A_hh() ; // forces the update at the current time step
619  for (int j=jmin; j<=jtime; j++) {
620  int indicator = (A_hh_evol.is_known(j)) ? 1 : 0 ;
621  fwrite_be(&indicator, sizeof(int), 1, fich) ;
622  if (indicator == 1) A_hh_evol[j].sauve(fich) ;
623  }
624 
625  // B_hh
626  B_hh() ; // forces the update at the current time step
627  for (int j=jmin; j<=jtime; j++) {
628  int indicator = (B_hh_evol.is_known(j)) ? 1 : 0 ;
629  fwrite_be(&indicator, sizeof(int), 1, fich) ;
630  if (indicator == 1) B_hh_evol[j].sauve(fich) ;
631  }
632 }
633 
634 
635 
636 
637 
638 
639 
640 }
void initialize_sources_copy() const
Copy the sources source_A_XXX_evol and source_B_XXX_evol to all time-steps.
void annule_domain(int l)
Sets the Tensor to zero in a given domain.
Definition: tensor.C:675
virtual const Scalar & psi() const
Conformal factor relating the physical metric to the conformal one: .
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.
const Sym_tensor get_aa() const
Returns .
const Scalar & psi4() const
Factor at the current time step (jtime ).
Cmp exp(const Cmp &)
Exponential.
Definition: cmp_math.C:273
void operator=(const Time_slice_conf &)
Assignment to another Time_slice_conf.
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:520
virtual const Vector & hdirac() const
Vector which vanishes in Dirac gauge.
virtual void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition: scalar.C:330
virtual void set_trh(const Scalar &trh_in)
Sets the trace, with respect to the flat metric ff , of .
virtual void del_deriv() const
Deletes all the derived quantities.
const Base_vect_spher & get_bvect_spher() const
Returns the orthonormal vectorial basis associated with the coordinates of the mapping.
Definition: map.h:795
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator<<).
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
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...
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:777
Flat metric for tensor calculation.
Definition: metric.h:261
Tensor field of valence 0 (or component of a tensorial field).
Definition: scalar.h:393
Evolution_std< Scalar > A_hata_evol
Potential A associated with the symmetric tensor .
Definition: time_slice.h:550
Base class for coordinate mappings.
Definition: map.h:682
virtual void update(const TyT &new_value, int j, double time_j)
Sets a new value at a given time step.
const Base_vect * get_triad() const
Returns the vectorial basis (triad) on which the metric is defined.
Definition: metric.h:309
void operator=(const Tslice_dirac_max &)
Assignment to another Tslice_dirac_max.
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
virtual void set_hh(const Sym_tensor &hh_in)
Sets the deviation of the conformal metric from the flat metric : .
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...
int get_etat() const
Returns the logical state ETATNONDEF (undefined), ETATZERO (null) or ETATQCQ (ordinary).
Definition: scalar.h:560
virtual void set_hh(const Sym_tensor &hh_in)
Sets the deviation of the conformal metric from the flat metric : .
Tensor field of valence 1.
Definition: vector.h:188
Vectorial bases (triads) with respect to which the tensorial components are defined.
Definition: base_vect.h:105
Vector * p_hdirac
Pointer on the vector (which vanishes in Dirac gauge), at the current time step (jtime).
Definition: time_slice.h:574
const Map & get_mp() const
Returns the mapping.
Definition: metric.h:202
Evolution_std< Scalar > B_hata_evol
Potential associated with the symmetric tensor .
Definition: time_slice.h:555
void set_dzpuis(int)
Modifies the dzpuis flag.
Definition: scalar.C:814
const Scalar & get_ln_psi() const
Returns .
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator<<).
Evolution_std< Scalar > source_A_hata_evol
The potential A of the source of equation for .
Definition: time_slice.h:1004
Evolution_std< Scalar > source_A_hh_evol
The A potential of the source of equation for .
Definition: time_slice.h:992
virtual const Scalar & A_hata() const
Returns the potential A of .
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 ~Tslice_dirac_max()
Destructor.
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
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:211
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:465
Spacelike time slice of a 3+1 spacetime with conformal decomposition.
Definition: time_slice.h:501
virtual const Scalar & B_hh() const
Returns the potential of .
void trace_from_det_one(const Sym_tensor_tt &htt, double precis=1.e-14, int it_max=100)
Assigns the derived member p_tt and computes the trace so that *this + the flat metric has a determin...
Evolution_std< Scalar > source_B_hh_evol
The potential of the source of equation for .
Definition: time_slice.h:998
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
Class for relativistic rotating stars in Dirac gauge and maximal slicing.
virtual void sauve(FILE *fich, bool partial_save) const
Total or partial saves in a binary file.
void set_khi_mu(const Scalar &khi_i, const Scalar &mu_i, int dzp=0, Param *par1=0x0, Param *par2=0x0, Param *par3=0x0)
Sets the component , as well as the angular potential (see member p_khi and p_mu )...
Transverse symmetric tensors of rank 2.
Definition: sym_tensor.h:611
Metric * p_tgamma
Pointer on the conformal metric at the current time step (jtime)
Definition: time_slice.h:563
Evolution_std< double > the_time
Time label of each slice.
Definition: time_slice.h:196
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...
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
const Metric_flat & ff
Pointer on the flat metric with respect to which the conformal decomposition is performed.
Definition: time_slice.h:510
Evolution_std< Scalar > B_hh_evol
The potential of .
Definition: time_slice.h:986
void set_longit_trans(const Vector &v, const Sym_tensor_trans &a)
Assigns the derived members p_longit_pot and p_transverse and updates the components accordingly...
virtual const Scalar & A_hh() const
Returns the potential A of .
Evolution_std< Scalar > trk_evol
Values at successive time steps of the trace K of the extrinsic curvature.
Definition: time_slice.h:227
Tslice_dirac_max(const Scalar &lapse_in, const Vector &shift_in, const Metric_flat &ff_in, const Scalar &psi_in, const Sym_tensor_trans &hh_in, const Sym_tensor &hata_in, int depth_in=3)
Constructor from conformal decomposition.
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 const Scalar & B_hata() const
Returns the potential of .
virtual const Sym_tensor & hh(Param *par_bc=0x0, Param *par_mat=0x0) const
Deviation of the conformal metric from the flat metric : .
Evolution_std< Sym_tensor > hata_evol
Values at successive time steps of the components .
Definition: time_slice.h:545
int get_type_r(int l) const
Returns the type of sampling in the radial direction in domain no.
Definition: grilles.h:491
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< Vector > beta_evol
Values at successive time steps of the shift vector .
Definition: time_slice.h:222
virtual void sauve(FILE *fich, bool partial_save) const
Total or partial saves in a binary file.
Evolution_std< Scalar > trh_evol
The trace, with respect to the flat metric ff , of .
Definition: time_slice.h:1013
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 void dec_dzpuis(int dec=1)
Decreases by dec units the value of dzpuis and changes accordingly the values of the Scalar in the co...
Spacelike time slice of a 3+1 spacetime with conformal decomposition in the maximal slicing and Dirac...
Definition: time_slice.h:971
Evolution_std< Sym_tensor > hh_evol
Values at successive time steps of the components .
Definition: time_slice.h:533
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
Evolution_std< Scalar > source_B_hata_evol
The potential of the source of equation for .
Definition: time_slice.h:1010
Transverse and traceless symmetric tensors of rank 2.
Definition: sym_tensor.h:933
virtual const Scalar & trh() const
Computes the trace h, with respect to the flat metric ff , of .
Evolution_std< Scalar > A_hh_evol
The A potential of .
Definition: time_slice.h:980