LORENE
blackhole_static.C
1 /*
2  * Method of class Black_hole to set metric quantities to a spherical,
3  * static, analytic solution
4  *
5  * (see file blackhole.h for documentation).
6  *
7  */
8 
9 /*
10  * Copyright (c) 2005-2007 Keisuke Taniguchi
11  *
12  * This file is part of LORENE.
13  *
14  * LORENE is free software; you can redistribute it and/or modify
15  * it under the terms of the GNU General Public License version 2
16  * as published by the Free Software Foundation.
17  *
18  * LORENE is distributed in the hope that it will be useful,
19  * but WITHOUT ANY WARRANTY; without even the implied warranty of
20  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
21  * GNU General Public License for more details.
22  *
23  * You should have received a copy of the GNU General Public License
24  * along with LORENE; if not, write to the Free Software
25  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
26  *
27  */
28 
29 char blackhole_static_C[] = "$Header: /cvsroot/Lorene/C++/Source/Black_hole/blackhole_static.C,v 1.3 2014/10/13 08:52:46 j_novak Exp $" ;
30 
31 /*
32  * $Id: blackhole_static.C,v 1.3 2014/10/13 08:52:46 j_novak Exp $
33  * $Log: blackhole_static.C,v $
34  * Revision 1.3 2014/10/13 08:52:46 j_novak
35  * Lorene classes and functions now belong to the namespace Lorene.
36  *
37  * Revision 1.2 2008/05/15 19:31:17 k_taniguchi
38  * Change of some parameters.
39  *
40  * Revision 1.1 2007/06/22 01:20:50 k_taniguchi
41  * *** empty log message ***
42  *
43  *
44  * $Header: /cvsroot/Lorene/C++/Source/Black_hole/blackhole_static.C,v 1.3 2014/10/13 08:52:46 j_novak Exp $
45  *
46  */
47 
48 // C++ headers
49 //#include <>
50 
51 // C headers
52 //#include <math.h>
53 
54 // Lorene headers
55 #include "blackhole.h"
56 #include "unites.h"
57 #include "utilitaires.h"
58 
59 namespace Lorene {
60 void Black_hole::static_bh(bool neumann, bool first) {
61 
62  // Fundamental constants and units
63  // -------------------------------
64  using namespace Unites ;
65 
66  double mass = ggrav * mass_bh ;
67 
68  Scalar rr(mp) ;
69  rr = mp.r ;
70  rr.std_spectral_base() ;
71 
72  //-------------------------------------//
73  // Metric quantities //
74  //-------------------------------------//
75 
76  Scalar st(mp) ;
77  st = mp.sint ;
78  st.std_spectral_base() ;
79  Scalar ct(mp) ;
80  ct = mp.cost ;
81  ct.std_spectral_base() ;
82  Scalar sp(mp) ;
83  sp = mp.sinp ;
84  sp.std_spectral_base() ;
85  Scalar cp(mp) ;
86  cp = mp.cosp ;
87  cp.std_spectral_base() ;
88 
89  Vector ll(mp, CON, mp.get_bvect_cart()) ;
90  ll.set_etat_qcq() ;
91  ll.set(1) = st * cp ;
92  ll.set(2) = st * sp ;
93  ll.set(3) = ct ;
94  ll.std_spectral_base() ;
95 
96  if (kerrschild) {
97 
98  // Lapconf function
99  // ----------------
100  lapconf = 1./sqrt(1.+2.*mass/rr) ;
103  lapconf.raccord(1) ;
104 
105  // Conformal factor
106  // ----------------
107  confo = 1. ;
109 
110  // Shift vector
111  // ------------
112  for (int i=1; i<=3; i++)
113  shift.set(i) = 2. * mass * lapconf % lapconf % ll(i) / rr ;
114 
115  shift.annule_domain(0) ;
117 
118  }
119  else { // Isotropic coordinates with the maximal slicing
120 
121  // Sets C/M^2 for each case of the lapse boundary condition
122  // --------------------------------------------------------
123  double cc ;
124 
125  if (neumann) { // Neumann boundary condition
126  if (first) { // First condition
127  // d(\alpha \psi)/dr = 0
128  // ---------------------
129  cc = 2. * (sqrt(13.) - 1.) / 3. ;
130  }
131  else { // Second condition
132  // d(\alpha \psi)/dr = (\alpha \psi)/(2 rah)
133  // -----------------------------------------
134  cc = 4. / 3. ;
135  }
136  }
137  else { // Dirichlet boundary condition
138  if (first) { // First condition
139  // (\alpha \psi) = 1/2
140  // -------------------
141  cout << "!!!!! WARNING: Not yet prepared !!!!!" << endl ;
142  abort() ;
143  }
144  else { // Second condition
145  // (\alpha \psi) = 1/sqrt(2.) \psi_KS
146  // ----------------------------------
147  cout << "!!!!! WARNING: Not yet prepared !!!!!" << endl ;
148  abort() ;
149  }
150  }
151 
152  Scalar r_are(mp) ;
153  r_are = r_coord(neumann, first) ;
154  r_are.std_spectral_base() ;
155 
156  // Lapconf function
157  // ----------------
158  lapconf = sqrt(1. - 2.*mass/r_are/rr
159  + cc*cc*pow(mass/r_are/rr,4.)) * sqrt(r_are) ;
160 
163  lapconf.raccord(1) ;
164 
165  // Conformal factor
166  // ----------------
167  confo = sqrt(r_are) ;
169  confo.annule_domain(0) ;
170  confo.raccord(1) ;
171 
172  // Lapse function
173  // --------------
174  lapse = lapconf / confo ;
176  lapse.annule_domain(0) ;
177  lapse.raccord(1) ;
178 
179  // Shift vector
180  // ------------
181  for (int i=1; i<=3; i++) {
182  shift.set(i) = mass * mass * cc * ll(i) / rr / rr
183  / pow(r_are,3.) ;
184  }
185 
187 
188  for (int i=1; i<=3; i++) {
189  shift.set(i).annule_domain(0) ;
190  shift.set(i).raccord(1) ;
191  }
192 
193  }
194 
195 }
196 }
Map & mp
Mapping associated with the black hole.
Definition: blackhole.h:80
void annule_domain(int l)
Sets the Tensor to zero in a given domain.
Definition: tensor.C:666
virtual void set_etat_qcq()
Sets the logical state of all components to ETATQCQ (ordinary state).
Definition: tensor.C:481
Cmp sqrt(const Cmp &)
Square root.
Definition: cmp_math.C:220
Lorene prototypes.
Definition: app_hor.h:64
double mass_bh
Gravitational mass of BH.
Definition: blackhole.h:88
Standard units of space, time and mass.
void static_bh(bool neumann, bool first)
Sets the metric quantities to a spherical, static blach-hole analytic solution.
Tensor field of valence 0 (or component of a tensorial field).
Definition: scalar.h:387
bool kerrschild
true for a Kerr-Schild background, false for a conformally flat background
Definition: blackhole.h:85
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
Tensor field of valence 1.
Definition: vector.h:188
virtual void std_spectral_base()
Sets the standard spectal bases of decomposition for each component.
Definition: vector.C:316
void raccord(int n)
Performs the matching of the nucleus with respect to the first shell.
Scalar confo
Conformal factor generated by the black hole.
Definition: blackhole.h:118
Coord sint
Definition: map.h:721
Coord sinp
Definition: map.h:723
Cmp pow(const Cmp &, int)
Power .
Definition: cmp_math.C:348
const Base_vect_cart & get_bvect_cart() const
Returns the Cartesian basis associated with the coordinates (x,y,z) of the mapping, i.e.
Definition: map.h:791
Scalar & set(const Itbl &ind)
Returns the value of a component (read/write version).
Definition: tensor.C:654
Scalar lapse
Lapse function generated by the black hole.
Definition: blackhole.h:106
Vector shift
Shift vector generated by the black hole.
Definition: blackhole.h:109
Coord cosp
Definition: map.h:724
Scalar & set(int)
Read/write access to a component.
Definition: vector.C:296
const Scalar r_coord(bool neumann, bool first) const
Expresses the areal radial coordinate by that in spatially isotropic coordinates. ...
Scalar lapconf
A function (lapse function * conformal factor) lapconf generated by the black hole.
Definition: blackhole.h:97
Coord r
r coordinate centered on the grid
Definition: map.h:718
Coord cost
Definition: map.h:722