27 char bin_ns_bh_C[] =
"$Header: /cvsroot/Lorene/C++/Source/Bin_ns_bh/bin_ns_bh.C,v 1.15 2014/10/13 08:52:42 j_novak Exp $" ;
94 #include "bin_ns_bh.h" 95 #include "utilitaires.h" 97 #include "graphique.h" 109 : ref_triad(0.,
"Absolute frame Cartesian basis"),
110 star(mp_ns, nzet, true, eos, irrot_ns, ref_triad),
122 :
ref_triad(0.,
"Absolute frame Cartesian basis"),
135 :
ref_triad(0.,
"Absolute frame Cartesian basis"),
137 hole(mp_bh, fich, old) {
153 Bin_ns_bh::~Bin_ns_bh(){
246 void Bin_ns_bh::sauve(FILE* fich)
const {
257 void Bin_ns_bh::init_auto () {
264 filtre_ns = 0.5 * (1 +
exp(-radius_ns*radius_ns/rlim_ns/rlim_ns)) ;
265 filtre_ns.std_base_scal() ;
271 filtre_bh = 0.5 * (1 +
exp(-radius_bh*radius_bh/rlim_bh/rlim_bh)) ;
272 filtre_bh.std_base_scal() ;
288 Cmp soustrait ((filtre_bh-0.5)*2*
exp(1.)) ;
296 double xa_abs, ya_abs, za_abs ;
297 double air, tet, phi ;
301 for (
int k=0 ; k<np ; k++)
302 for (
int j=0 ; j<nt ; j++) {
304 xa_abs = xa_hole(1,k,j,0) ;
305 ya_abs = ya_hole(1,k,j,0) ;
306 za_abs = za_hole(1,k,j,0) ;
309 for (
int l=1 ; l<nz ; l++)
311 hole.set_n_auto().
set(l,k,j,i) -= (val_star+val_hole)*soustrait(l,k,j,i) ;
321 void Bin_ns_bh::affecte(
const Bin_ns_bh& so) {
348 for (
int k=0 ; k<np ; k++)
349 for (
int j=0 ; j<nt ; j++) {
350 double rcourant = map_et_so->
val_r(
star.
nzet-1, 1, tet(0,k,j,0), phi(0,k,j,0)) ;
351 if (rcourant > rmax) {
358 double old_r = map_et->
val_r(
star.
nzet-1, 1, tet(0,kmax,jmax,0), phi(0,kmax,jmax,0)) ;
370 double precis_secant = 1.e-14 ;
371 double alpha_r = 1. ;
372 double reg_map = 1. ;
375 par_adapt.
add_int(nitermax, 0) ;
377 par_adapt.
add_int(nz_search, 2) ;
378 par_adapt.
add_int(adapt_flag, 3) ;
385 par_adapt.
add_tbl(ent_limit, 0) ;
387 Map_et mp_prev = *map_et ;
391 ent_limit.
set(l) =
star.
ent()(l, kmax, jmax, nr-1) ;
422 for (
int i=0 ; i<3 ; i++)
486 ost <<
"Neutron star - black hole binary system" << endl ;
487 ost <<
"=======================================" << endl ;
489 "Orbital angular velocity : " <<
omega * f_unit <<
" rad/s" << endl ;
491 "Absolute coordinate X of the rotation axis : " <<
x_axe / km
493 ost << endl <<
"Neutron star : " << endl ;
494 ost <<
"============ " << endl ;
495 ost <<
star << endl ;
497 ost <<
"Black hole : " << endl ;
498 ost <<
"========== " << endl ;
499 ost <<
"Coordinate radius of the throat : " <<
hole.
get_rayon() / km <<
" km" << endl ;
500 ost <<
"Absolute abscidia of the throat center : " << (
hole.
get_mp()).get_ori_x() / km
Coord xa
Absolute x coordinate.
void set_omega(double)
Sets the orbital angular velocity [{ f_unit}].
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Cmp exp(const Cmp &)
Exponential.
void set_der_0x0() const
Sets to { 0x0} all the pointers on derived quantities.
const Map & get_mp() const
Returns the mapping.
Radial mapping of rather general form.
void add_int(const int &n, int position=0)
Adds the address of a new int to the list.
Neutron star - black hole binary system.
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
void fait_shift_auto()
Computes shift_auto from w_shift and khi_shift according to Shibata's prescription [Prog...
Cmp sqrt(const Cmp &)
Square root.
Bin_ns_bh(Map &mp_ns, int nzet, const Eos &eos, bool irrot_ns, Map_af &mp_bh)
Standard constructor.
void set_rot_phi(double phi0)
Sets a new rotation angle.
void set_std_base()
Set the standard spectal basis of decomposition for each component.
double * p_total_ener
Total energy of the system.
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.
virtual void kinematics(double omega, double x_axe)
Computes the quantities bsn and pot_centri .
virtual double val_r(int l, double xi, double theta, double pphi) const
Returns the value of the radial coordinate r for a given in a given domain.
double & set(int i)
Read/write of a particular element (index i) (1D case)
void set_rayon(double ray)
Sets the radius of the horizon to ray .
Base class for coordinate mappings.
double get_ori_x() const
Returns the x coordinate of the origin.
virtual double val_r(int l, double xi, double theta, double pphi) const
Returns the value of the radial coordinate r for a given in a given domain.
double * p_ham_constr
Relative error on the Hamiltonian constraint.
const Tenseur & get_n_auto() const
Returns the part of N generated by the hole.
const Tenseur & get_logn_auto() const
Returns the logarithm of the part of the lapse N generated principaly by the star.
void update_metric_der_comp(const Bhole &comp)
Computes the derivative of metric functions related to the companion black hole.
bool is_irrotational() const
Returns true for an irrotational motion, false for a corotating one.
Tenseur psi_auto
Part of generated by the hole.
Tenseur & set_n_auto()
Read/write the lapse { N} generated principaly by the star.
double * p_mass_adm
Total ADM mass of the system.
virtual void sauve(FILE *) const
Save in a file.
double get_rot_phi() const
Returns the angle between the x –axis and X –axis.
Tenseur n_auto
Part of the lapse { N} generated principaly by the star.
Tenseur & set_confpsi_auto()
Read/write the conformal factor $$ generated principaly by the star.
Coord tet
coordinate centered on the grid
Tenseur shift_auto
Part of generated by the hole.
void set_ori(double xa0, double ya0, double za0)
Sets a new origin.
const Tenseur & get_psi_auto() const
Returns the part of generated by the hole.
void del_deriv() const
Destructor.
Coord phi
coordinate centered on the grid
const Base_vect_cart ref_triad
Cartesian triad of the absolute reference frame.
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
double * p_virial_fus
Virial theorem error by J.L.Friedman, K.Uryu, and M.Shibata.
double * p_virial
Virial theorem error.
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...
double * p_mass_kom
Total Komar mass of the system.
friend ostream & operator<<(ostream &, const Bin_ns_bh &)
Save in a file.
Map & set_mp()
Read/write of the mapping.
double x_axe
Absolute X coordinate of the rotation axis.
Et_bin_nsbh star
The neutron star.
Tbl * p_mom_constr
Relative error on the momentum constraint.
Map_af & mp
Affine mapping.
const Tenseur & get_n_auto() const
Returns the part of the lapse { N} generated principaly by the star.
void add_tbl(const Tbl &ti, int position=0)
Adds the address of a new Tbl to the list.
Map & mp
Mapping associated with the star.
int get_nzone() const
Returns the number of domains.
virtual void equation_of_state()
Computes the proper baryon and energy density, as well as pressure from the enthalpy.
virtual void homothetie(double lambda)
Sets a new radial scale.
ostream & operator>>(ostream &) const
Operator >> (function called by the operator <<).
virtual void reevaluate_symy(const Map *mp_prev, int nzet, Cmp &uu) const
Recomputes the values of a Cmp at the collocation points after a change in the mapping.
double get_rayon() const
Returns the radius of the horizon.
Tenseur khi_shift
Scalar used in the decomposition of shift_auto , following Shibata's prescription [Prog...
void dec2_dzpuis()
Decreases by 2 the value of dzpuis and changes accordingly the values of the Cmp in the external comp...
double omega
Angular velocity with respect to an asymptotically inertial observer.
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.
void sauve(FILE *fich) const
Write on a file.
double rayon
Radius of the horizon in LORENE's units.
void update_metric(const Bhole &comp)
Computes metric coefficients from known potentials, when the companion is a black hole...
const Map_af & get_mp() const
Returns the mapping (readonly).
void std_base_scal()
Sets the spectral bases of the Valeur va to the standard ones for a scalar.
int nzet
Number of domains of *mp occupied by the star.
int fread_be(int *aa, int size, int nb, FILE *fich)
Reads integer(s) from a binary file according to the big endian convention.
int get_nr(int l) const
Returns the number of points in the radial direction ( ) in domain no. l.
Coord ya
Absolute y coordinate.
virtual void adapt(const Cmp &ent, const Param &par, int nbr_filtre=0)
Adaptation of the mapping to a given scalar field.
Tbl & set(int l)
Read/write of the value in a given domain.
double get_omega() const
Returns the angular velocity.
Tenseur confpsi_auto
Part of the conformal factor $$ generated principaly by the star.
void allocate_all()
Sets the logical state to ETATQCQ (ordinary state) and performs the memory allocation of all the elem...
Map_af & set_mp()
Read/write of the mapping.
void operator=(const Bin_ns_bh &)
Assignment to another Bin_ns_bh.
Tenseur ent
Log-enthalpy (relativistic case) or specific enthalpy (Newtonian case)
void raccord(int n)
Performs the matching of the nucleus with respect to the first shell.
Tbl * p_angu_mom
Total angular momentum of the system.
Coord za
Absolute z coordinate.
void set_omega(double ome)
Sets the angular velocity to ome .
void add_double(const double &x, int position=0)
Adds the the address of a new double to the list.
void fait_tkij(int bound_nn=-1, double lim_nn=0)
Computation of the extrinsic curvature tensor for both { star} and { bhole}.
double * p_virial_gb
Virial theorem error by E.Gourgoulhon and S.Bonazzola.
Bhole hole
The black hole.
Tenseur n_auto
Part of N generated by the hole.
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
void set_x_axe(double)
Sets the absolute coordinate X of the rotation axis [{ r_unit}].
void homothetie_interne(double lambda)
Sets a new radial scale at the bondary between the nucleus and the first shell.
void import(const Cmp &ci)
Assignment to another Cmp defined on a different mapping.
double separation() const
Return the separation.
void convert_absolute(double xx, double yy, double zz, double &rr, double &theta, double &pphi) const
Determines the coordinates corresponding to given absolute Cartesian coordinates (X...
Tensor handling *** DEPRECATED : use class Tensor instead ***.
void add_int_mod(int &n, int position=0)
Adds the address of a new modifiable int to the list.
const Tenseur & get_beta_auto() const
Returns the logarithm of the part of the product AN generated principaly by the star.
Coord r
r coordinate centered on the grid
Tenseur d_psi
Gradient of (in the irrotational case) (Cartesian components with respect to ref_triad ) ...
virtual void hydro_euler()
Computes the hydrodynamical quantities relative to the Eulerian observer from those in the fluid fram...