LORENE
bin_bhns.C
1/*
2 * Methods of class Bin_bhns
3 *
4 * (see file bin_bhns.h for documentation).
5 *
6 */
7
8/*
9 * Copyright (c) 2005-2007 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_C[] = "$Header: /cvsroot/Lorene/C++/Source/Bin_bhns/bin_bhns.C,v 1.4 2014/10/13 08:52:41 j_novak Exp $" ;
29
30/*
31 * $Id: bin_bhns.C,v 1.4 2014/10/13 08:52:41 j_novak Exp $
32 * $Log: bin_bhns.C,v $
33 * Revision 1.4 2014/10/13 08:52:41 j_novak
34 * Lorene classes and functions now belong to the namespace Lorene.
35 *
36 * Revision 1.3 2014/10/06 15:13:00 j_novak
37 * Modified #include directives to use c++ syntax.
38 *
39 * Revision 1.2 2008/05/15 18:58:01 k_taniguchi
40 * Change of some parameters and introduction of new
41 * global quantities.
42 *
43 * Revision 1.1 2007/06/22 01:08:53 k_taniguchi
44 * *** empty log message ***
45 *
46 *
47 * $Header: /cvsroot/Lorene/C++/Source/Bin_bhns/bin_bhns.C,v 1.4 2014/10/13 08:52:41 j_novak Exp $
48 *
49 */
50
51// C++ headers
52//#include <>
53
54// C headers
55#include <cmath>
56
57// Lorene headers
58#include "bin_bhns.h"
59#include "eos.h"
60#include "utilitaires.h"
61#include "unites.h"
62
63 //---------------------//
64 // Constructor //
65 //---------------------//
66
67// Standard constructor
68// --------------------
69namespace Lorene {
70Bin_bhns::Bin_bhns(Map& mp_bh, Map& mp_ns, int nzet_i, const Eos& eos_i,
71 bool irrot_ns, bool kerrschild_i,
72 bool bc_lapconf_nd, bool bc_lapconf_fs, bool irrot_bh,
73 double massbh)
74 : ref_triad(0., "Absolute frame Cartesian basis"),
75 hole(mp_bh,kerrschild_i,bc_lapconf_nd,bc_lapconf_fs,irrot_bh, massbh),
76 star(mp_ns, nzet_i, eos_i, irrot_ns)
77{
78
79 omega = 0. ;
80 separ = 0. ;
81 x_rot = 0. ;
82 y_rot = 0. ;
83
84 // Pointers of derived quantities are initialized to zero
85 set_der_0x0() ;
86
87}
88
89// Copy constructor
90// ----------------
92 : ref_triad(0., "Absolute frame Cartesian basis"),
93 hole(bibi.hole),
94 star(bibi.star),
95 omega(bibi.omega),
96 separ(bibi.separ),
97 x_rot(bibi.x_rot),
98 y_rot(bibi.y_rot)
99{
100
101 // Pointers of derived quantities are initialized to zero
102 set_der_0x0() ;
103
104}
105
106// Constructor from a file
107// -----------------------
108Bin_bhns::Bin_bhns(Map& mp_bh, Map& mp_ns, const Eos& eos_i, FILE* fich)
109 : ref_triad(0., "Absolute frame Cartesian basis"),
110 hole(mp_bh, fich),
111 star(mp_ns, eos_i, fich)
112{
113
114 // omega, separ, and x_rot are read from the file
115 fread_be(&omega, sizeof(double), 1, fich) ;
116 fread_be(&separ, sizeof(double), 1, fich) ;
117 fread_be(&x_rot, sizeof(double), 1, fich) ;
118 fread_be(&y_rot, sizeof(double), 1, fich) ;
119
120 // Pointers of derived quantities are initialized to zero
121 set_der_0x0() ;
122
123}
124
125
126 //--------------------//
127 // Destructor //
128 //--------------------//
129
131{
132
133 del_deriv() ;
134
135}
136
137
138 //------------------------------------------//
139 // Management of derived quantities //
140 //------------------------------------------//
141
143
144 if (p_mass_adm_bhns_surf != 0x0) delete p_mass_adm_bhns_surf ;
145 if (p_mass_adm_bhns_vol != 0x0) delete p_mass_adm_bhns_vol ;
146 if (p_mass_kom_bhns_surf != 0x0) delete p_mass_kom_bhns_surf ;
147 if (p_mass_kom_bhns_vol != 0x0) delete p_mass_kom_bhns_vol ;
148 if (p_line_mom_bhns != 0x0) delete p_line_mom_bhns ;
149 if (p_angu_mom_bhns != 0x0) delete p_angu_mom_bhns ;
150 if (p_virial_bhns_surf != 0x0) delete p_virial_bhns_surf ;
151 if (p_virial_bhns_vol != 0x0) delete p_virial_bhns_vol ;
152 if (p_xa_barycenter != 0x0) delete p_xa_barycenter ;
153 if (p_ya_barycenter != 0x0) delete p_ya_barycenter ;
154 if (p_omega_two_points != 0x0) delete p_omega_two_points ;
155 // if (p_ham_constr_bhns != 0x0) delete p_ham_constr_bhns ;
156 // if (p_mom_constr_bhns != 0x0) delete p_mom_constr_bhns ;
157
158 set_der_0x0() ;
159
160}
161
163
165 p_mass_adm_bhns_vol = 0x0 ;
167 p_mass_kom_bhns_vol = 0x0 ;
168 p_line_mom_bhns = 0x0 ;
169 p_angu_mom_bhns = 0x0 ;
170 p_virial_bhns_surf = 0x0 ;
171 p_virial_bhns_vol = 0x0 ;
172 p_xa_barycenter = 0x0 ;
173 p_ya_barycenter = 0x0 ;
174 p_omega_two_points = 0x0 ;
175 // p_ham_constr_bhns = 0x0 ;
176 // p_mom_constr_bhns = 0x0 ;
177
178}
179
180
181 //--------------------//
182 // Assignment //
183 //--------------------//
184
185// Assignment to anothe Bin_bhns
186// -----------------------------
187void Bin_bhns::operator=(const Bin_bhns& bibi) {
188
189 assert( bibi.ref_triad == ref_triad ) ;
190
191 hole = bibi.hole ;
192 star = bibi.star ;
193
194 omega = bibi.omega ;
195 separ = bibi.separ ;
196 x_rot = bibi.x_rot ;
197 y_rot = bibi.y_rot ;
198
199 del_deriv() ; // Deletes all derived quantities
200
201}
202
203
204 //-----------------//
205 // Outputs //
206 //-----------------//
207
208// Save in a file
209// --------------
210void Bin_bhns::sauve(FILE* fich) const {
211
212 hole.sauve(fich) ;
213 star.sauve(fich) ;
214
215 fwrite_be(&omega, sizeof(double), 1, fich) ;
216 fwrite_be(&separ, sizeof(double), 1, fich) ;
217 fwrite_be(&x_rot, sizeof(double), 1, fich) ;
218 fwrite_be(&y_rot, sizeof(double), 1, fich) ;
219
220}
221
222// Printing
223// --------
224ostream& operator<<(ostream& ost, const Bin_bhns& bibi) {
225
226 bibi >> ost ;
227 return ost ;
228
229}
230
231ostream& Bin_bhns::operator>>(ostream& ost) const {
232
233 using namespace Unites ;
234
235 ost << endl ;
236 ost << "BH-NS binary system" << endl ;
237 ost << "===================" << endl ;
238
239 ost << endl
240 << "Orbital angular velocity : "
241 << omega * f_unit << " [rad/s]" << endl ;
242 ost << "Coordinate separation between BH and NS : "
243 << separ / km << " [km]" << endl ;
244 ost << "Position of the rotation axis : "
245 << x_rot / km << " [km]" << endl ;
246
247 ost << endl << "Black hole : " ;
248 ost << endl << "========== " << endl ;
249
250 int nt = (hole.get_mp()).get_mg()->get_nt(1) ;
251
252 ost << "Irreducible mass of BH : "
253 << hole.mass_irr_bhns() / msol << " [Mo]" << endl ;
254 ost << "Mass in the background : "
255 << hole.get_mass_bh() / msol << " [Mo]" << endl ;
256 ost << "Radius of the apparent horison : "
257 << hole.rad_ah() / km << " [km]" << endl ;
258 ost << "Lapse function on the AH : "
259 << (hole.get_lapse_tot()).val_grid_point(1,0,nt-1,0) << endl ;
260 ost << "Conformal factor on the AH : "
261 << (hole.get_confo_tot()).val_grid_point(1,0,nt-1,0) << endl ;
262 ost << "shift(1) on the AH : "
263 << (hole.get_shift_tot()(1)).val_grid_point(1,0,nt-1,0) << endl ;
264 ost << "shift(2) on the AH : "
265 << (hole.get_shift_tot()(2)).val_grid_point(1,0,nt-1,0) << endl ;
266 ost << "shift(3) on the AH : "
267 << (hole.get_shift_tot()(3)).val_grid_point(1,0,nt-1,0) << endl ;
268
269 ost << endl << "Neutron star : " ;
270 ost << endl << "============ " << endl ;
271
272 ost << "Baryon mass of NS in isolation : "
273 << star.mass_b()/msol << " [Mo]" << endl ;
274 ost << "Gravitational mass of NS : "
275 << star.mass_g_bhns()/msol << " [Mo]" << endl ;
276 ost << "Baryon mass of NS in BH bg. : "
277 << star.mass_b_bhns(hole.is_kerrschild(),hole.get_mass_bh(),separ)/msol
278 << " [Mo]" << endl ;
279 ost << "Coordinate radius R_eq_tow : "
280 << star.ray_eq_pi() / km << " [km]" << endl ;
281 ost << "Coordinate radius R_eq_opp : "
282 << star.ray_eq() / km << " [km]" << endl ;
283 ost << "Coordinate radius R_eq_orb : "
284 << star.ray_eq_pis2() / km << " [km]" << endl ;
285 ost << "Coordinate radius R_eq_orb_opp : "
286 << star.ray_eq_3pis2() / km << " [km]" << endl ;
287 ost << "Coordinate radius R_pole : "
288 << star.ray_pole() / km << " [km]" << endl ;
289 ost << "Central enthalpy H_ent : "
290 << (star.get_ent()).val_grid_point(0,0,0,0) << endl ;
291 ost << "Lapse function at the center of NS : "
292 << (star.get_lapse_tot()).val_grid_point(0,0,0,0) << endl ;
293 ost << "Conformal factor at the center of NS : "
294 << (star.get_confo_tot()).val_grid_point(0,0,0,0) << endl ;
295 ost << "shift(1) at the center of NS : "
296 << (star.get_shift_tot()(1)).val_grid_point(0,0,0,0) << endl ;
297 ost << "shift(2) at the center of NS : "
298 << (star.get_shift_tot()(2)).val_grid_point(0,0,0,0) << endl ;
299 ost << "shift(3) at the center of NS : "
300 << (star.get_shift_tot()(3)).val_grid_point(0,0,0,0) << endl ;
301
302
303 return ost ;
304
305}
306
307// Display in polytropic unites
308// ----------------------------
309void Bin_bhns::display_poly(ostream& ost) const {
310
311 using namespace Unites ;
312
313 const Eos* p_eos = &( star.get_eos() ) ;
314 const Eos_poly* p_eos_poly = dynamic_cast<const Eos_poly*>( p_eos ) ;
315
316 if (p_eos_poly != 0x0) {
317
318 double kappa = p_eos_poly->get_kap() ;
319 double gamma = p_eos_poly->get_gam() ;
320 double kap_ns2 = pow( kappa, 0.5 /(gamma-1) ) ;
321
322 // Polytropic unit of length in terms of r_unit :
323 double r_poly = kap_ns2 / sqrt(ggrav) ;
324
325 // Polytropic unit of time in terms of t_unit :
326 double t_poly = r_poly ;
327
328 // Polytropic unit of mass in terms of m_unit :
329 double m_poly = r_poly / ggrav ;
330
331 // Polytropic unit of angular momentum in terms of j_unit :
332 // double j_poly = r_poly * r_poly / ggrav ;
333
334 double r0 = 0.5 * ( star.ray_eq() + star.ray_eq_pi() ) ;
335 // = 1 in Baumgarte et al.
336 // double d_ns = separ + star.ray_eq() - r0 ;
337 // Orbital separation of Baumgarte et al.
338
339 double y_ns = star.get_mp().get_ori_y() ;
340 double d_coord = sqrt(separ*separ + y_ns*y_ns) ;
341
342 ost.precision(16) ;
343 ost << endl << "Quantities in polytropic units : " ;
344 ost << endl << "==============================" << endl ;
345 ost << " ( r_poly = " << r_poly / km << " km )" << endl ;
346 ost << " d_separ : " << separ / r_poly << endl ;
347 ost << " d_coord : " << d_coord / r_poly << endl ;
348 ost << " Omega : " << omega * t_poly << endl ;
349 ost << " Omega M_irr : "
350 << omega * t_poly * hole.mass_irr_bhns() / m_poly << endl ;
351 ost << " Omega_spin : " << hole.get_omega_spin() * t_poly << endl ;
352 ost << " M_irr : " << hole.mass_irr_bhns() / m_poly << endl ;
353 ost << " M_bh : " << hole.get_mass_bh() / m_poly << endl ;
354 ost << " R_ah : " << hole.rad_ah() / r_poly << endl ;
355 ost << " M_bar(NS)iso : " << star.mass_b() / m_poly
356 << endl ;
357 ost << " M_bar(NS) : "
358 << star.mass_b_bhns(hole.is_kerrschild(),hole.get_mass_bh(),separ)
359 / m_poly << endl ;
360 ost << " M_g(NS)iso : " << star.mass_g_bhns() / m_poly << endl ;
361 ost << " R_0(NS) : " << r0 / r_poly << endl ;
362
363 }
364
365}
366}
Class for computing a black hole - neutron star binary system with comparable mass ()
Definition bin_bhns.h:63
void display_poly(ostream &) const
Display in polytropic units.
Definition bin_bhns.C:309
void del_deriv() const
Deletes all the derived quantities.
Definition bin_bhns.C:142
double * p_mass_adm_bhns_surf
Total ADM mass of the system calculated by the surface integral at infinity.
Definition bin_bhns.h:97
virtual void sauve(FILE *) const
Save in a file.
Definition bin_bhns.C:210
Bin_bhns(Map &mp_bh, Map &mp_ns, int nzet, const Eos &eos, bool irrot_ns, bool kerrschild, bool bc_lapse_nd, bool bc_lapse_fs, bool irrot_bh, double mass_bh)
Relative error on the Hamiltonian constraint.
Definition bin_bhns.C:70
double y_rot
Absolute Y coordinate of the rotation axis.
Definition bin_bhns.h:89
const Base_vect_cart ref_triad
Cartesian triad of the absolute reference frame.
Definition bin_bhns.h:69
double * p_mass_kom_bhns_surf
Total Komar mass of the system calculated by the surface integral at infinity.
Definition bin_bhns.h:107
Hole_bhns hole
Black hole.
Definition bin_bhns.h:72
void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition bin_bhns.C:162
double omega
Angular velocity with respect to an asymptotically inertial observer.
Definition bin_bhns.h:80
double * p_mass_adm_bhns_vol
Total ADM mass of the system calculated by the volume integral and the surface integral at the appare...
Definition bin_bhns.h:102
double * p_virial_bhns_vol
Virial theorem error calculated by the ADM mass and the Komar mass of the volume integral.
Definition bin_bhns.h:128
double * p_omega_two_points
Orbital angular velocity derived from another method.
Definition bin_bhns.h:137
double separ
Absolute orbital separation between two centers of BH and NS.
Definition bin_bhns.h:83
Tbl * p_angu_mom_bhns
Total angular momentum of the system.
Definition bin_bhns.h:118
double * p_ya_barycenter
Absolute coordinate Y of the barycenter of the baryon density.
Definition bin_bhns.h:134
double x_rot
Absolute X coordinate of the rotation axis.
Definition bin_bhns.h:86
virtual ~Bin_bhns()
Destructor.
Definition bin_bhns.C:130
void operator=(const Bin_bhns &)
Assignment to another Bin_bhns.
Definition bin_bhns.C:187
double * p_virial_bhns_surf
Virial theorem error calculated by the ADM mass and the Komar mass of the surface integral at infinit...
Definition bin_bhns.h:123
double * p_mass_kom_bhns_vol
Total Komar mass of the system calculated by the volume integral and the surface integral at the appa...
Definition bin_bhns.h:112
Star_bhns star
Neutron star.
Definition bin_bhns.h:75
double * p_xa_barycenter
Absolute coordinate X of the barycenter of the baryon density.
Definition bin_bhns.h:131
ostream & operator>>(ostream &) const
Operator >> (function called by the operator <<)
Definition bin_bhns.C:231
Tbl * p_line_mom_bhns
Total linear momentum of the system.
Definition bin_bhns.h:115
bool is_kerrschild() const
Returns true for a Kerr-Schild background, false for a Conformally flat one.
Definition blackhole.h:218
virtual double rad_ah() const
Radius of the apparent horizon.
const Map & get_mp() const
Returns the mapping.
Definition blackhole.h:213
double get_mass_bh() const
Returns the gravitational mass of BH [{\tt m_unit}].
Definition blackhole.h:221
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
double get_omega_spin() const
Returns the spin angular velocity of the black hole [{\tt f_unit}].
Definition hole_bhns.h:362
const Vector & get_shift_tot() const
Returns the total shift vector.
Definition hole_bhns.h:416
virtual double mass_irr_bhns() const
Irreducible mass of the black hole.
const Scalar & get_confo_tot() const
Returns the total conformal factor.
Definition hole_bhns.h:447
const Scalar & get_lapse_tot() const
Returns the total lapse function.
Definition hole_bhns.h:386
virtual void sauve(FILE *) const
Save in a file.
Definition hole_bhns.C:603
Base class for coordinate mappings.
Definition map.h:670
double get_ori_y() const
Returns the y coordinate of the origin.
Definition map.h:770
virtual void sauve(FILE *) const
Save in a file.
Definition star_bhns.C:470
const Scalar & get_lapse_tot() const
Returns the total lapse function.
Definition star_bhns.h:353
const Scalar & get_confo_tot() const
Returns the total conformal factor.
Definition star_bhns.h:391
const Vector & get_shift_tot() const
Returns the part total shift vector.
Definition star_bhns.h:372
virtual double mass_b() const
Baryon mass.
const Map & get_mp() const
Returns the mapping.
Definition star.h:355
double ray_eq_3pis2() const
Coordinate radius at , [r_unit].
const Eos & get_eos() const
Returns the equation of state.
Definition star.h:361
double ray_eq() const
Coordinate radius at , [r_unit].
double ray_eq_pis2() const
Coordinate radius at , [r_unit].
const Scalar & get_ent() const
Returns the enthalpy field.
Definition star.h:364
double ray_eq_pi() const
Coordinate radius at , [r_unit].
double ray_pole() const
Coordinate radius at [r_unit].
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.