LORENE
et_bin_nsbh.C
1 /*
2  * Methods for the class Et_bin_nsbh
3  *
4  * (see file et_bin_nsbh.h for documentation).
5  *
6  */
7 
8 /*
9  * Copyright (c) 2003 Keisuke Taniguchi
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 char et_bin_nsbh_C[] = "$Header: /cvsroot/Lorene/C++/Source/Etoile/et_bin_nsbh.C,v 1.12 2014/10/13 08:52:56 j_novak Exp $" ;
29 
30 /*
31  * $Id: et_bin_nsbh.C,v 1.12 2014/10/13 08:52:56 j_novak Exp $
32  * $Log: et_bin_nsbh.C,v $
33  * Revision 1.12 2014/10/13 08:52:56 j_novak
34  * Lorene classes and functions now belong to the namespace Lorene.
35  *
36  * Revision 1.11 2014/10/06 15:13:08 j_novak
37  * Modified #include directives to use c++ syntax.
38  *
39  * Revision 1.10 2006/09/25 10:01:49 p_grandclement
40  * Addition of N-dimensional Tbl
41  *
42  * Revision 1.9 2006/06/01 12:47:53 p_grandclement
43  * update of the Bin_ns_bh project
44  *
45  * Revision 1.8 2005/08/29 15:10:16 p_grandclement
46  * Addition of things needed :
47  * 1) For BBH with different masses
48  * 2) Provisory files for the mixted binaries (Bh and NS) : THIS IS NOT
49  * WORKING YET !!!
50  *
51  * Revision 1.7 2004/06/09 06:39:35 k_taniguchi
52  * Introduce set_n_auto() and set_confpsi_auto().
53  *
54  * Revision 1.6 2004/03/25 10:29:04 j_novak
55  * All LORENE's units are now defined in the namespace Unites (in file unites.h).
56  *
57  * Revision 1.5 2003/12/16 05:30:18 k_taniguchi
58  * Some changes in the constructor from file and the saving to file.
59  *
60  * Revision 1.4 2003/11/28 04:30:22 k_taniguchi
61  * Add the output operator.
62  *
63  * Revision 1.3 2003/10/24 12:30:58 k_taniguchi
64  * Add some derivative terms
65  *
66  * Revision 1.2 2003/10/22 08:12:46 k_taniguchi
67  * Supress some terms in Et_bin_nsbh::sauve.
68  *
69  * Revision 1.1 2003/10/21 11:50:38 k_taniguchi
70  * Methods for class Et_bin_nsbh
71  *
72  *
73  * $Header: /cvsroot/Lorene/C++/Source/Etoile/et_bin_nsbh.C,v 1.12 2014/10/13 08:52:56 j_novak Exp $
74  *
75  */
76 
77 // C headers
78 #include <cmath>
79 
80 // Lorene headers
81 #include "et_bin_nsbh.h"
82 #include "etoile.h"
83 #include "eos.h"
84 #include "utilitaires.h"
85 #include "param.h"
86 #include "unites.h"
87 
88  //--------------//
89  // Constructors //
90  //--------------//
91 
92 // Standard constructor
93 // --------------------
94 namespace Lorene {
95 Et_bin_nsbh::Et_bin_nsbh(Map& mp_i, int nzet_i, bool relat, const Eos& eos_i,
96  bool irrot, const Base_vect& ref_triad_i)
97  : Etoile_bin(mp_i, nzet_i, relat, eos_i, irrot, ref_triad_i),
98  n_auto(mp_i),
99  n_comp(mp_i),
100  d_n_auto(mp_i, 1, COV, ref_triad_i),
101  d_n_comp(mp_i, 1, COV, ref_triad_i),
102  confpsi(mp_i),
103  confpsi_auto(mp_i),
104  confpsi_comp(mp_i),
105  d_confpsi_auto(mp_i, 1, COV, ref_triad_i),
106  d_confpsi_comp(mp_i, 1, COV, ref_triad_i),
107  taij_auto(mp_i, 2, CON, ref_triad_i) ,
108  taij_comp(mp_i, 2, CON, ref_triad_i),
109  taij_tot(mp_i, 2, CON, ref_triad_i) ,
110  tkij_auto(mp_i, 2, CON, ref_triad_i),
111  tkij_tot(mp_i, 2, CON, ref_triad_i),
112  ssjm1_lapse(mp_i),
113  ssjm1_confpsi(mp_i) {
114 
115  // Pointers of derived quantities initialized to zero :
116  set_der_0x0() ;
117 
118  // The metric is initialized to the flat one :
119  n_auto = 0.5 ;
120  n_auto.set_std_base() ;
121  n_comp = 0.5 ;
122  n_comp.set_std_base() ;
123  d_n_auto = 0 ;
124  d_n_comp = 0 ;
125  confpsi = 1 ;
127  confpsi_auto = 0.5 ;
129  confpsi_comp = 0.5 ;
131  d_confpsi_auto = 0 ;
132  d_confpsi_comp = 0 ;
133 
139 
141  ssjm1_lapse = 0. ;
143  ssjm1_confpsi = 0. ;
144 
145 }
146 
147 // Copy constructor
148 // ----------------
150  : Etoile_bin(et),
151  n_auto(et.n_auto),
152  n_comp(et.n_comp),
153  d_n_auto(et.d_n_auto),
154  d_n_comp(et.d_n_comp),
155  confpsi(et.confpsi),
156  confpsi_auto(et.confpsi_auto),
157  confpsi_comp(et.confpsi_comp),
158  d_confpsi_auto(et.d_confpsi_auto),
159  d_confpsi_comp(et.d_confpsi_comp),
160  taij_auto(et.taij_auto),
161  taij_comp(et.taij_comp),
162  taij_tot(et.taij_tot),
163  tkij_auto(et.tkij_auto),
164  tkij_tot(et.tkij_tot),
165  ssjm1_lapse(et.ssjm1_lapse),
166  ssjm1_confpsi(et.ssjm1_confpsi) {
167 
168  set_der_0x0() ;
169 }
170 
171 // Constructor from a file
172 // -----------------------
173 Et_bin_nsbh::Et_bin_nsbh(Map& mp_i, const Eos& eos_i,
174  const Base_vect& ref_triad_i, FILE* fich, bool old)
175  : Etoile_bin(mp_i, eos_i, ref_triad_i, fich),
176  n_auto(mp_i),
177  n_comp(mp_i),
178  d_n_auto(mp_i, 1, COV, ref_triad_i),
179  d_n_comp(mp_i, 1, COV, ref_triad_i ),
180  confpsi(mp_i),
181  confpsi_auto(mp_i),
182  confpsi_comp(mp_i),
183  d_confpsi_auto(mp_i, 1, COV, ref_triad_i),
184  d_confpsi_comp(mp_i, 1, COV, ref_triad_i),
185  taij_auto(mp_i, 2, CON, ref_triad_i),
186  taij_comp(mp_i, 2, CON, ref_triad_i),
187  taij_tot(mp_i, 2, CON, ref_triad_i),
188  tkij_auto(mp_i, 2, CON, ref_triad_i),
189  tkij_tot(mp_i, 2, CON, ref_triad_i),
190  ssjm1_lapse(mp_i),
191  ssjm1_confpsi(mp_i) {
192 
193  // Construct from data in Etoile_bin
194  Cmp n_from_file (mp_i, *(mp_i.get_mg()), fich) ;
195  n_auto.set_etat_qcq() ;
196  n_auto.set() = n_from_file ;
197 
198  Cmp psi_from_file (mp_i, *(mp_i.get_mg()), fich) ;
200  confpsi_auto.set() = psi_from_file ;
201 
202  if (!old) {
203  Tenseur shift_auto_file (mp, ref_triad_i, fich) ;
204  shift_auto = shift_auto_file ;
205  }
206 
207  // All other fields are initialized to zero or some constants :
208  // ----------------------------------------------------------
209  n_comp = 0.5 ;
210  n_comp.set_std_base() ;
211  d_n_auto = 0 ;
212  d_n_comp = 0 ;
213  confpsi = 1 ;
215  confpsi_comp = 0.5 ;
217  d_confpsi_auto = 0 ;
218  d_confpsi_comp = 0 ;
224 
225  // Read of the saved fields :
226  Cmp ssjm1_lapse_file(mp_i, *(mp_i.get_mg()), fich) ;
227  ssjm1_lapse = ssjm1_lapse_file ;
228 
229  Cmp ssjm1_confpsi_file(mp_i, *(mp_i.get_mg()), fich) ;
230  ssjm1_confpsi = ssjm1_confpsi_file ;
231 
232  // Pointers of derived quantities initialized to zero
233  // --------------------------------------------------
234  set_der_0x0() ;
235 }
236 
237  //------------//
238  // Destructor //
239  //------------//
240 
241 Et_bin_nsbh::~Et_bin_nsbh(){
242 
243  del_deriv() ;
244 
245 }
246 
247  //--------------//
248  // Assignment //
249  //--------------//
250 
251 // Assignment to another Et_bin_nsbh
252 // ---------------------------------
254 
255  // Assignment of quantities common to the derived classes of Etoile_bin
257 
258  n_auto = et.n_auto ;
259  n_comp = et.n_comp ;
260  d_n_auto = et.d_n_auto ;
261  d_n_comp = et.d_n_comp ;
262  confpsi = et.confpsi ;
267  taij_auto = et.taij_auto ;
268  taij_comp = et.taij_comp ;
269  taij_tot = et.taij_tot ;
270  tkij_auto = et.tkij_auto ;
271  tkij_tot = et.tkij_tot ;
272  ssjm1_lapse = et.ssjm1_lapse ;
274 
275  del_deriv() ; // Deletes all derived quantities
276 }
277 
279 
280  del_deriv() ; // sets to 0x0 all the derived quantities
281  return n_auto ;
282 
283 }
284 
286 
287  del_deriv() ; // sets to 0x0 all the derived quantities
288  return n_comp ;
289 
290 }
291 
293 
294  del_deriv() ; // sets to 0x0 all the derived quantities
295  return confpsi_auto ;
296 
297 }
298 
300 
301  del_deriv() ; // sets to 0x0 all the derived quantities
302  return confpsi_comp ;
303 
304 }
305 
306  //--------------//
307  // Outputs //
308  //--------------//
309 
310 // Save in a file
311 // --------------
312 void Et_bin_nsbh::sauve(FILE* fich) const {
313 
314  Etoile_bin::sauve(fich) ;
315 
316  n_auto().sauve(fich) ;
317  confpsi_auto().sauve(fich) ;
318  shift_auto.sauve(fich) ;
319 
320  ssjm1_lapse.sauve(fich) ;
321  ssjm1_confpsi.sauve(fich) ;
322 }
323 
324 
325 // Printing
326 // --------
327 
328 ostream& Et_bin_nsbh::operator>>(ostream& ost) const {
329 
330  using namespace Unites ;
331 
332  Etoile::operator>>(ost) ;
333 
334  ost << endl ;
335  ost << "Neutron star in a binary system" << endl ;
336  ost << "-------------------------------" << endl ;
337 
338  if (irrotational) {
339  ost << "irrotational configuration" << endl ;
340  }
341  else {
342  ost << "corotating configuration" << endl ;
343  }
344 
345  ost << "Absolute abscidia of the stellar center: " <<
346  mp.get_ori_x() / km << " km" << endl ;
347 
348  ost << "Absolute abscidia of the barycenter of the baryon density : " <<
349  xa_barycenter() / km << " km" << endl ;
350 
351  double r_0 = 0.5 * ( ray_eq() + ray_eq_pi() ) ;
352  double d_ns = fabs( mp.get_ori_x() ) + ray_eq_pi() - r_0 ;
353  double d_tilde = 2 * d_ns / r_0 ;
354 
355  ost << "d_tilde : " << d_tilde << endl ;
356 
357  ost << "Orientation with respect to the absolute frame : " <<
358  mp.get_rot_phi() << " rad" << endl ;
359 
360  ost << "Central value of gam_euler : "
361  << gam_euler()(0, 0, 0, 0) << endl ;
362 
363  ost << "Central u_euler (U^X, U^Y, U^Z) [c] : "
364  << u_euler(0)(0, 0, 0, 0) << " "
365  << u_euler(1)(0, 0, 0, 0) << " "
366  << u_euler(2)(0, 0, 0, 0) << endl ;
367 
368  if (irrotational) {
369  ost << "Central d_psi (X, Y, Z) [c] : "
370  << d_psi(0)(0, 0, 0, 0) << " "
371  << d_psi(1)(0, 0, 0, 0) << " "
372  << d_psi(2)(0, 0, 0, 0) << endl ;
373 
374  ost << "Central vel. / co-orb. (W^X, W^Y, W^Z) [c] : "
375  << wit_w(0)(0, 0, 0, 0) << " "
376  << wit_w(1)(0, 0, 0, 0) << " "
377  << wit_w(2)(0, 0, 0, 0) << endl ;
378 
379  ost << "Max vel. / co-orb. (W^X, W^Y, W^Z) [c] : "
380  << max(max(wit_w(0))) << " "
381  << max(max(wit_w(1))) << " "
382  << max(max(wit_w(2))) << endl ;
383 
384  ost << "Min vel. / co-orb. (W^X, W^Y, W^Z) [c] : "
385  << min(min(wit_w(0))) << " "
386  << min(min(wit_w(1))) << " "
387  << min(min(wit_w(2))) << endl ;
388 
389  double r_surf = mp.val_r(0,1.,M_PI/4,M_PI/4) ;
390 
391  ost << "Velocity at (r_surf,pi/4,pi/4) / co-orb. [c] : "
392  << wit_w(0).val_point(r_surf,M_PI/4,M_PI/4) << " "
393  << wit_w(1).val_point(r_surf,M_PI/4,M_PI/4) << " "
394  << wit_w(2).val_point(r_surf,M_PI/4,M_PI/4) << endl ;
395 
396  ost << "Central value of loggam : "
397  << loggam()(0, 0, 0, 0) << endl ;
398  }
399 
400  ost << "Central value of lapse(N) auto : "
401  << n_auto()(0, 0, 0, 0) << endl ;
402 
403  ost << "Central value of confpsi auto : "
404  << confpsi_auto()(0, 0, 0, 0) << endl ;
405 
406  ost << "Central value of shift (N^X, N^Y, N^Z) [c] : "
407  << shift(0)(0, 0, 0, 0) << " "
408  << shift(1)(0, 0, 0, 0) << " "
409  << shift(2)(0, 0, 0, 0) << endl ;
410 
411  ost << " ... shift_auto part of it [c] : "
412  << shift_auto(0)(0, 0, 0, 0) << " "
413  << shift_auto(1)(0, 0, 0, 0) << " "
414  << shift_auto(2)(0, 0, 0, 0) << endl ;
415 
416  ost << " ... shift_comp part of it [c] : "
417  << shift_comp(0)(0, 0, 0, 0) << " "
418  << shift_comp(1)(0, 0, 0, 0) << " "
419  << shift_comp(2)(0, 0, 0, 0) << endl ;
420 
421  ost << " ... w_shift (NB: components in the star Cartesian frame) [c] : "
422  << endl
423  << w_shift(0)(0, 0, 0, 0) << " "
424  << w_shift(1)(0, 0, 0, 0) << " "
425  << w_shift(2)(0, 0, 0, 0) << endl ;
426 
427  ost << "Central value of khi_shift [km c] : "
428  << khi_shift()(0, 0, 0, 0) / km << endl ;
429 
430  ost << endl << "Central value of (B^X, B^Y, B^Z)/N [c] : "
431  << bsn(0)(0, 0, 0, 0) << " "
432  << bsn(1)(0, 0, 0, 0) << " "
433  << bsn(2)(0, 0, 0, 0) << endl ;
434 
435  ost << endl <<
436  "Central (d/dX,d/dY,d/dZ)(logn_auto) [km^{-1}] : "
437  << d_n_auto(0)(0, 0, 0, 0) * km << " "
438  << d_n_auto(1)(0, 0, 0, 0) * km << " "
439  << d_n_auto(2)(0, 0, 0, 0) * km << endl ;
440 
441  ost << endl <<
442  "Central (d/dX,d/dY,d/dZ)(beta_auto) [km^{-1}] : "
443  << d_confpsi_auto(0)(0, 0, 0, 0) * km << " "
444  << d_confpsi_auto(1)(0, 0, 0, 0) * km << " "
445  << d_confpsi_auto(2)(0, 0, 0, 0) * km << endl ;
446 
447  return ost ;
448 }
449 }
Lorene::Et_bin_nsbh::set_n_auto
Tenseur & set_n_auto()
Read/write the lapse {\it N} generated principaly by the star.
Definition: et_bin_nsbh.C:278
Lorene::Etoile::u_euler
Tenseur u_euler
Fluid 3-velocity with respect to the Eulerian observer.
Definition: etoile.h:474
Lorene::Tenseur::set_etat_zero
void set_etat_zero()
Sets the logical state to ETATZERO (zero state).
Definition: tenseur.C:645
Lorene::Map::get_ori_x
double get_ori_x() const
Returns the x coordinate of the origin.
Definition: map.h:768
Lorene::Tenseur
Tensor handling *** DEPRECATED : use class Tensor instead ***.
Definition: tenseur.h:301
Lorene::Etoile_bin::shift_auto
Tenseur shift_auto
Part of the shift vector generated principaly by the star.
Definition: etoile.h:889
Lorene::Et_bin_nsbh::confpsi
Tenseur confpsi
Total conformal factor $\Psi$.
Definition: et_bin_nsbh.h:101
Lorene::Etoile_bin::w_shift
Tenseur w_shift
Vector used in the decomposition of shift_auto , following Shibata's prescription [Prog.
Definition: etoile.h:908
Lorene::Tenseur::set
Cmp & set()
Read/write for a scalar (see also operator=(const Cmp&) ).
Definition: tenseur.C:824
Lorene::Et_bin_nsbh::d_confpsi_auto
Tenseur d_confpsi_auto
Gradient of {\tt confpsi_auto} (Cartesian components with respect to {\tt ref_triad})
Definition: et_bin_nsbh.h:114
Lorene::Etoile_bin::bsn
Tenseur bsn
3-vector shift, divided by N , of the rotating coordinates, .
Definition: etoile.h:950
Lorene::Etoile_bin::operator=
void operator=(const Etoile_bin &)
Assignment to another Etoile_bin.
Definition: etoile_bin.C:482
Lorene::Map::get_rot_phi
double get_rot_phi() const
Returns the angle between the x –axis and X –axis.
Definition: map.h:775
Lorene::Et_bin_nsbh
Class for a star in a NS-BH binary system.
Definition: et_bin_nsbh.h:79
Lorene
Lorene prototypes.
Definition: app_hor.h:64
Lorene::Et_bin_nsbh::taij_auto
Tenseur_sym taij_auto
Part of the extrinsic curvature tensor $\tilde A^{ij} = 2 N K^{ij}$ generated by {\tt shift_auto}.
Definition: et_bin_nsbh.h:126
Lorene::Map::get_mg
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:765
Lorene::Etoile_bin::d_psi
Tenseur d_psi
Gradient of (in the irrotational case) (Cartesian components with respect to ref_triad )
Definition: etoile.h:838
Lorene::Etoile::mp
Map & mp
Mapping associated with the star.
Definition: etoile.h:429
Lorene::Tenseur::sauve
void sauve(FILE *) const
Save in a file.
Definition: tenseur.C:1325
Lorene::Eos
Equation of state base class.
Definition: eos.h:190
Lorene::Map::val_r
virtual double val_r(int l, double xi, double theta, double pphi) const =0
Returns the value of the radial coordinate r for a given in a given domain.
Lorene::max
Tbl max(const Cmp &)
Maximum values of a Cmp in each domain.
Definition: cmp_math.C:435
Lorene::Et_bin_nsbh::ssjm1_confpsi
Cmp ssjm1_confpsi
Effective source at the previous step for the resolution of the Poisson equation for {\tt confpsi_aut...
Definition: et_bin_nsbh.h:164
Lorene::Et_bin_nsbh::set_confpsi_comp
Tenseur & set_confpsi_comp()
Read/write the conformal factor $\Psi$ generated principaly by the companion star.
Definition: et_bin_nsbh.C:299
Lorene::Et_bin_nsbh::taij_tot
Tenseur_sym taij_tot
Total extrinsic curvature tensor $\tilde A^{ij} = 2 N K^{ij}$ generated by {\tt shift_auto} and {\tt ...
Definition: et_bin_nsbh.h:140
Lorene::Et_bin_nsbh::n_auto
Tenseur n_auto
Part of the lapse {\it N} generated principaly by the star.
Definition: et_bin_nsbh.h:85
Lorene::Etoile::ray_eq_pi
double ray_eq_pi() const
Coordinate radius at , [r_unit].
Definition: etoile_global.C:256
Lorene::Etoile_bin
Class for stars in binary system.
Definition: etoile.h:814
Lorene::Et_bin_nsbh::tkij_tot
Tenseur_sym tkij_tot
Total extrinsic curvature tensor $K^{ij}$ generated by {\tt shift_auto} and {\tt shift_comp}.
Definition: et_bin_nsbh.h:152
Unites
Standard units of space, time and mass.
Lorene::Etoile_bin::loggam
Tenseur loggam
Logarithm of the Lorentz factor between the fluid and the co-orbiting observer.
Definition: etoile.h:849
Lorene::Et_bin_nsbh::sauve
virtual void sauve(FILE *) const
Save in a file.
Definition: et_bin_nsbh.C:312
Lorene::Et_bin_nsbh::confpsi_comp
Tenseur confpsi_comp
Part of the conformal factor $\Psi$ generated principaly by the companion star.
Definition: et_bin_nsbh.h:109
Lorene::Et_bin_nsbh::set_n_comp
Tenseur & set_n_comp()
Read/write the lapse {\it N} generated principaly by the companion star.
Definition: et_bin_nsbh.C:285
Lorene::Et_bin_nsbh::ssjm1_lapse
Cmp ssjm1_lapse
Effective source at the previous step for the resolution of the Poisson equation for {\tt n_auto} by ...
Definition: et_bin_nsbh.h:158
Lorene::Et_bin_nsbh::confpsi_auto
Tenseur confpsi_auto
Part of the conformal factor $\Psi$ generated principaly by the star.
Definition: et_bin_nsbh.h:104
Lorene::Cmp
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition: cmp.h:446
Lorene::Etoile_bin::sauve
virtual void sauve(FILE *) const
Save in a file.
Definition: etoile_bin.C:559
Lorene::min
Tbl min(const Cmp &)
Minimum values of a Cmp in each domain.
Definition: cmp_math.C:458
Lorene::Et_bin_nsbh::d_confpsi_comp
Tenseur d_confpsi_comp
Gradient of {\tt confpsi_comp} (Cartesian components with respect to {\tt ref_triad})
Definition: et_bin_nsbh.h:119
Lorene::Et_bin_nsbh::d_n_comp
Tenseur d_n_comp
Gradient of {\tt n_comp} (Cartesian components with respect to {\tt ref_triad})
Definition: et_bin_nsbh.h:98
Lorene::Etoile::shift
Tenseur shift
Total shift vector.
Definition: etoile.h:512
Lorene::Et_bin_nsbh::n_comp
Tenseur n_comp
Part of the lapse {\it N} generated principaly by the companion star.
Definition: et_bin_nsbh.h:88
Lorene::Etoile::operator>>
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition: etoile.C:511
Lorene::Et_bin_nsbh::operator>>
virtual ostream & operator>>(ostream &) const
Save in a file.
Definition: et_bin_nsbh.C:328
Lorene::Tenseur::set_etat_qcq
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: tenseur.C:636
Lorene::Etoile_bin::wit_w
Tenseur wit_w
Spatial projection of the fluid 3-velocity with respect to the co-orbiting observer.
Definition: etoile.h:844
Lorene::Etoile_bin::set_der_0x0
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: etoile_bin.C:459
Lorene::Etoile_bin::xa_barycenter
virtual double xa_barycenter() const
Absolute coordinate X of the barycenter of the baryon density, defined according to the formula.
Definition: et_bin_global.C:170
Lorene::Cmp::sauve
void sauve(FILE *) const
Save in a file.
Definition: cmp.C:561
Lorene::Et_bin_nsbh::set_confpsi_auto
Tenseur & set_confpsi_auto()
Read/write the conformal factor $\Psi$ generated principaly by the star.
Definition: et_bin_nsbh.C:292
Lorene::Etoile::gam_euler
Tenseur gam_euler
Lorentz factor between the fluid and Eulerian observers.
Definition: etoile.h:471
Lorene::Etoile::ray_eq
double ray_eq() const
Coordinate radius at , [r_unit].
Definition: etoile_global.C:120
Lorene::Et_bin_nsbh::operator=
void operator=(const Et_bin_nsbh &)
Destructor.
Definition: et_bin_nsbh.C:253
Lorene::Et_bin_nsbh::d_n_auto
Tenseur d_n_auto
Gradient of {\tt n_auto} (Cartesian components with respect to {\tt ref_triad})
Definition: et_bin_nsbh.h:93
Lorene::Etoile_bin::irrotational
bool irrotational
true for an irrotational star, false for a corotating one
Definition: etoile.h:822
Lorene::Et_bin_nsbh::taij_comp
Tenseur_sym taij_comp
Part of the extrinsic curvature tensor $\tilde A^{ij} = 2 N K^{ij}$ generated by {\tt shift_comp}.
Definition: et_bin_nsbh.h:133
Lorene::Etoile_bin::khi_shift
Tenseur khi_shift
Scalar used in the decomposition of shift_auto , following Shibata's prescription [Prog.
Definition: etoile.h:918
Lorene::Et_bin_nsbh::tkij_auto
Tenseur_sym tkij_auto
Part of the extrinsic curvature tensor $K^{ij}$ generated by {\tt shift_auto}.
Definition: et_bin_nsbh.h:146
Lorene::Cmp::set_etat_qcq
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: cmp.C:304
Lorene::Map
Base class for coordinate mappings.
Definition: map.h:670
Lorene::Etoile_bin::del_deriv
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: etoile_bin.C:447
Lorene::Base_vect
Vectorial bases (triads) with respect to which the tensorial components are defined.
Definition: base_vect.h:105
Lorene::Tenseur::set_std_base
void set_std_base()
Set the standard spectal basis of decomposition for each component.
Definition: tenseur.C:1170
Lorene::Et_bin_nsbh::Et_bin_nsbh
Et_bin_nsbh(Map &mp_i, int nzet_i, bool relat, const Eos &eos_i, bool irrot, const Base_vect &ref_triad_i)
Standard constructor.
Definition: et_bin_nsbh.C:95
Lorene::Etoile_bin::shift_comp
Tenseur shift_comp
Part of the shift vector generated principaly by the companion star.
Definition: etoile.h:895