LORENE
blackhole.C
1 /*
2  * Methods of class Black_hole
3  *
4  * (see file blackhole.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 
28 char blackhole_C[] = "$Header: /cvsroot/Lorene/C++/Source/Black_hole/blackhole.C,v 1.5 2014/10/13 08:52:45 j_novak Exp $" ;
29 
30 /*
31  * $Id: blackhole.C,v 1.5 2014/10/13 08:52:45 j_novak Exp $
32  * $Log: blackhole.C,v $
33  * Revision 1.5 2014/10/13 08:52:45 j_novak
34  * Lorene classes and functions now belong to the namespace Lorene.
35  *
36  * Revision 1.4 2014/10/06 15:13:02 j_novak
37  * Modified #include directives to use c++ syntax.
38  *
39  * Revision 1.3 2008/07/03 14:55:26 k_taniguchi
40  * Addition of the angular momentum.
41  *
42  * Revision 1.2 2008/05/15 19:25:07 k_taniguchi
43  * Change of some parameters.
44  *
45  * Revision 1.1 2007/06/22 01:18:47 k_taniguchi
46  * *** empty log message ***
47  *
48  *
49  * $Header: /cvsroot/Lorene/C++/Source/Black_hole/blackhole.C,v 1.5 2014/10/13 08:52:45 j_novak Exp $
50  *
51  */
52 
53 // C++ headers
54 //#include <>
55 
56 // C headers
57 #include <cmath>
58 
59 // Lorene headers
60 #include "blackhole.h"
61 #include "utilitaires.h"
62 #include "unites.h"
63 
64  //----------------------//
65  // Constructors //
66  //----------------------//
67 
68 // Standard constructor
69 // --------------------
70 namespace Lorene {
71 Black_hole::Black_hole(Map& mp_i, bool kerrschild_i, double massbh)
72  : mp(mp_i),
73  kerrschild(kerrschild_i),
74  mass_bh(massbh),
75  lapconf(mp_i),
76  lapconf_rs(mp_i),
77  lapconf_bh(mp_i),
78  lapse(mp_i),
79  shift(mp_i, CON, mp_i.get_bvect_cart()),
80  shift_rs(mp_i, CON, mp_i.get_bvect_cart()),
81  shift_bh(mp_i, CON, mp_i.get_bvect_cart()),
82  confo(mp_i),
83  taij(mp_i, CON, mp_i.get_bvect_cart()),
84  taij_rs(mp_i, CON, mp_i.get_bvect_cart()),
85  taij_quad(mp_i),
86  taij_quad_rs(mp_i),
87  flat(mp_i, mp_i.get_bvect_cart()) {
88 
89  // Pointers of derived quantities are initialized to zero
90  set_der_0x0() ;
91 
92  // The metric quantities are initialized to the flat one
93  lapconf = 1. ;
95  lapconf_rs = 0. ;
97  lapconf_bh = 1. ;
99  lapse = 1. ;
101  shift.set_etat_zero() ;
104  confo = 1. ;
106 
107  taij.set_etat_zero() ;
109  taij_quad = 0. ;
111  taij_quad_rs = 0. ;
113 
114 }
115 
116 // Copy constructor
117 // ----------------
119  : mp(bh.mp),
120  kerrschild(bh.kerrschild),
121  mass_bh(bh.mass_bh),
122  lapconf(bh.lapconf),
123  lapconf_rs(bh.lapconf_rs),
124  lapconf_bh(bh.lapconf_bh),
125  lapse(bh.lapse),
126  shift(bh.shift),
127  shift_rs(bh.shift_rs),
128  shift_bh(bh.shift_bh),
129  confo(bh.confo),
130  taij(bh.taij),
131  taij_rs(bh.taij_rs),
132  taij_quad(bh.taij_quad),
133  taij_quad_rs(bh.taij_quad_rs),
134  flat(bh.flat) {
135 
136  set_der_0x0() ;
137 
138 }
139 
140 // Constructor from a file
141 // -----------------------
142 Black_hole::Black_hole(Map& mp_i, FILE* fich)
143  : mp(mp_i),
144  lapconf(mp_i),
145  lapconf_rs(mp_i, *(mp_i.get_mg()), fich),
146  lapconf_bh(mp_i),
147  lapse(mp_i),
148  shift(mp_i, CON, mp_i.get_bvect_cart()),
149  shift_rs(mp_i, mp_i.get_bvect_cart(), fich),
150  shift_bh(mp_i, CON, mp_i.get_bvect_cart()),
151  confo(mp_i, *(mp_i.get_mg()), fich),
152  taij(mp_i, CON, mp_i.get_bvect_cart()),
153  taij_rs(mp_i, CON, mp_i.get_bvect_cart()),
154  taij_quad(mp_i),
155  taij_quad_rs(mp_i),
156  flat(mp_i, mp_i.get_bvect_cart()) {
157 
158  // Background
159  // ----------
160  fread(&kerrschild, sizeof(bool), 1, fich) ;
161  fread(&mass_bh, sizeof(double), 1, fich) ;
162 
163  // All other fields are initialized to zero or some values
164  // -------------------------------------------------------
165  lapconf = lapconf_rs ;
167  lapconf_bh = 0. ;
169 
170  lapse = lapconf / confo ;
171 
172  shift = shift_rs ;
175 
176  taij.set_etat_zero() ;
178  taij_quad = 0. ;
180  taij_quad_rs = 0. ;
182 
183  // Pointers of derived quantities are initialized to zero
184  // ------------------------------------------------------
185  set_der_0x0() ;
186 
187 }
188 
189 
190  //--------------------//
191  // Destructor //
192  //--------------------//
193 
195 
196  del_deriv() ;
197 
198 }
199 
200 
201  //------------------------------------------//
202  // Management of derived quantities //
203  //------------------------------------------//
204 
205 void Black_hole::del_deriv() const {
206 
207  if (p_mass_irr != 0x0) delete p_mass_irr ;
208  if (p_mass_adm != 0x0) delete p_mass_adm ;
209  if (p_mass_kom != 0x0) delete p_mass_kom ;
210  if (p_rad_ah != 0x0) delete p_rad_ah ;
211  if (p_spin_am_bh != 0x0) delete p_spin_am_bh ;
212  if (p_angu_mom_bh != 0x0) delete p_angu_mom_bh ;
213 
214  set_der_0x0() ;
215 
216 }
217 
219 
220  p_mass_irr = 0x0 ;
221  p_mass_adm = 0x0 ;
222  p_mass_kom = 0x0 ;
223  p_rad_ah = 0x0 ;
224  p_spin_am_bh = 0x0 ;
225  p_angu_mom_bh = 0x0 ;
226 
227 }
228 
229 
230  //--------------------//
231  // Assignment //
232  //--------------------//
233 
234 // Assignment to another Black_hole
235 // --------------------------------
237 
238  assert( &(bh.mp) == &mp ) ; // Same mapping
239 
240  kerrschild = bh.kerrschild ;
241  mass_bh = bh.mass_bh ;
242  lapconf = bh.lapconf ;
243  lapconf_rs = bh.lapconf_rs ;
244  lapconf_bh = bh.lapconf_bh ;
245  lapse = bh.lapse ;
246  shift = bh.shift ;
247  shift_rs = bh.shift_rs ;
248  shift_bh = bh.shift_bh ;
249  confo = bh.confo ;
250  taij = bh.taij ;
251  taij_rs = bh.taij_rs ;
252  taij_quad = bh.taij_quad ;
254  flat = bh.flat ;
255 
256  del_deriv() ; // Deletes all derived quantities
257 
258 }
259 
260 
261  //-----------------//
262  // Outputs //
263  //-----------------//
264 
265 // Save in a file
266 // --------------
267 void Black_hole::sauve(FILE* fich) const {
268 
269  lapconf_rs.sauve(fich) ;
270  shift_rs.sauve(fich) ;
271  confo.sauve(fich) ;
272 
273  fwrite(&kerrschild, sizeof(bool), 1, fich) ;
274  fwrite(&mass_bh, sizeof(double), 1, fich) ;
275 
276 }
277 
278 // Printing
279 // --------
280 ostream& operator<<(ostream& ost, const Black_hole& bh) {
281 
282  bh >> ost ;
283  return ost ;
284 
285 }
286 
287 ostream& Black_hole::operator>>(ostream& ost) const {
288 
289  using namespace Unites ;
290 
291  const Mg3d* mg = mp.get_mg() ;
292  int nt = mg->get_nt(1) ;
293 
294  ost << endl ;
295  if (kerrschild) {
296  ost << "Kerr-Schild background" << endl ;
297  ost << "----------------------" << endl ;
298  }
299  else {
300  ost << "Conformally flat background" << endl ;
301  ost << "---------------------------" << endl ;
302  }
303 
304  ost << "lapconf on the AH : "
305  << lapconf.val_grid_point(1,0,nt-1,0) << endl ;
306  ost << "lapse on the AH : "
307  << lapse.val_grid_point(1,0,nt-1,0) << endl ;
308  ost << "shift(1) on the AH : "
309  << shift(1).val_grid_point(1,0,nt-1,0) << endl ;
310  ost << "shift(2) on the AH : "
311  << shift(2).val_grid_point(1,0,nt-1,0) << endl ;
312  ost << "shift(3) on the AH : "
313  << shift(3).val_grid_point(1,0,nt-1,0) << endl ;
314  ost << "confo on the AH : "
315  << confo.val_grid_point(1,0,nt-1,0) << endl ;
316  ost << "Gravitational mass : "
317  << mass_bh / msol << " M_sol" << endl ;
318  ost << "Irreducible mass : "
319  << mass_irr() / msol << " M_sol" << endl ;
320  ost << "ADM mass : "
321  << mass_adm() / msol << " M_sol" << endl ;
322  ost << "Komar mass : "
323  << mass_kom() / msol << " M_sol" << endl ;
324 
325  double irr_gm, adm_gm, kom_gm ;
326  irr_gm = mass_irr() / mass_bh - 1. ;
327  adm_gm = mass_adm() / mass_bh - 1. ;
328  kom_gm = mass_kom() / mass_bh - 1. ;
329  ost << "Diff. (Mirr-Mg)/Mg : " << irr_gm << endl ;
330  ost << "Diff. (Madm-Mg)/Mg : " << adm_gm << endl ;
331  ost << "Diff. (Mkom-Mg)/Mg : " << kom_gm << endl ;
332 
333  return ost ;
334 
335 }
336 }
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
virtual double mass_kom() const
Komar mass.
Vector shift_rs
Part of the shift vector from the numerical computation.
Definition: blackhole.h:112
Black_hole(Map &mp_i, bool Kerr_schild, double massbh)
Standard constructor.
Definition: blackhole.C:71
Scalar taij_quad
Part of the scalar generated by .
Definition: blackhole.h:135
Scalar taij_quad_rs
Part of the scalar.
Definition: blackhole.h:138
Sym_tensor taij
Trace of the extrinsic curvature.
Definition: blackhole.h:127
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
Scalar lapconf
A function (lapse function * conformal factor) lapconf generated by the black hole.
Definition: blackhole.h:97
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
Scalar lapconf_rs
Part of lapconf from the numerical computation.
Definition: blackhole.h:100
void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: blackhole.C:218
Tbl * p_angu_mom_bh
Spin angular momentum.
Definition: blackhole.h:162
virtual ~Black_hole()
Destructor.
Definition: blackhole.C:194
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
Scalar lapconf_bh
Part of lapconf from the analytic background.
Definition: blackhole.h:103
double * p_mass_kom
ADM mass.
Definition: blackhole.h:155
void operator=(const Black_hole &)
Assignment to another Black_hole.
Definition: blackhole.C:236
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
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<)
Definition: blackhole.C:287
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
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:765
Multi-domain grid.
Definition: grilles.h:273
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
Definition: grilles.h:457
virtual void sauve(FILE *) const
Save in a file.
Definition: scalar.C:686
virtual void std_spectral_base()
Sets the spectral bases of the Valeur va to the standard ones for a scalar field.
Definition: scalar.C:784
double val_grid_point(int l, int k, int j, int i) const
Returns the value of the field at a specified grid point.
Definition: scalar.h:637
virtual void std_spectral_base()
Sets the standard spectal bases of decomposition for each component.
Definition: vector.C:316
virtual void sauve(FILE *) const
Save in a binary file.
Definition: tensor.C:906
virtual void set_etat_zero()
Sets the logical state of all components to ETATZERO (zero state).
Definition: tensor.C:497
Lorene prototypes.
Definition: app_hor.h:64
Standard units of space, time and mass.