LORENE
bin_bhns_extr.C
1/*
2 * Methods of class Bin_bhns_extr
3 *
4 * (see file bin_bhns_extr.h for documentation).
5 *
6 */
7
8/*
9 * Copyright (c) 2004-2005 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
28char bin_bhns_extr_C[] = "$Header: /cvsroot/Lorene/C++/Source/Bin_bhns_extr/bin_bhns_extr.C,v 1.10 2014/10/13 08:52:41 j_novak Exp $" ;
29
30/*
31 * $Id: bin_bhns_extr.C,v 1.10 2014/10/13 08:52:41 j_novak Exp $
32 * $Log: bin_bhns_extr.C,v $
33 * Revision 1.10 2014/10/13 08:52:41 j_novak
34 * Lorene classes and functions now belong to the namespace Lorene.
35 *
36 * Revision 1.9 2014/10/06 15:13:00 j_novak
37 * Modified #include directives to use c++ syntax.
38 *
39 * Revision 1.8 2005/02/28 23:06:16 k_taniguchi
40 * Modification of some arguments to include the case of the conformally
41 * flat background metric.
42 *
43 * Revision 1.7 2004/12/13 21:08:59 k_taniguchi
44 * Addition of some outputs in display_poly.
45 *
46 * Revision 1.6 2004/12/03 20:17:52 k_taniguchi
47 * Correction of "display_poly">
48 *
49 * Revision 1.5 2004/12/02 22:43:35 k_taniguchi
50 * Modification of "display_poly".
51 *
52 * Revision 1.4 2004/12/02 17:36:42 k_taniguchi
53 * Modification of "display_poly()".
54 *
55 * Revision 1.3 2004/12/01 16:37:11 k_taniguchi
56 * Modification of the output again.
57 *
58 * Revision 1.2 2004/12/01 16:27:09 k_taniguchi
59 * Modification of the output.
60 *
61 * Revision 1.1 2004/11/30 20:45:49 k_taniguchi
62 * *** empty log message ***
63 *
64 *
65 * $Header: /cvsroot/Lorene/C++/Source/Bin_bhns_extr/bin_bhns_extr.C,v 1.10 2014/10/13 08:52:41 j_novak Exp $
66 *
67 */
68
69// C headers
70#include <cmath>
71
72// Lorene headers
73#include "bin_bhns_extr.h"
74#include "eos.h"
75#include "utilitaires.h"
76#include "unites.h"
77
78 //--------------------------------//
79 // Constructors //
80 //--------------------------------//
81
82// Standard constructor
83// --------------------
84namespace Lorene {
85Bin_bhns_extr::Bin_bhns_extr(Map& mp, int nzet, const Eos& eos, bool irrot,
86 bool relat, bool kerrs, bool multi)
87 : ref_triad(0., "Absolute frame Cartesian basis"),
88 star(mp, nzet, relat, eos, irrot, ref_triad, kerrs, multi)
89{
90
91 omega = 0. ;
92 separ = 0. ;
93 mass_bh = 0. ;
94
95 // Pointers of derived quantities initialized to zero :
96 set_der_0x0() ;
97
98}
99
100// Copy constructor
101// ----------------
103 : ref_triad(0., "Absolute frame Cartesian basis"),
104 star(bibi.star),
105 omega(bibi.omega),
106 separ(bibi.separ),
107 mass_bh(bibi.mass_bh)
108{
109
110 // Pointers of derived quantities initialized to zero :
111 set_der_0x0() ;
112
113}
114
115// Constructor from a file
116// -----------------------
117Bin_bhns_extr::Bin_bhns_extr(Map& mp, const Eos& eos, FILE* fich)
118 : ref_triad(0., "Absolute frame Cartesian basis"),
119 star(mp, eos, ref_triad, fich)
120{
121
122 // omega and x_axe are read in the file :
123 fread_be(&omega, sizeof(double), 1, fich) ;
124 fread_be(&separ, sizeof(double), 1, fich) ;
125 fread_be(&mass_bh, sizeof(double), 1, fich) ;
126
127 // Pointers of derived quantities initialized to zero :
128 set_der_0x0() ;
129
130}
131
132 //------------------------------//
133 // Destructor //
134 //------------------------------//
135
140
141 //----------------------------------------------------//
142 // Management of derived quantities //
143 //----------------------------------------------------//
144
146
147 if (p_xa_barycenter_extr != 0x0) delete p_xa_barycenter_extr ;
148 if (p_ya_barycenter_extr != 0x0) delete p_ya_barycenter_extr ;
149 if (p_mass_b_extr != 0x0) delete p_mass_b_extr ;
150
151 set_der_0x0() ;
152
153}
154
156
159 p_mass_b_extr = 0x0 ;
160
161}
162
163 //------------------------------//
164 // Assignment //
165 //------------------------------//
166
167// Assignment to another Bin_bhns_extr
168// -----------------------------------
170
171 assert( bibi.ref_triad == ref_triad ) ;
172
173 star = bibi.star ;
174
175 omega = bibi.omega ;
176 separ = bibi.separ ;
177 mass_bh = bibi.mass_bh ;
178
179 // ref_triad remains unchanged
180
181 del_deriv() ; // Deletes all derived quantities
182
183}
184
185 //---------------------------//
186 // Outputs //
187 //---------------------------//
188
189// Save in a file
190// --------------
191void Bin_bhns_extr::sauve(FILE* fich) const {
192
193 star.sauve(fich) ;
194
195 fwrite_be(&omega, sizeof(double), 1, fich) ;
196 fwrite_be(&separ, sizeof(double), 1, fich) ;
197 fwrite_be(&mass_bh, sizeof(double), 1, fich) ;
198
199}
200
201// Printing
202// --------
203ostream& operator<<(ostream& ost, const Bin_bhns_extr& bibi) {
204
205 bibi >> ost ;
206 return ost ;
207
208}
209
210ostream& Bin_bhns_extr::operator>>(ostream& ost) const {
211
212 using namespace Unites ;
213
214 ost << endl ;
215 ost << "Binary BH-NS system" << endl ;
216 ost << "===================" << endl ;
217 ost << endl ;
218 if (star.in_kerrschild()) {
219 ost << "Kerr-Schild background metric" << endl ;
220 ost << "-----------------------------" << endl ;
221 }
222 else {
223 ost << "Conformally flat background metric" << endl ;
224 ost << "----------------------------------" << endl ;
225 }
226 if (star.with_multipole()) {
227 ost << "Multipole falloff boundary condition" << endl ;
228 ost << "------------------------------------" << endl ;
229 }
230 else {
231 ost << "1/r falloff boundary condition" << endl ;
232 ost << "------------------------------" << endl ;
233 }
234 ost << endl
235 << "Orbital angular velocity : "
236 << omega * f_unit << " rad/s" << endl ;
237 ost << "Coordinate separation between BH and NS : "
238 << separ / km << " km" << endl ;
239
240 ost << endl << "Neutron star : " ;
241 ost << endl << "============ " << endl ;
242 ost << endl ;
243 if (star.is_relativistic()) {
244 ost << "Relativistic star" << endl ;
245 ost << "-----------------" << endl ;
246 }
247 else {
248 ost << "WARNING : BH-NS binary should be relativistic !!!" << endl ;
249 ost << "-------------------------------------------------" << endl ;
250 }
251 ost << "Number of domains occupied by the star : " << star.get_nzet()
252 << endl ;
253 ost << "Equation of state : " << endl ;
254 ost << star.get_eos() << endl ;
255
256 ost << endl
257 << "Enthalpy at the coordinate origin : "
258 << star.get_ent()()(0,0,0,0) << " c^2" << endl ;
259 ost << "Proper baryon density at the coordinate origin : "
260 << star.get_nbar()()(0,0,0,0) << " x 0.1 fm^-3" << endl ;
261 ost << "Proper energy density at the coordinate origin : "
262 << star.get_ener()()(0,0,0,0) << " rho_nuc c^2" << endl ;
263 ost << "Pressure at the coordinate origin : "
264 << star.get_press()()(0,0,0,0) << " rho_nuc c^2" << endl ;
265 ost << endl
266 << "Lapse N at the coordinate origin : "
267 << star.get_nnn()()(0,0,0,0) << endl ;
268 ost << "Conformal factor A^2 at the coordinate origin : "
269 << star.get_a_car()()(0,0,0,0) << endl ;
270
271 ost << endl
272 << "Equatorial radius (to BH) a_to : "
273 << star.ray_eq_pi()/km << " km" << endl ;
274 ost << "Equatorial radius (opp. to BH) a_opp : "
275 << star.ray_eq()/km << " km" << endl ;
276
277 ost << endl
278 << "Baryon mass in isolation : " << star.mass_b() / msol
279 << " Mo" << endl
280 << "Gravitational mass in isolation : " << star.mass_g() / msol
281 << " Mo" << endl
282 << "Baryon mass in a binary system : " << mass_b_extr() / msol
283 << " Mo" << endl ;
284
285 ost << endl ;
286 ost << "Star in a binary system" << endl ;
287 ost << "-----------------------" << endl ;
288
289 if (star.is_irrotational()) {
290 ost << "irrotational configuration" << endl ;
291 }
292 else {
293 ost << "corotating configuration" << endl ;
294 }
295
296 ost << "Absolute abscidia of the stellar center: "
297 << star.get_mp().get_ori_x()/km << " km" << endl ;
298 ost << "Orientation with respect to the absolute frame : "
299 << star.get_mp().get_rot_phi() << " rad" << endl ;
300
301 double r0 = 0.5 * ( star.ray_eq() + star.ray_eq_pi() ) ;
302 // = 1 in Baumgarte et al.
303 double d_ns = separ + star.ray_eq() - r0 ;
304 // Orbital separation of Baumgarte et al.
305 double d_tilde = d_ns / r0 ; // Normalized separation of Baumgarte et al.
306
307 ost << endl << "Comparison with those by Baumgarte et al. :" << endl ;
308 ost << " Radius r0 : " << r0 / km << " km " << endl ;
309 ost << " Separation d : " << d_ns / km << " km " << endl ;
310 ost << " Normalized sep. (d/r0) : " << d_tilde << endl ;
311
312 ost << endl << "Black hole : " ;
313 ost << endl << "========== " << endl ;
314 ost << "Gravitational mass of BH : "
315 << mass_bh / msol << " M_sol" << endl ;
316
317 return ost ;
318
319}
320
321// Display in polytropic units
322// ---------------------------
323void Bin_bhns_extr::display_poly(ostream& ost) const {
324
325 using namespace Unites ;
326
327 const Eos* p_eos = &( star.get_eos() ) ;
328 const Eos_poly* p_eos_poly = dynamic_cast<const Eos_poly*>( p_eos ) ;
329
330 if (p_eos_poly != 0x0) {
331
332 double kappa = p_eos_poly->get_kap() ;
333 double gamma = p_eos_poly->get_gam() ;
334 double kap_ns2 = pow( kappa, 0.5 /(gamma-1) ) ;
335
336 // Polytropic unit of length in terms of r_unit :
337 double r_poly = kap_ns2 / sqrt(ggrav) ;
338
339 // Polytropic unit of time in terms of t_unit :
340 double t_poly = r_poly ;
341
342 // Polytropic unit of mass in terms of m_unit :
343 double m_poly = r_poly / ggrav ;
344
345 // Polytropic unit of angular momentum in terms of j_unit :
346 // double j_poly = r_poly * r_poly / ggrav ;
347
348 double r0 = 0.5 * ( star.ray_eq() + star.ray_eq_pi() ) ;
349 // = 1 in Baumgarte et al.
350 double d_ns = separ + star.ray_eq() - r0 ;
351 // Orbital separation of Baumgarte et al.
352
353 ost.precision(16) ;
354 ost << endl << "Quantities in polytropic units : " ;
355 ost << endl << "==============================" << endl ;
356 ost << " ( r_poly = " << r_poly / km << " km )" << endl ;
357 ost << " d_e_max : " << separ / r_poly << endl ;
358 ost << " d_G_x : " << xa_barycenter_extr() / r_poly << endl
359 << " d_G_y : " << ya_barycenter_extr() / r_poly << endl ;
360 ost << " d_bh/M_bh : " << d_ns/r_poly / (mass_bh/m_poly) << endl ;
361 ost << " Omega : " << omega * t_poly << endl ;
362 ost << " Omega M_bh : " << omega * t_poly * mass_bh / m_poly
363 << endl ;
364 ost << " M_bar(NS) : " << mass_b_extr() / m_poly << endl ;
365 ost << " M_bar(NS_0) : " << star.mass_b() / m_poly << endl ;
366 ost << " R_0(NS) : "
367 << 0.5 * (star.ray_eq() + star.ray_eq_pi()) / r_poly << endl ;
368 ost << " M_grv(BH) : " << mass_bh / m_poly << endl ;
369
370 }
371
372}
373}
Class for computing a Black hole - Neutron star binary system with an extreme mass ratio.
const Base_vect_cart ref_triad
Cartesian triad of the absolute reference frame.
double * p_ya_barycenter_extr
Absolute coordinate Y of the barycenter of the baryon density in the Kerr-Schild background metric.
double * p_mass_b_extr
Baryon mass of the neutron star in the Kerr-Schild background metric or in the conformally flat one.
double ya_barycenter_extr() const
in the Kerr-Schild background metric
double mass_b_extr() const
Baryon mass of the neutron star in the Kerr-Schild background metric or in the conformally flat.
Et_bin_bhns_extr star
Neutron star.
void operator=(const Bin_bhns_extr &)
Assignment to another Bin_bhns_extr.
Bin_bhns_extr(Map &mp, int nzet, const Eos &eos, bool irrot, bool relat, bool kerrs, bool multi)
Standard constructor.
void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
double separ
Absolute orbital separation between two centers of BH and NS.
double mass_bh
Gravitational mass of BH.
double xa_barycenter_extr() const
Absolute coordinate X of the barycenter of the baryon density in the Kerr-Schild background metric or...
void display_poly(ostream &) const
Display in polytropic units.
void sauve(FILE *) const
Save in a file.
ostream & operator>>(ostream &) const
Operator >> (function called by the operator <<)
~Bin_bhns_extr()
Destructor.
double omega
Angular velocity with respect to an asymptotically inertial observer.
double * p_xa_barycenter_extr
Absolute coordinate X of the barycenter of the baryon density in the Kerr-Schild background metric or...
void del_deriv() const
Deletes all the derived quantities.
Polytropic equation of state (relativistic case).
Definition eos.h:757
double get_gam() const
Returns the adiabatic index (cf. Eq. (3))
Definition eos_poly.C:256
double get_kap() const
Returns the pressure coefficient (cf.
Definition eos_poly.C:260
Equation of state base class.
Definition eos.h:190
bool in_kerrschild() const
Returns true for the Kerr-Schild background metric, false for the conformally flat one.
virtual void sauve(FILE *) const
Save in a file.
bool with_multipole() const
Returns true for the multipole falloff boundary condition, false for the one.
virtual double mass_b() const
Baryon mass.
bool is_irrotational() const
Returns true for an irrotational motion, false for a corotating one.
Definition etoile.h:1092
virtual double mass_g() const
Gravitational mass.
const Tenseur & get_nnn() const
Returns the total lapse function N.
Definition etoile.h:727
int get_nzet() const
Returns the number of domains occupied by the star.
Definition etoile.h:662
double ray_eq_pi() const
Coordinate radius at , [r_unit].
double ray_eq() const
Coordinate radius at , [r_unit].
const Map & get_mp() const
Returns the mapping.
Definition etoile.h:659
bool is_relativistic() const
Returns true for a relativistic star, false for a Newtonian one.
Definition etoile.h:667
const Tenseur & get_ent() const
Returns the enthalpy field.
Definition etoile.h:673
const Eos & get_eos() const
Returns the equation of state.
Definition etoile.h:670
const Tenseur & get_a_car() const
Returns the total conformal factor .
Definition etoile.h:733
const Tenseur & get_ener() const
Returns the proper total energy density.
Definition etoile.h:679
const Tenseur & get_press() const
Returns the fluid pressure.
Definition etoile.h:682
const Tenseur & get_nbar() const
Returns the proper baryon density.
Definition etoile.h:676
Base class for coordinate mappings.
Definition map.h:670
double get_ori_x() const
Returns the x coordinate of the origin.
Definition map.h:768
double get_rot_phi() const
Returns the angle between the x –axis and X –axis.
Definition map.h:775
Cmp sqrt(const Cmp &)
Square root.
Definition cmp_math.C:220
Cmp pow(const Cmp &, int)
Power .
Definition cmp_math.C:348
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
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 prototypes.
Definition app_hor.h:64
Standard units of space, time and mass.