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 $" ;
81 #include "et_bin_nsbh.h" 84 #include "utilitaires.h" 97 :
Etoile_bin(mp_i, nzet_i, relat, eos_i, irrot, ref_triad_i),
100 d_n_auto(mp_i, 1, COV, ref_triad_i),
101 d_n_comp(mp_i, 1, COV, ref_triad_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),
113 ssjm1_confpsi(mp_i) {
174 const Base_vect& ref_triad_i, FILE* fich,
bool old)
178 d_n_auto(mp_i, 1, COV, ref_triad_i),
179 d_n_comp(mp_i, 1, COV, ref_triad_i ),
187 taij_tot(mp_i, 2, CON, ref_triad_i),
189 tkij_tot(mp_i, 2, CON, ref_triad_i),
194 Cmp n_from_file (mp_i, *(mp_i.
get_mg()), fich) ;
198 Cmp psi_from_file (mp_i, *(mp_i.
get_mg()), fich) ;
203 Tenseur shift_auto_file (
mp, ref_triad_i, fich) ;
226 Cmp ssjm1_lapse_file(mp_i, *(mp_i.
get_mg()), fich) ;
229 Cmp ssjm1_confpsi_file(mp_i, *(mp_i.
get_mg()), fich) ;
241 Et_bin_nsbh::~Et_bin_nsbh(){
335 ost <<
"Neutron star in a binary system" << endl ;
336 ost <<
"-------------------------------" << endl ;
339 ost <<
"irrotational configuration" << endl ;
342 ost <<
"corotating configuration" << endl ;
345 ost <<
"Absolute abscidia of the stellar center: " <<
348 ost <<
"Absolute abscidia of the barycenter of the baryon density : " <<
353 double d_tilde = 2 * d_ns / r_0 ;
355 ost <<
"d_tilde : " << d_tilde << endl ;
357 ost <<
"Orientation with respect to the absolute frame : " <<
360 ost <<
"Central value of gam_euler : " 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 ;
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 ;
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 ;
379 ost <<
"Max vel. / co-orb. (W^X, W^Y, W^Z) [c] : " 384 ost <<
"Min vel. / co-orb. (W^X, W^Y, W^Z) [c] : " 389 double r_surf =
mp.
val_r(0,1.,M_PI/4,M_PI/4) ;
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 ;
396 ost <<
"Central value of loggam : " 397 <<
loggam()(0, 0, 0, 0) << endl ;
400 ost <<
"Central value of lapse(N) auto : " 401 <<
n_auto()(0, 0, 0, 0) << endl ;
403 ost <<
"Central value of confpsi auto : " 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 ;
411 ost <<
" ... shift_auto part of it [c] : " 416 ost <<
" ... shift_comp part of it [c] : " 421 ost <<
" ... w_shift (NB: components in the star Cartesian frame) [c] : " 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 ;
427 ost <<
"Central value of khi_shift [km c] : " 428 <<
khi_shift()(0, 0, 0, 0) / km << endl ;
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 ;
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 ;
442 "Central (d/dX,d/dY,d/dZ)(beta_auto) [km^{-1}] : " Tenseur & set_confpsi_comp()
Read/write the conformal factor $$ generated principaly by the companion star.
virtual ostream & operator>>(ostream &) const
Save in a file.
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Tenseur shift_comp
Part of the shift vector generated principaly by the companion star.
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
virtual void sauve(FILE *) const
Save in a file.
void set_std_base()
Set the standard spectal basis of decomposition for each component.
Standard units of space, time and mass.
Equation of state base class.
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
double ray_eq() const
Coordinate radius at , [r_unit].
Base class for coordinate mappings.
double get_ori_x() const
Returns the x coordinate of the origin.
Tenseur confpsi_comp
Part of the conformal factor $$ generated principaly by the companion star.
Tenseur wit_w
Spatial projection of the fluid 3-velocity with respect to the co-orbiting observer.
Tenseur_sym taij_comp
Part of the extrinsic curvature tensor $ A^{ij} = 2 N K^{ij}$ generated by { shift_comp}.
Tenseur & set_n_comp()
Read/write the lapse { N} generated principaly by the companion star.
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.
Tenseur & set_n_auto()
Read/write the lapse { N} generated principaly by the star.
Class for stars in binary system.
virtual double xa_barycenter() const
Absolute coordinate X of the barycenter of the baryon density, defined according to the formula wher...
virtual void sauve(FILE *) const
Save in a file.
double get_rot_phi() const
Returns the angle between the x –axis and X –axis.
Tbl min(const Cmp &)
Minimum values of a Cmp in each domain.
Tenseur n_auto
Part of the lapse { N} generated principaly by the star.
Vectorial bases (triads) with respect to which the tensorial components are defined.
Tenseur shift
Total shift vector.
Tenseur confpsi
Total conformal factor $$.
Tenseur & set_confpsi_auto()
Read/write the conformal factor $$ generated principaly by the star.
Tenseur d_confpsi_comp
Gradient of { confpsi_comp} (Cartesian components with respect to { ref_triad})
Tenseur shift_auto
Part of the shift vector generated principaly by the star.
bool irrotational
true for an irrotational star, false for a corotating one
Tenseur_sym taij_tot
Total extrinsic curvature tensor $ A^{ij} = 2 N K^{ij}$ generated by { shift_auto} and { shift_comp}...
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.
Cmp & set()
Read/write for a scalar (see also operator=(const Cmp&) ).
Tenseur w_shift
Vector used in the decomposition of shift_auto , following Shibata's prescription [Prog...
Tenseur u_euler
Fluid 3-velocity with respect to the Eulerian observer.
virtual void del_deriv() const
Deletes all the derived quantities.
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
void sauve(FILE *) const
Save in a file.
Tenseur bsn
3-vector shift, divided by N , of the rotating coordinates, .
Tenseur gam_euler
Lorentz factor between the fluid and Eulerian observers.
Etoile_bin(Map &mp_i, int nzet_i, bool relat, const Eos &eos_i, bool irrot, const Base_vect &ref_triad_i)
Standard constructor.
Map & mp
Mapping associated with the star.
Tbl max(const Cmp &)
Maximum values of a Cmp in each domain.
Tenseur khi_shift
Scalar used in the decomposition of shift_auto , following Shibata's prescription [Prog...
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Tenseur d_confpsi_auto
Gradient of { confpsi_auto} (Cartesian components with respect to { ref_triad})
Tenseur d_n_auto
Gradient of { n_auto} (Cartesian components with respect to { ref_triad})
void operator=(const Etoile_bin &)
Assignment to another Etoile_bin.
Tenseur confpsi_auto
Part of the conformal factor $$ generated principaly by the star.
Class for a star in a NS-BH binary system.
double ray_eq_pi() const
Coordinate radius at , [r_unit].
Cmp ssjm1_confpsi
Effective source at the previous step for the resolution of the Poisson equation for { confpsi_auto} ...
Tenseur_sym taij_auto
Part of the extrinsic curvature tensor $ A^{ij} = 2 N K^{ij}$ generated by { shift_auto}.
Tenseur loggam
Logarithm of the Lorentz factor between the fluid and the co-orbiting observer.
Cmp ssjm1_lapse
Effective source at the previous step for the resolution of the Poisson equation for { n_auto} by mea...
Tenseur_sym tkij_tot
Total extrinsic curvature tensor $K^{ij}$ generated by { shift_auto} and { shift_comp}.
void operator=(const Et_bin_nsbh &)
Destructor.
Tenseur_sym tkij_auto
Part of the extrinsic curvature tensor $K^{ij}$ generated by { shift_auto}.
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
void sauve(FILE *) const
Save in a file.
void set_etat_zero()
Sets the logical state to ETATZERO (zero state).
Tenseur n_comp
Part of the lapse { N} generated principaly by the companion star.
Tenseur d_n_comp
Gradient of { n_comp} (Cartesian components with respect to { ref_triad})
Tensor handling *** DEPRECATED : use class Tensor instead ***.
Tenseur d_psi
Gradient of (in the irrotational case) (Cartesian components with respect to ref_triad ) ...