LORENE
single_hor.C
1 /*
2  * Methods of class single_hor
3  *
4  * (see file isol_hor.h for documentation).
5  *
6  */
7 
8 /*
9  * Copyright (c) 2004 Jose Luis Jaramillo
10  * Francois Limousin
11  *
12  * This file is part of LORENE.
13  *
14  * LORENE is free software; you can redistribute it and/or modify
15  * it under the terms of the GNU General Public License version 2
16  * as published by the Free Software Foundation.
17  *
18  * LORENE is distributed in the hope that it will be useful,
19  * but WITHOUT ANY WARRANTY; without even the implied warranty of
20  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
21  * GNU General Public License for more details.
22  *
23  * You should have received a copy of the GNU General Public License
24  * along with LORENE; if not, write to the Free Software
25  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
26  *
27  */
28 
29 char single_hor_C[] = "$Header: /cvsroot/Lorene/C++/Source/Isol_hor/single_hor.C,v 1.3 2014/10/13 08:53:01 j_novak Exp $" ;
30 
31 /*
32  * $Id: single_hor.C,v 1.3 2014/10/13 08:53:01 j_novak Exp $
33  * $Log: single_hor.C,v $
34  * Revision 1.3 2014/10/13 08:53:01 j_novak
35  * Lorene classes and functions now belong to the namespace Lorene.
36  *
37  * Revision 1.2 2014/10/06 15:13:11 j_novak
38  * Modified #include directives to use c++ syntax.
39  *
40  * Revision 1.1 2007/04/13 15:28:35 f_limousin
41  * Lots of improvements, generalisation to an arbitrary state of
42  * rotation, implementation of the spatial metric given by Samaya.
43  *
44  *
45  * $Header: /cvsroot/Lorene/C++/Source/Isol_hor/single_hor.C,v 1.3 2014/10/13 08:53:01 j_novak Exp $
46  *
47  */
48 
49 // C headers
50 #include <cstdlib>
51 #include <cassert>
52 
53 // Lorene headers
54 #include "param.h"
55 #include "utilitaires.h"
56 #include "time_slice.h"
57 #include "isol_hor.h"
58 #include "tensor.h"
59 #include "metric.h"
60 #include "evolution.h"
61 //#include "graphique.h"
62 
63 //--------------//
64 // Constructors //
65 //--------------//
66 // Standard constructor
67 // --------------------
68 
69 namespace Lorene {
71  mp(mpi), nz(mpi.get_mg()->get_nzone()), radius ((mpi.get_alpha())[0]),
72  omega(0),regul(0),
73  n_auto(mpi), n_comp(mpi), nn(mpi),
74  psi_auto(mpi), psi_comp(mpi), psi(mpi),
75  dn(mpi, COV, mpi.get_bvect_spher()),
76  dpsi(mpi, COV, mpi.get_bvect_spher()),
77  beta_auto(mpi, CON, mpi.get_bvect_spher()),
78  beta_comp(mpi, CON, mpi.get_bvect_spher()),
79  beta(mpi, CON, mpi.get_bvect_spher()),
80  aa_auto(mpi, CON, mpi.get_bvect_spher()),
81  aa_comp(mpi, CON, mpi.get_bvect_spher()),
82  aa(mpi, CON, mpi.get_bvect_spher()),
83  tgam(mpi.flat_met_spher()),
84  ff(mpi.flat_met_spher()),
85  hh(mpi, CON, mpi.get_bvect_spher()),
86  gamt_point(mpi, CON, mpi.get_bvect_spher()),
87  trK(mpi), trK_point(mpi), decouple(mpi){
88 
89  hh.set_etat_zero() ;
90  set_der_0x0() ;
91 }
92 
93 // Copy constructor
94 // ----------------
95 
96 Single_hor::Single_hor(const Single_hor& singlehor_in)
97  : mp(singlehor_in.mp),
98  nz(singlehor_in.nz),
99  radius(singlehor_in.radius),
100  omega(singlehor_in.omega),
101  regul(singlehor_in.regul),
102  n_auto(singlehor_in.n_auto),
103  n_comp(singlehor_in.n_comp),
104  nn(singlehor_in.nn),
105  psi_auto(singlehor_in.psi_auto),
106  psi_comp(singlehor_in.psi_comp),
107  psi(singlehor_in.psi),
108  dn(singlehor_in.dn),
109  dpsi(singlehor_in.dpsi),
110  beta_auto(singlehor_in.beta_auto),
111  beta_comp(singlehor_in.beta_comp),
112  beta(singlehor_in.beta),
113  aa_auto(singlehor_in.aa_auto),
114  aa_comp(singlehor_in.aa_comp),
115  aa(singlehor_in.aa),
116  tgam(singlehor_in.tgam),
117  ff(singlehor_in.ff),
118  hh(singlehor_in.hh),
119  gamt_point(singlehor_in.gamt_point),
120  trK(singlehor_in.trK),
121  trK_point(singlehor_in.trK_point),
122  decouple(singlehor_in.decouple){
123 
124  set_der_0x0() ;
125 }
126 
127 // Constructor from a file
128 // -----------------------
129 
131  : mp(mpi), nz(mpi.get_mg()->get_nzone()), radius ((mpi.get_alpha())[0]),
132  omega(0),regul(0),
133  n_auto(mpi, *(mpi.get_mg()), fich), n_comp(mpi),
134  nn(mpi),
135  psi_auto(mpi, *(mpi.get_mg()), fich), psi_comp(mpi),
136  psi(mpi),
137  dn(mpi, COV, mpi.get_bvect_spher()),
138  dpsi(mpi, COV, mpi.get_bvect_spher()),
139  beta_auto(mpi, mpi.get_bvect_spher(), fich),
140  beta_comp(mpi, CON, mpi.get_bvect_spher()),
141  beta(mpi, CON, mpi.get_bvect_spher()),
142  aa_auto(mpi, CON, mpi.get_bvect_spher()),
143  aa_comp(mpi, CON, mpi.get_bvect_spher()),
144  aa(mpi, CON, mpi.get_bvect_spher()),
145  tgam(mpi.flat_met_spher()),
146  ff(mpi.flat_met_spher()),
147  hh(mpi, CON, mpi.get_bvect_spher()),
148  gamt_point(mpi, CON, mpi.get_bvect_spher()),
149  trK(mpi), trK_point(mpi), decouple(mpi){
150 
151  fread_be(&omega, sizeof(double), 1, fich) ;
152 
153  // tgam, gamt_point, trK, trK_point
154 
155  Sym_tensor met_file (mp, mp.get_bvect_spher(), fich) ;
156  tgam = met_file ;
157 
158  Sym_tensor gamt_point_file (mp, mp.get_bvect_spher(), fich) ;
159  gamt_point = gamt_point_file ;
160 
161  Scalar trK_file (mp, *(mp.get_mg()), fich) ;
162  trK = trK_file ;
163 
164  Scalar trK_point_file (mp, *(mp.get_mg()), fich) ;
165  trK_point = trK_point_file ;
166 
167  set_der_0x0() ;
168  hh = tgam.con() - ff.con() ;
169 
170 }
171 
172  //--------------//
173  // Destructor //
174  //--------------//
175 
177 
179 }
180 
181 
182  //-----------------------//
183  // Mutators / assignment //
184  //-----------------------//
185 
186 void Single_hor::operator=(const Single_hor& singlehor_in) {
187 
188  mp = singlehor_in.mp ;
189  nz = singlehor_in.nz ;
190  radius = singlehor_in.radius ;
191  omega = singlehor_in.omega ;
192  regul = singlehor_in.regul ;
193  n_auto = singlehor_in.n_auto ;
194  n_comp = singlehor_in.n_comp ;
195  nn = singlehor_in.nn ;
196  psi_auto = singlehor_in.psi_auto ;
197  psi_comp = singlehor_in.psi_comp ;
198  psi = singlehor_in.psi ;
199  dn = singlehor_in.dn ;
200  dpsi = singlehor_in.dpsi ;
201  beta_auto = singlehor_in.beta_auto ;
202  beta_comp = singlehor_in.beta_comp ;
203  beta = singlehor_in.beta ;
204  aa_auto = singlehor_in.aa_auto ;
205  aa_comp = singlehor_in.aa_comp ;
206  aa = singlehor_in.aa ;
207  tgam = singlehor_in.tgam ;
208  ff = singlehor_in.ff ;
209  hh = singlehor_in.hh ;
210  gamt_point = singlehor_in.gamt_point ;
211  trK = singlehor_in.trK ;
212  trK_point = singlehor_in.trK_point ;
213  decouple = singlehor_in.decouple ;
214 }
215 
216 
217  //---------------------//
218  // Memory management //
219  //---------------------//
220 
221 void Single_hor::del_deriv() const {
222 
223  if (p_psi4 != 0x0) delete p_psi4 ;
224  if (p_gam != 0x0) delete p_gam ;
225  if (p_k_dd != 0x0) delete p_k_dd ;
226 
227  set_der_0x0() ;
228 }
229 
230 
232 
233  p_psi4 = 0x0 ;
234  p_gam = 0x0 ;
235  p_k_dd = 0x0 ;
236 
237 }
238 
239  //--------------------------//
240  // Save in a file //
241  //--------------------------//
242 
243 
244 void Single_hor::sauve(FILE* fich) const {
245 
246  n_auto.sauve(fich) ;
247  psi_auto.sauve(fich) ;
248  beta_auto.sauve(fich) ;
249 
250  fwrite_be (&omega, sizeof(double), 1, fich) ;
251 
252  tgam.con().sauve(fich) ;
253  gamt_point.sauve(fich) ;
254  trK.sauve(fich) ;
255  trK_point.sauve(fich) ;
256 }
257 
258 // Accessors
259 // ---------
260 
262 
263  return n_auto ;
264 }
265 
267 
268  return n_comp ;
269 }
270 const Scalar& Single_hor::get_nn() const {
271 
272  return nn ;
273 }
274 
276 
277  return psi_auto ;
278 }
279 
281 
282  return psi_comp ;
283 }
284 const Scalar& Single_hor::get_psi() const {
285 
286  return psi ;
287 }
288 const Scalar& Single_hor::get_psi4() const {
289 
290  if (p_psi4 == 0x0) {
291 
292  p_psi4 = new Scalar( pow( psi, 4.) ) ;
294  }
295 
296  return *p_psi4 ;
297 }
298 
299 const Vector& Single_hor::get_dn() const {
300 
301  return dn ;
302 }
303 
304 const Vector& Single_hor::get_dpsi() const {
305 
306  return dpsi ;
307 }
308 
310 
311  return beta_auto ;
312 }
313 
315 
316  return beta_comp ;
317 
318 }
319 const Vector& Single_hor::get_beta() const {
320 
321  return beta ;
322 }
323 
325 
326  return aa_auto ;
327 }
328 
330 
331  return aa_comp ;
332 
333 }
335 
336  return aa ;
337 
338 }
339 const Metric& Single_hor::get_gam() const {
340 
341  if (p_gam == 0x0) {
342  p_gam = new Metric( get_psi4()*tgam.cov() ) ;
343  }
344 
345  return *p_gam ;
346 }
347 
349 
350  if (p_k_dd == 0x0) {
351 
352  Sym_tensor temp (aa/get_psi4()+1./3.*trK*get_gam().con()) ;
353 
354  p_k_dd = new Sym_tensor( temp.up_down(get_gam()) ) ;
356  }
357 
358  return *p_k_dd ;
359 }
360 
361 
362 void Single_hor::set_psi_auto(const Scalar& psi_in) {
363 
364  psi_auto = psi_in ;
365 }
366 void Single_hor::set_n_auto(const Scalar& n_in) {
367 
368  n_auto = n_in ;
369 }
370 void Single_hor::set_beta_auto(const Scalar& beta_in) {
371 
372  beta_auto = beta_in ;
373 }
374 void Single_hor::set_aa_auto(const Scalar& aa_in) {
375 
376  aa_auto = aa_in ;
377 }
378 void Single_hor::set_aa_comp(const Scalar& aa_in) {
379 
380  aa_comp = aa_in ;
381 }
382 void Single_hor::set_aa(const Scalar& aa_in) {
383 
384  aa = aa_in ;
385 }
386 
387 
388 // Import the lapse from the companion (Bhole case)
389 
391 
392  Scalar temp (mp) ;
393  temp.import(comp.n_auto) ;
394  temp.std_spectral_base() ;
395  n_comp = temp ;
396  nn = temp + n_auto ;
397 
398  Vector dn_comp (mp, COV, mp.get_bvect_cart()) ;
399  dn_comp.set_etat_qcq() ;
400  Vector auxi (comp.n_auto.derive_cov(comp.ff)) ;
401  auxi.dec_dzpuis(2) ;
402  auxi.change_triad(auxi.get_mp().get_bvect_cart()) ;
403  for (int i=1 ; i<=3 ; i++){
404  if (auxi(i).get_etat() != ETATZERO)
405  auxi.set(i).raccord(3) ;
406  }
407 
408  auxi.change_triad(mp.get_bvect_cart()) ;
409  assert ( *(auxi.get_triad()) == *(dn_comp.get_triad())) ;
410 
411  for (int i=1 ; i<=3 ; i++){
412  dn_comp.set(i).import(auxi(i)) ;
413  dn_comp.set(i).set_spectral_va().set_base(auxi(i).get_spectral_va().
414  get_base()) ;
415  }
416  dn_comp.inc_dzpuis(2) ;
417  dn_comp.change_triad(mp.get_bvect_spher()) ;
418 
419  dn = n_auto.derive_cov(ff) + dn_comp ;
420 }
421 
422 // Import the conformal factor from the companion (Bhole case)
423 
425 
426  Scalar temp (mp) ;
427  temp.import(comp.psi_auto) ;
428  temp.std_spectral_base() ;
429  psi_comp = temp ;
430  psi = temp + psi_auto ;
431 
432  Vector dpsi_comp (mp, COV, mp.get_bvect_cart()) ;
433  dpsi_comp.set_etat_qcq() ;
434  Vector auxi (comp.psi_auto.derive_cov(comp.ff)) ;
435  auxi.dec_dzpuis(2) ;
436  auxi.change_triad(auxi.get_mp().get_bvect_cart()) ;
437  for (int i=1 ; i<=3 ; i++){
438  if (auxi(i).get_etat() != ETATZERO)
439  auxi.set(i).raccord(3) ;
440  }
441 
442  auxi.change_triad(mp.get_bvect_cart()) ;
443  assert ( *(auxi.get_triad()) == *(dpsi_comp.get_triad())) ;
444 
445  for (int i=1 ; i<=3 ; i++){
446  dpsi_comp.set(i).import(auxi(i)) ;
447  dpsi_comp.set(i).set_spectral_va().set_base(auxi(i).get_spectral_va().
448  get_base()) ;
449  }
450  dpsi_comp.inc_dzpuis(2) ;
451  dpsi_comp.change_triad(mp.get_bvect_spher()) ;
452  /*
453  Vector dpsi_comp_zec (psi_comp().derive_cov(ff)) ;
454  for (int i=1 ; i<=3 ; i++)
455  for (int l=nz-1 ; l<=nz-1 ; l++) {
456  if (dpsi_comp.set(i).get_etat() == ETATQCQ)
457  dpsi_comp.set(i).set_domain(l) = dpsi_comp_zec(i).domain(l) ;
458  }
459  */
460 
461  dpsi = psi_auto.derive_cov(ff) + dpsi_comp ;
462 
463 }
464 
466 
467  Vector tmp_vect (mp, CON, mp.get_bvect_cart()) ;
468  Vector shift_comp (comp.beta_auto) ;
469  shift_comp.change_triad(comp.mp.get_bvect_cart()) ;
470  shift_comp.change_triad(mp.get_bvect_cart()) ;
471  assert (*(shift_comp.get_triad()) == *(tmp_vect.get_triad())) ;
472 
473  tmp_vect.set(1).import(shift_comp(1)) ;
474  tmp_vect.set(2).import(shift_comp(2)) ;
475  tmp_vect.set(3).import(shift_comp(3)) ;
476  tmp_vect.std_spectral_base() ;
477  tmp_vect.change_triad(mp.get_bvect_spher()) ;
478 
479  beta_comp = tmp_vect ;
480  beta = beta_auto + beta_comp ;
481 }
482 
483 //Initialisation to Schwartzchild
485 
486  Scalar auxi(mp) ;
487 
488  // Initialisation of the lapse different of zero on the horizon
489  // at the first step
490  auxi = 0.5 - 0.5/mp.r ;
491  auxi.annule(0, 0);
492  auxi.set_dzpuis(0) ;
493 
494  Scalar temp(mp) ;
495  temp = auxi;
496  temp.std_spectral_base() ;
497  temp.raccord(1) ;
498  n_auto = temp ;
499 
500  temp = 0.5 ;
501  temp.std_spectral_base() ;
502  n_comp = temp ;
503  nn = n_auto ;
504 
505  auxi = 0.5 + radius/mp.r ;
506  auxi.annule(0, 0);
507  auxi.set_dzpuis(0) ;
508  temp = auxi;
509  temp.std_spectral_base() ;
510  temp.raccord(1) ;
511  psi_auto = temp ;
512 
513  temp = 0.5 ;
514  temp.std_spectral_base() ;
515  psi_comp = temp ;
516  psi = psi_auto + psi_comp ;
517 
518  dn = nn.derive_cov(ff) ;
519  dpsi = psi.derive_cov(ff) ;
520 
521  Vector temp_vect1(mp, CON, mp.get_bvect_spher()) ;
522  temp_vect1.set(1) = 0.0/mp.r/mp.r ;
523  temp_vect1.set(2) = 0. ;
524  temp_vect1.set(3) = 0. ;
525  temp_vect1.std_spectral_base() ;
526 
527  Vector temp_vect2(mp, CON, mp.get_bvect_spher()) ;
528  temp_vect2.set_etat_zero() ;
529 
530  beta_auto = temp_vect1 ;
531  beta_comp = temp_vect2 ;
532  beta = temp_vect1 ;
533 
534 }
535 
537 
538  Metric flat (mp.flat_met_spher()) ;
539  tgam = flat ;
540 
542  trK.set_etat_zero() ;
544 
545 }
546 
547 
549 
550  Scalar auxi(mp) ;
551 
552  auxi = (1-radius/mp.r)/(1+radius/mp.r) ;
553  auxi.annule(0, 0);
554  auxi.set_outer_boundary((*mp.get_mg()).get_nzone(), 1.) ;
555  auxi.set_dzpuis(0) ;
556 
557  Scalar temp(mp) ;
558  temp = auxi;
559  temp.std_spectral_base() ;
560  temp.raccord(1) ;
561  n_auto = temp;
562 
563  temp.set_etat_zero() ;
564  n_comp = temp ;
565  nn = temp ;
566 
567  auxi = 1 + radius/mp.r ;
568  auxi.annule(0, 0);
569  auxi.set_outer_boundary((*mp.get_mg()).get_nzone(), 1.) ;
570  auxi.set_dzpuis(0) ;
571 
572  temp = auxi;
573  temp.std_spectral_base() ;
574  temp.raccord(1) ;
575  psi_auto = temp ;
576  temp.set_etat_zero() ;
577  psi_comp = temp ;
578  psi = temp ;
579 
580  dn = nn.derive_cov(ff) ;
581  dpsi = psi.derive_cov(ff) ;
582 
583  Vector temp_vect(mp, CON, mp.get_bvect_spher()) ;
584  temp_vect.set_etat_zero() ;
585  beta_auto = temp_vect ;
586  beta_comp = temp_vect ;
587  beta = temp_vect ;
588 
589 }
590 
591 
592 
593 }
Lorene::Single_hor::dn
Vector dn
Covariant derivative of the lapse with respect to the flat metric .
Definition: isol_hor.h:937
Lorene::Scalar::raccord
void raccord(int n)
Performs the matching of the nucleus with respect to the first shell.
Definition: scalar_raccord.C:62
Lorene::Single_hor::~Single_hor
virtual ~Single_hor()
Destructor.
Definition: single_hor.C:176
Lorene::Metric_flat::con
virtual const Sym_tensor & con() const
Read-only access to the contravariant representation.
Definition: metric_flat.C:153
Lorene::Single_hor::get_dn
const Vector & get_dn() const
Covariant derivative of the lapse function .
Definition: single_hor.C:299
Lorene::Single_hor::nn
Scalar nn
Lapse function .
Definition: isol_hor.h:921
Lorene::Single_hor::get_beta
const Vector & get_beta() const
Shift function .
Definition: single_hor.C:319
Lorene::Tensor_sym::sauve
virtual void sauve(FILE *) const
Save in a binary file.
Definition: tensor_sym.C:372
Lorene::Scalar::sauve
virtual void sauve(FILE *) const
Save in a file.
Definition: scalar.C:686
Lorene::Scalar::set_outer_boundary
void set_outer_boundary(int l, double x)
Sets the value of the Scalar at the outer boundary of a given domain.
Definition: scalar_manip.C:315
Lorene::Single_hor::get_n_comp
const Scalar & get_n_comp() const
Lapse function .
Definition: single_hor.C:266
Lorene::Single_hor::ff
Metric_flat ff
3 metric flat
Definition: isol_hor.h:980
Lorene::Single_hor::beta_comp
Vector beta_comp
Shift function .
Definition: isol_hor.h:947
Lorene::Single_hor::psi
Scalar psi
Conformal factor .
Definition: isol_hor.h:930
Lorene::Single_hor::get_dpsi
const Vector & get_dpsi() const
Covariant derivative with respect to the flat metric of the conformal factor .
Definition: single_hor.C:304
Lorene::Single_hor::aa
Sym_tensor aa
Components of the conformal representation of the traceless part of the extrinsic curvature:
Definition: isol_hor.h:971
Lorene::Single_hor::set_beta_auto
void set_beta_auto(const Scalar &shift_in)
Sets the shift.
Definition: single_hor.C:370
Lorene::Sym_tensor
Class intended to describe valence-2 symmetric tensors.
Definition: sym_tensor.h:223
Lorene::Metric
Metric for tensor calculation.
Definition: metric.h:90
Lorene::Single_hor::set_aa_auto
void set_aa_auto(const Scalar &aa_auto_in)
Sets aa_auto.
Definition: single_hor.C:374
Lorene::Single_hor::n_auto
Scalar n_auto
Lapse function .
Definition: isol_hor.h:915
Lorene::Single_hor::get_aa_auto
const Sym_tensor & get_aa_auto() const
Conformal representation of the traceless part of the extrinsic curvature:
Definition: single_hor.C:324
Lorene::Single_hor::set_psi_auto
void set_psi_auto(const Scalar &psi_in)
Sets the conformal factor relating the physical metric to the conformal one: .
Definition: single_hor.C:362
Lorene::Single_hor::set_aa_comp
void set_aa_comp(const Scalar &aa_comp_in)
Sets aa_comp.
Definition: single_hor.C:378
Lorene
Lorene prototypes.
Definition: app_hor.h:64
Lorene::Tensor::set_etat_zero
virtual void set_etat_zero()
Sets the logical state of all components to ETATZERO (zero state).
Definition: tensor.C:497
Lorene::Single_hor::get_psi_comp
const Scalar & get_psi_comp() const
Conformal factor .
Definition: single_hor.C:280
Lorene::Single_hor::beta
Vector beta
Shift function .
Definition: isol_hor.h:950
Lorene::Map::get_bvect_spher
const Base_vect_spher & get_bvect_spher() const
Returns the orthonormal vectorial basis associated with the coordinates of the mapping.
Definition: map.h:783
Lorene::Map::get_mg
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:765
Lorene::Map::flat_met_spher
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:321
Lorene::Single_hor::aa_comp
Sym_tensor aa_comp
Components of the conformal representation of the traceless part of the extrinsic curvature:
Definition: isol_hor.h:965
Lorene::Single_hor::init_met_trK
void init_met_trK()
Sets the 3-metric tilde to the flat metric and gamt_point, trK and trK_point to zero.
Definition: single_hor.C:536
Lorene::Single_hor::sauve
virtual void sauve(FILE *fich) const
Total or partial saves in a binary file.
Definition: single_hor.C:244
Lorene::Tensor::dec_dzpuis
virtual void dec_dzpuis(int dec=1)
Decreases by dec units the value of dzpuis and changes accordingly the values in the compactified ext...
Definition: tensor.C:808
Lorene::Scalar::derive_cov
const Vector & derive_cov(const Metric &gam) const
Returns the gradient (1-form = covariant vector) of *this
Definition: scalar_deriv.C:390
Lorene::Metric::cov
virtual const Sym_tensor & cov() const
Read-only access to the covariant representation.
Definition: metric.C:280
Lorene::Map::r
Coord r
r coordinate centered on the grid
Definition: map.h:718
Lorene::Single_hor::trK_point
Scalar trK_point
Time derivative of the trace of the extrinsic curvature.
Definition: isol_hor.h:992
Lorene::Tensor::up_down
Tensor up_down(const Metric &gam) const
Computes a new tensor by raising or lowering all the indices of *this .
Definition: tensor_calculus.C:305
Lorene::Scalar::import
void import(const Scalar &ci)
Assignment to another Scalar defined on a different mapping.
Definition: scalar_import.C:68
Lorene::Scalar
Tensor field of valence 0 (or component of a tensorial field).
Definition: scalar.h:387
Lorene::Single_hor::aa_auto
Sym_tensor aa_auto
Components of the conformal representation of the traceless part of the extrinsic curvature:
Definition: isol_hor.h:959
Lorene::Tensor::get_mp
const Map & get_mp() const
Returns the mapping.
Definition: tensor.h:861
Lorene::pow
Cmp pow(const Cmp &, int)
Power .
Definition: cmp_math.C:348
Lorene::Vector::change_triad
virtual void change_triad(const Base_vect &)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
Definition: vector_change_triad.C:75
Lorene::Single_hor::regul
double regul
Intensity of the correction on the shift vector.
Definition: isol_hor.h:912
Lorene::Single_hor::get_psi4
const Scalar & get_psi4() const
Conformal factor .
Definition: single_hor.C:288
Lorene::Single_hor::dpsi
Vector dpsi
Covariant derivative of the conformal factor .
Definition: isol_hor.h:941
Lorene::Single_hor::get_beta_auto
const Vector & get_beta_auto() const
Shift function .
Definition: single_hor.C:309
Lorene::Single_hor::del_deriv
void del_deriv() const
Deletes all the derived quantities.
Definition: single_hor.C:221
Lorene::Scalar::set_spectral_va
Valeur & set_spectral_va()
Returns va (read/write version)
Definition: scalar.h:604
Lorene::Tensor::set_etat_qcq
virtual void set_etat_qcq()
Sets the logical state of all components to ETATQCQ (ordinary state).
Definition: tensor.C:481
Lorene::Single_hor::mp
Map_af & mp
Affine mapping.
Definition: isol_hor.h:900
Lorene::Single_hor::init_bhole
void init_bhole()
Sets the values of the fields to :
Definition: single_hor.C:484
Lorene::Single_hor::get_psi_auto
const Scalar & get_psi_auto() const
Conformal factor
Definition: single_hor.C:275
Lorene::Single_hor::get_nn
const Scalar & get_nn() const
Lapse function .
Definition: single_hor.C:270
Lorene::Single_hor::psi_auto
Scalar psi_auto
Conformal factor .
Definition: isol_hor.h:924
Lorene::Vector::set
Scalar & set(int)
Read/write access to a component.
Definition: vector.C:296
Lorene::Metric::con
virtual const Sym_tensor & con() const
Read-only access to the contravariant representation.
Definition: metric.C:290
Lorene::Single_hor::get_beta_comp
const Vector & get_beta_comp() const
Shift function .
Definition: single_hor.C:314
Lorene::Single_hor::nz
int nz
Number of zones.
Definition: isol_hor.h:903
Lorene::Single_hor::p_gam
Metric * p_gam
Spatial metric .
Definition: isol_hor.h:953
Lorene::Single_hor::p_k_dd
Sym_tensor * p_k_dd
Components of the extrinsic curvature:
Definition: isol_hor.h:974
Lorene::Single_hor::radius
double radius
Radius of the horizon in LORENE's units.
Definition: isol_hor.h:906
Lorene::Single_hor::operator=
void operator=(const Single_hor &)
Assignment to another Single_hor.
Definition: single_hor.C:186
Lorene::Tensor::std_spectral_base
virtual void std_spectral_base()
Sets the standard spectal bases of decomposition for each component.
Definition: tensor.C:926
Lorene::fread_be
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:69
Lorene::Single_hor::gamt_point
Sym_tensor gamt_point
Time derivative of the 3-metric tilde.
Definition: isol_hor.h:986
Lorene::Map_af
Affine radial mapping.
Definition: map.h:2027
Lorene::Valeur::set_base
void set_base(const Base_val &)
Sets the bases for spectral expansions (member base )
Definition: valeur.C:810
Lorene::Map::get_bvect_cart
const Base_vect_cart & get_bvect_cart() const
Returns the Cartesian basis associated with the coordinates (x,y,z) of the mapping,...
Definition: map.h:791
Lorene::Scalar::annule
virtual void annule(int l_min, int l_max)
Sets the Scalar to zero in several domains.
Definition: scalar.C:391
Lorene::Single_hor::get_n_auto
const Scalar & get_n_auto() const
Lapse function .
Definition: single_hor.C:261
Lorene::Scalar::set_etat_zero
virtual void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition: scalar.C:324
Lorene::Scalar::set_dzpuis
void set_dzpuis(int)
Modifies the dzpuis flag.
Definition: scalar.C:808
Lorene::Vector
Tensor field of valence 1.
Definition: vector.h:188
Lorene::Tensor::get_triad
const Base_vect * get_triad() const
Returns the vectorial basis (triad) on which the components are defined.
Definition: tensor.h:866
Lorene::Single_hor::psi_comp_import
void psi_comp_import(const Single_hor &comp)
Imports the part of due to the companion hole comp .
Definition: single_hor.C:424
Lorene::Single_hor::hh
Sym_tensor hh
Deviation metric.
Definition: isol_hor.h:983
Lorene::Single_hor::init_bhole_seul
void init_bhole_seul()
Initiates for a single black hole.
Definition: single_hor.C:548
Lorene::Single_hor::set_n_auto
void set_n_auto(const Scalar &nn_in)
Sets the lapse.
Definition: single_hor.C:366
Lorene::Tensor::inc_dzpuis
virtual void inc_dzpuis(int inc=1)
Increases by inc units the value of dzpuis and changes accordingly the values in the compactified ext...
Definition: tensor.C:816
Lorene::Single_hor::get_k_dd
const Sym_tensor & get_k_dd() const
k_dd
Definition: single_hor.C:348
Lorene::Scalar::std_spectral_base
virtual void std_spectral_base()
Sets the spectral bases of the Valeur va to the standard ones for a scalar field.
Definition: scalar.C:784
Lorene::Single_hor::decouple
Scalar decouple
Function used to construct from the total .
Definition: isol_hor.h:1002
Lorene::Single_hor::get_aa_comp
const Sym_tensor & get_aa_comp() const
Conformal representation of the traceless part of the extrinsic curvature:
Definition: single_hor.C:329
Lorene::Single_hor::psi_comp
Scalar psi_comp
Conformal factor .
Definition: isol_hor.h:927
Lorene::Single_hor::p_psi4
Scalar * p_psi4
Conformal factor .
Definition: isol_hor.h:933
Lorene::Single_hor::get_gam
const Metric & get_gam() const
metric
Definition: single_hor.C:339
Lorene::Single_hor::set_aa
void set_aa(const Scalar &aa_in)
Sets aa.
Definition: single_hor.C:382
Lorene::fwrite_be
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:70
Lorene::Single_hor::get_psi
const Scalar & get_psi() const
Conformal factor .
Definition: single_hor.C:284
Lorene::Single_hor::omega
double omega
Angular velocity in LORENE's units.
Definition: isol_hor.h:909
Lorene::Single_hor::get_aa
const Sym_tensor & get_aa() const
Conformal representation of the traceless part of the extrinsic curvature:
Definition: single_hor.C:334
Lorene::Single_hor
Binary black holes system.
Definition: isol_hor.h:894
Lorene::Single_hor::beta_auto
Vector beta_auto
Shift function .
Definition: isol_hor.h:944
Lorene::Single_hor::n_comp
Scalar n_comp
Lapse function .
Definition: isol_hor.h:918
Lorene::Single_hor::n_comp_import
void n_comp_import(const Single_hor &comp)
Imports the part of N due to the companion hole comp .
Definition: single_hor.C:390
Lorene::Single_hor::set_der_0x0
void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: single_hor.C:231
Lorene::Single_hor::beta_comp_import
void beta_comp_import(const Single_hor &comp)
Imports the part of due to the companion hole comp.
Definition: single_hor.C:465
Lorene::Single_hor::Single_hor
Single_hor(Map_af &mpi)
Standard constructor.
Definition: single_hor.C:70
Lorene::Single_hor::trK
Scalar trK
Trace of the extrinsic curvature.
Definition: isol_hor.h:989
Lorene::Single_hor::tgam
Metric tgam
3 metric tilde
Definition: isol_hor.h:977
Lorene::Tensor::sauve
virtual void sauve(FILE *) const
Save in a binary file.
Definition: tensor.C:906
Lorene::Vector::std_spectral_base
virtual void std_spectral_base()
Sets the standard spectal bases of decomposition for each component.
Definition: vector.C:316