84 #include "et_bin_nsbh.h" 87 #include "utilitaires.h" 100 :
Etoile_bin(mp_i, nzet_i, relat, eos_i, irrot, ref_triad_i),
103 d_n_auto(mp_i, 1, COV, ref_triad_i),
104 d_n_comp(mp_i, 1, COV, ref_triad_i),
108 d_confpsi_auto(mp_i, 1, COV, ref_triad_i),
109 d_confpsi_comp(mp_i, 1, COV, ref_triad_i),
110 taij_auto(mp_i, 2, CON, ref_triad_i) ,
111 taij_comp(mp_i, 2, CON, ref_triad_i),
112 taij_tot(mp_i, 2, CON, ref_triad_i) ,
113 tkij_auto(mp_i, 2, CON, ref_triad_i),
114 tkij_tot(mp_i, 2, CON, ref_triad_i),
116 ssjm1_confpsi(mp_i) {
156 d_n_auto(et.d_n_auto),
157 d_n_comp(et.d_n_comp),
159 confpsi_auto(et.confpsi_auto),
160 confpsi_comp(et.confpsi_comp),
161 d_confpsi_auto(et.d_confpsi_auto),
162 d_confpsi_comp(et.d_confpsi_comp),
163 taij_auto(et.taij_auto),
164 taij_comp(et.taij_comp),
165 taij_tot(et.taij_tot),
166 tkij_auto(et.tkij_auto),
167 tkij_tot(et.tkij_tot),
168 ssjm1_lapse(et.ssjm1_lapse),
169 ssjm1_confpsi(et.ssjm1_confpsi) {
177 const Base_vect& ref_triad_i, FILE* fich,
bool old)
181 d_n_auto(mp_i, 1, COV, ref_triad_i),
182 d_n_comp(mp_i, 1, COV, ref_triad_i ),
186 d_confpsi_auto(mp_i, 1, COV, ref_triad_i),
187 d_confpsi_comp(mp_i, 1, COV, ref_triad_i),
188 taij_auto(mp_i, 2, CON, ref_triad_i),
189 taij_comp(mp_i, 2, CON, ref_triad_i),
190 taij_tot(mp_i, 2, CON, ref_triad_i),
191 tkij_auto(mp_i, 2, CON, ref_triad_i),
192 tkij_tot(mp_i, 2, CON, ref_triad_i),
194 ssjm1_confpsi(mp_i) {
197 Cmp n_from_file (mp_i, *(mp_i.
get_mg()), fich) ;
201 Cmp psi_from_file (mp_i, *(mp_i.
get_mg()), fich) ;
206 Tenseur shift_auto_file (
mp, ref_triad_i, fich) ;
229 Cmp ssjm1_lapse_file(mp_i, *(mp_i.
get_mg()), fich) ;
232 Cmp ssjm1_confpsi_file(mp_i, *(mp_i.
get_mg()), fich) ;
244 Et_bin_nsbh::~Et_bin_nsbh(){
338 ost <<
"Neutron star in a binary system" << endl ;
339 ost <<
"-------------------------------" << endl ;
342 ost <<
"irrotational configuration" << endl ;
345 ost <<
"corotating configuration" << endl ;
348 ost <<
"Absolute abscidia of the stellar center: " <<
351 ost <<
"Absolute abscidia of the barycenter of the baryon density : " <<
356 double d_tilde = 2 * d_ns / r_0 ;
358 ost <<
"d_tilde : " << d_tilde << endl ;
360 ost <<
"Orientation with respect to the absolute frame : " <<
363 ost <<
"Central value of gam_euler : " 366 ost <<
"Central u_euler (U^X, U^Y, U^Z) [c] : " 367 <<
u_euler(0)(0, 0, 0, 0) <<
" " 368 <<
u_euler(1)(0, 0, 0, 0) <<
" " 369 <<
u_euler(2)(0, 0, 0, 0) << endl ;
372 ost <<
"Central d_psi (X, Y, Z) [c] : " 373 <<
d_psi(0)(0, 0, 0, 0) <<
" " 374 <<
d_psi(1)(0, 0, 0, 0) <<
" " 375 <<
d_psi(2)(0, 0, 0, 0) << endl ;
377 ost <<
"Central vel. / co-orb. (W^X, W^Y, W^Z) [c] : " 378 <<
wit_w(0)(0, 0, 0, 0) <<
" " 379 <<
wit_w(1)(0, 0, 0, 0) <<
" " 380 <<
wit_w(2)(0, 0, 0, 0) << endl ;
382 ost <<
"Max vel. / co-orb. (W^X, W^Y, W^Z) [c] : " 387 ost <<
"Min vel. / co-orb. (W^X, W^Y, W^Z) [c] : " 392 double r_surf =
mp.
val_r(0,1.,M_PI/4,M_PI/4) ;
394 ost <<
"Velocity at (r_surf,pi/4,pi/4) / co-orb. [c] : " 395 <<
wit_w(0).val_point(r_surf,M_PI/4,M_PI/4) <<
" " 396 <<
wit_w(1).val_point(r_surf,M_PI/4,M_PI/4) <<
" " 397 <<
wit_w(2).val_point(r_surf,M_PI/4,M_PI/4) << endl ;
399 ost <<
"Central value of loggam : " 400 <<
loggam()(0, 0, 0, 0) << endl ;
403 ost <<
"Central value of lapse(N) auto : " 404 <<
n_auto()(0, 0, 0, 0) << endl ;
406 ost <<
"Central value of confpsi auto : " 409 ost <<
"Central value of shift (N^X, N^Y, N^Z) [c] : " 410 <<
shift(0)(0, 0, 0, 0) <<
" " 411 <<
shift(1)(0, 0, 0, 0) <<
" " 412 <<
shift(2)(0, 0, 0, 0) << endl ;
414 ost <<
" ... shift_auto part of it [c] : " 419 ost <<
" ... shift_comp part of it [c] : " 424 ost <<
" ... w_shift (NB: components in the star Cartesian frame) [c] : " 426 <<
w_shift(0)(0, 0, 0, 0) <<
" " 427 <<
w_shift(1)(0, 0, 0, 0) <<
" " 428 <<
w_shift(2)(0, 0, 0, 0) << endl ;
430 ost <<
"Central value of khi_shift [km c] : " 431 <<
khi_shift()(0, 0, 0, 0) / km << endl ;
433 ost << endl <<
"Central value of (B^X, B^Y, B^Z)/N [c] : " 434 <<
bsn(0)(0, 0, 0, 0) <<
" " 435 <<
bsn(1)(0, 0, 0, 0) <<
" " 436 <<
bsn(2)(0, 0, 0, 0) << endl ;
439 "Central (d/dX,d/dY,d/dZ)(logn_auto) [km^{-1}] : " 440 <<
d_n_auto(0)(0, 0, 0, 0) * km <<
" " 441 <<
d_n_auto(1)(0, 0, 0, 0) * km <<
" " 442 <<
d_n_auto(2)(0, 0, 0, 0) * km << endl ;
445 "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.
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 ) ...