LORENE
blackhole.h
1/*
2 * Definition of Lorene class Black_hole
3 *
4 */
5
6/*
7 * Copyright (c) 2005-2007 Keisuke Taniguchi
8 *
9 * This file is part of LORENE.
10 *
11 * LORENE is free software; you can redistribute it and/or modify
12 * it under the terms of the GNU General Public License version 2
13 * as published by the Free Software Foundation.
14 *
15 * LORENE is distributed in the hope that it will be useful,
16 * but WITHOUT ANY WARRANTY; without even the implied warranty of
17 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18 * GNU General Public License for more details.
19 *
20 * You should have received a copy of the GNU General Public License
21 * along with LORENE; if not, write to the Free Software
22 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
23 *
24 */
25
26#ifndef __BLACKHOLE_H_
27#define __BLACKHOLE_H_
28
29/*
30 * $Id: blackhole.h,v 1.4 2014/10/13 08:52:32 j_novak Exp $
31 * $Log: blackhole.h,v $
32 * Revision 1.4 2014/10/13 08:52:32 j_novak
33 * Lorene classes and functions now belong to the namespace Lorene.
34 *
35 * Revision 1.3 2008/07/02 20:41:40 k_taniguchi
36 * Addition of the routines to compute angular momentum
37 * and modification of the argument of equilibrium_spher.
38 *
39 * Revision 1.2 2008/05/15 18:53:37 k_taniguchi
40 * Change of some parameters and introduction of some
41 * computational routines.
42 *
43 * Revision 1.1 2007/06/22 01:02:57 k_taniguchi
44 * *** empty log message ***
45 *
46 *
47 * $Header: /cvsroot/Lorene/C++/Include/blackhole.h,v 1.4 2014/10/13 08:52:32 j_novak Exp $
48 *
49 */
50
51// External classes which appear in the declaration of class Black_hole:
52//class YYY ;
53
54// Headers Lorene
55#include "metric.h"
56
57 //-------------------------------//
58 // Base class Black_hole //
59 //-------------------------------//
60
61namespace Lorene {
75
76 // Data :
77 // -----
78 protected:
81
86
88 double mass_bh ;
89
90 // Metric quantities
91 // -----------------
92
97 Scalar lapconf ; // lapconf = lapconf_rs + lapconf_bh
98
101
104
107
109 Vector shift ; // shift = shift_rs + shift_bh
110
113
116
119
121 // Scalar trace_k ;
122
128
131
136
139
144
146 // Metric tgij ;
147
148 // Derived data :
149 // ------------
150 protected:
151 mutable double* p_mass_irr ;
152
153 mutable double* p_mass_adm ;
154
155 mutable double* p_mass_kom ;
156
157 mutable double* p_rad_ah ;
158
159 mutable double* p_spin_am_bh ;
160
163
164 // Constructors - Destructor
165 // -------------------------
166 public:
167
172 Black_hole(Map& mp_i, bool Kerr_schild, double massbh) ;
173
174 Black_hole(const Black_hole& ) ;
175
182 Black_hole(Map& mp_i, FILE* fich) ;
183
184 virtual ~Black_hole() ;
185
186
187 // Memory management
188 // -----------------
189 protected:
191 virtual void del_deriv() const ;
192
194 void set_der_0x0() const ;
195
196
197 // Mutators / assignment
198 // ---------------------
199 public:
201 void operator=(const Black_hole&) ;
202
204 Map& set_mp() {return mp; } ;
205
207 double& set_mass_bh() {return mass_bh; } ;
208
209 // Accessors
210 // ---------
211 public:
213 const Map& get_mp() const {return mp; } ;
214
218 bool is_kerrschild() const {return kerrschild; } ;
219
221 double get_mass_bh() const {return mass_bh; } ;
222
224 const Scalar& get_lapconf() const {return lapconf; } ;
225
227 const Scalar& get_lapconf_rs() const {return lapconf_rs; } ;
228
230 const Scalar& get_lapse() const {return lapse; } ;
231
233 const Vector& get_shift() const {return shift; } ;
234
238 const Vector& get_shift_rs() const {return shift_rs; } ;
239
241 const Scalar& get_confo() const {return confo; } ;
242
243 // Outputs
244 // -------
245 public:
246 virtual void sauve(FILE *) const ;
247
249 friend ostream& operator<<(ostream& , const Black_hole& ) ;
250
251 protected:
253 virtual ostream& operator>>(ostream& ) const ;
254
255 // Global quantities
256 // -----------------
257 public:
259 virtual double mass_irr() const ;
260
262 virtual double mass_adm() const ;
263
265 virtual double mass_kom() const ;
266
268 virtual double rad_ah() const ;
269
271 double spin_am_bh(bool bclapconf_nd, bool bclapconf_fs,
272 const Tbl& xi_i, const double& phi_i,
273 const double& theta_i, const int& nrk_phi,
274 const int& nrk_theta) const ;
275
283 const Tbl& angu_mom_bh() const ;
284
285 // Computational routines
286 // ----------------------
287 public:
291 const Valeur bc_lapconf(bool neumann, bool first) const ;
292
296 const Valeur bc_shift_x(double omega_r) const ;
297
301 const Valeur bc_shift_y(double omega_r) const ;
302
306 const Valeur bc_shift_z() const ;
307
311 const Valeur bc_confo() const ;
312
316 void extr_curv_bh() ;
317
335 void equilibrium_spher(bool neumann, bool first, double spin_omega,
336 double precis = 1.e-14,
337 double precis_shift = 1.e-8) ;
338
347 void static_bh(bool neumann, bool first) ;
348
357 double rah_iso(bool neumann, bool first) const ;
358
367 const Scalar r_coord(bool neumann, bool first) const ;
368
380 Tbl runge_kutta_phi_bh(const Tbl& xi_i, const double& phi_i,
381 const int& nrk) const ;
382
396 Tbl runge_kutta_theta_bh(const Tbl& xi_i, const double& theta_i,
397 const double& phi, const int& nrk) const ;
398
415 Vector killing_vect_bh(const Tbl& xi_i, const double& phi_i,
416 const double& theta_i, const int& nrk_phi,
417 const int& nrk_theta) const ;
418
419};
420ostream& operator<<(ostream& , const Black_hole& ) ;
421
422}
423#endif
Base class for black holes.
Definition blackhole.h:74
Vector shift_bh
Part of the shift vector from the analytic background.
Definition blackhole.h:115
double rah_iso(bool neumann, bool first) const
Computes a radius of apparent horizon (excised surface) in isotropic coordinates.
virtual double mass_kom() const
Komar mass.
Vector shift_rs
Part of the shift vector from the numerical computation.
Definition blackhole.h:112
Tbl runge_kutta_phi_bh(const Tbl &xi_i, const double &phi_i, const int &nrk) const
Compute a forth-order Runge-Kutta integration to the phi direction for the solution of the Killing ve...
double spin_am_bh(bool bclapconf_nd, bool bclapconf_fs, const Tbl &xi_i, const double &phi_i, const double &theta_i, const int &nrk_phi, const int &nrk_theta) const
Spin angular momentum.
const Valeur bc_confo() const
Boundary condition on the apparent horizon of the black hole for the conformal factor: 2-D Valeur.
double & set_mass_bh()
Read/write of the gravitational mass of BH [{\tt m_unit}].
Definition blackhole.h:207
Scalar taij_quad
Part of the scalar generated by .
Definition blackhole.h:135
const Scalar & get_lapse() const
Returns the lapse function generated by the black hole.
Definition blackhole.h:230
const Valeur bc_shift_x(double omega_r) const
Boundary condition on the apparent horizon of the black hole for the shift vector of the \fx\f direct...
Scalar taij_quad_rs
Part of the scalar.
Definition blackhole.h:138
Sym_tensor taij
Trace of the extrinsic curvature.
Definition blackhole.h:127
const Scalar r_coord(bool neumann, bool first) const
Expresses the areal radial coordinate by that in spatially isotropic coordinates.
const Vector & get_shift() const
Returns the shift vector generated by the black hole.
Definition blackhole.h:233
Map & set_mp()
Read/write of the mapping.
Definition blackhole.h:204
virtual double mass_irr() const
Irreducible mass of the black hole.
double * p_spin_am_bh
Radius of the apparent horizon.
Definition blackhole.h:159
Sym_tensor taij_rs
Part of the extrinsic curvature tensor.
Definition blackhole.h:130
Map & mp
Mapping associated with the black hole.
Definition blackhole.h:80
bool is_kerrschild() const
Returns true for a Kerr-Schild background, false for a Conformally flat one.
Definition blackhole.h:218
Scalar lapconf
A function (lapse function * conformal factor) lapconf generated by the black hole.
Definition blackhole.h:97
const Valeur bc_lapconf(bool neumann, bool first) const
Boundary condition on the apparent horizon of the black hole for the lapse function: 2-D Valeur.
Vector shift
Shift vector generated by the black hole.
Definition blackhole.h:109
Metric_flat flat
Flat metric defined on the mapping (Spherical components with respect to the mapping of the black hol...
Definition blackhole.h:143
virtual double rad_ah() const
Radius of the apparent horizon.
const Scalar & get_lapconf_rs() const
Returns the part of lapconf from the numerical computation.
Definition blackhole.h:227
Scalar lapconf_rs
Part of lapconf from the numerical computation.
Definition blackhole.h:100
const Scalar & get_confo() const
Returns the conformal factor generated by the black hole.
Definition blackhole.h:241
void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition blackhole.C:218
void static_bh(bool neumann, bool first)
Sets the metric quantities to a spherical, static blach-hole analytic solution.
Tbl runge_kutta_theta_bh(const Tbl &xi_i, const double &theta_i, const double &phi, const int &nrk) const
Compute a forth-order Runge-Kutta integration to the theta direction for the solution of the Killing ...
friend ostream & operator<<(ostream &, const Black_hole &)
Display.
Definition blackhole.C:280
void equilibrium_spher(bool neumann, bool first, double spin_omega, double precis=1.e-14, double precis_shift=1.e-8)
Computes a spherical, static black-hole by giving boundary conditions on the apparent horizon.
const Map & get_mp() const
Returns the mapping.
Definition blackhole.h:213
Tbl * p_angu_mom_bh
Spin angular momentum.
Definition blackhole.h:162
virtual ~Black_hole()
Destructor.
Definition blackhole.C:194
const Scalar & get_lapconf() const
Returns lapconf generated by the black hole.
Definition blackhole.h:224
Vector killing_vect_bh(const Tbl &xi_i, const double &phi_i, const double &theta_i, const int &nrk_phi, const int &nrk_theta) const
Compute the Killing vector of a black hole normalized so that its affine length is 2 M_PI.
virtual double mass_adm() const
ADM mass.
double * p_rad_ah
Komar mass.
Definition blackhole.h:157
virtual void sauve(FILE *) const
Save in a file.
Definition blackhole.C:267
const Valeur bc_shift_z() const
Boundary condition on the apparent horizon of the black hole for the shift vector of the \fz\f direct...
Scalar lapconf_bh
Part of lapconf from the analytic background.
Definition blackhole.h:103
const Valeur bc_shift_y(double omega_r) const
Boundary condition on the apparent horizon of the black hole for the shift vector of the \fy\f direct...
double * p_mass_kom
ADM mass.
Definition blackhole.h:155
void operator=(const Black_hole &)
Assignment to another Black_hole.
Definition blackhole.C:236
const Vector & get_shift_rs() const
Returns the part of the shift vector from the numerical computation.
Definition blackhole.h:238
const Tbl & angu_mom_bh() const
Total angular momentum.
Scalar lapse
Lapse function generated by the black hole.
Definition blackhole.h:106
bool kerrschild
true for a Kerr-Schild background, false for a conformally flat background
Definition blackhole.h:85
Scalar confo
Conformal factor generated by the black hole.
Definition blackhole.h:118
double mass_bh
Gravitational mass of BH.
Definition blackhole.h:88
double * p_mass_adm
Irreducible mass of the black hole.
Definition blackhole.h:153
void extr_curv_bh()
Computes taij , taij_quad from shift , lapse , confo .
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<)
Definition blackhole.C:287
double get_mass_bh() const
Returns the gravitational mass of BH [{\tt m_unit}].
Definition blackhole.h:221
double * p_mass_irr
Conformal metric .
Definition blackhole.h:151
virtual void del_deriv() const
Deletes all the derived quantities.
Definition blackhole.C:205
Base class for coordinate mappings.
Definition map.h:670
Flat metric for tensor calculation.
Definition metric.h:261
Tensor field of valence 0 (or component of a tensorial field).
Definition scalar.h:387
Class intended to describe valence-2 symmetric tensors.
Definition sym_tensor.h:223
Basic array class.
Definition tbl.h:161
Values and coefficients of a (real-value) function.
Definition valeur.h:287
Tensor field of valence 1.
Definition vector.h:188
Lorene prototypes.
Definition app_hor.h:64