LORENE
hole_bhns_upmetr.C
1/*
2 * Method of class Hole_bhns to compute metric quantities from
3 * the companion neutron-star and total metric quantities
4 *
5 * (see file hole_bhns.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
29char hole_bhns_upmetr_C[] = "$Header: /cvsroot/Lorene/C++/Source/Hole_bhns/hole_bhns_upmetr.C,v 1.3 2014/10/13 08:53:00 j_novak Exp $" ;
30
31/*
32 * $Id: hole_bhns_upmetr.C,v 1.3 2014/10/13 08:53:00 j_novak Exp $
33 * $Log: hole_bhns_upmetr.C,v $
34 * Revision 1.3 2014/10/13 08:53:00 j_novak
35 * Lorene classes and functions now belong to the namespace Lorene.
36 *
37 * Revision 1.2 2008/05/15 19:08:15 k_taniguchi
38 * Change of some parameters.
39 *
40 * Revision 1.1 2007/06/22 01:25:31 k_taniguchi
41 * *** empty log message ***
42 *
43 *
44 * $Header: /cvsroot/Lorene/C++/Source/Hole_bhns/hole_bhns_upmetr.C,v 1.3 2014/10/13 08:53:00 j_novak Exp $
45 *
46 */
47
48// C++ headers
49//#include <>
50
51// C headers
52//#include <>
53
54// Lorene headers
55#include "hole_bhns.h"
56#include "star_bhns.h"
57#include "utilitaires.h"
58#include "unites.h"
59
60namespace Lorene {
62 const Hole_bhns& hole_prev,
63 double relax) {
64
65 // Fundamental constants and units
66 // -------------------------------
67 using namespace Unites ;
68
69 //-----------------------------------------------------
70 // Computation of quantities coming from the companion
71 //-----------------------------------------------------
72
73 const Map& mp_ns (star.get_mp()) ;
74
75 double mass = ggrav * mass_bh ;
76
77 // Lapconf function
78 // ----------------
79
80 if ( (star.get_lapconf_auto()).get_etat() == ETATZERO ) {
82 }
83 else {
87 }
88
89 // Shift vector
90 // ------------
91
92 if ( (star.get_shift_auto())(2).get_etat() == ETATZERO ) {
93 assert( (star.get_shift_auto())(1).get_etat() == ETATZERO ) ;
94 assert( (star.get_shift_auto())(3).get_etat() == ETATZERO ) ;
95
97 }
98 else {
101
102 Vector comp_shift(star.get_shift_auto()) ;
103 comp_shift.change_triad(mp_ns.get_bvect_cart()) ;
104 comp_shift.change_triad(mp.get_bvect_cart()) ;
105
106 assert( *(shift_comp.get_triad()) == *(comp_shift.get_triad()) ) ;
107
108 (shift_comp.set(1)).import( comp_shift(1) ) ;
109 (shift_comp.set(2)).import( comp_shift(2) ) ;
110 (shift_comp.set(3)).import( comp_shift(3) ) ;
111
113 }
114
115 // Conformal factor
116 // ----------------
117
118 if ( (star.get_confo_auto()).get_etat() == ETATZERO ) {
120 }
121 else {
125 }
126
127 //----------------------------------------------------
128 // Relaxation on lapconf_comp, shift_comp, confo_comp
129 //----------------------------------------------------
130
131 double relax_jm1 = 1. - relax ;
132
133 lapconf_comp = relax * lapconf_comp
134 + relax_jm1 * (hole_prev.lapconf_comp) ;
135
136 shift_comp = relax * shift_comp + relax_jm1 * (hole_prev.shift_comp) ;
137
138 confo_comp = relax * confo_comp + relax_jm1 * (hole_prev.confo_comp) ;
139
140
141 // Coordinates
142 // -----------
143 Scalar rr(mp) ;
144 rr = mp.r ;
145 rr.std_spectral_base() ;
146 Scalar st(mp) ;
147 st = mp.sint ;
148 st.std_spectral_base() ;
149 Scalar ct(mp) ;
150 ct = mp.cost ;
151 ct.std_spectral_base() ;
152 Scalar sp(mp) ;
153 sp = mp.sinp ;
154 sp.std_spectral_base() ;
155 Scalar cp(mp) ;
156 cp = mp.cosp ;
157 cp.std_spectral_base() ;
158
159 Vector ll(mp, CON, mp.get_bvect_cart()) ;
160 ll.set_etat_qcq() ;
161 ll.set(1) = st % cp ;
162 ll.set(2) = st % sp ;
163 ll.set(3) = ct ;
164 ll.std_spectral_base() ;
165
166
167 if (kerrschild) {
168
169 //----------------------------------------------
170 // Metric quantities from the analytic solution
171 //----------------------------------------------
172
173 lapconf_auto_bh = 1. / sqrt(1.+2.*mass/rr) ;
177
178 confo_auto_bh = 1. ;
180
181 shift_auto_bh = 2. * lapconf_auto_bh*lapconf_auto_bh*mass * ll / rr ;
184
185 //---------------------------------
186 // Derivative of metric quantities
187 //---------------------------------
188
190 for (int i=1; i<=3; i++)
192
194
196 for (int i=1; i<=3; i++) {
197 d_lapconf_auto_bh.set(i) = pow(lapconf_auto_bh,3.) * mass * ll(i)
198 / rr / rr ;
199 }
203
206
208 for (int i=1; i<=3; i++) {
209 for (int j=1; j<=3; j++) {
210 d_shift_auto_rs.set(i,j) = shift_auto_rs(j).deriv(i) ;
211 }
212 }
213
215
217 for (int i=1; i<=3; i++) {
218 for (int j=1; j<=3; j++) {
220 *lapconf_auto_bh*mass
221 * (flat.con()(i,j)
222 - 2.*lapconf_auto_bh*lapconf_auto_bh*(1.+mass/rr)
223 * ll(i) * ll(j))
224 / rr / rr ;
225 }
226 }
230
233
235 for (int i=1; i<=3; i++)
237
239
241 for (int i=1; i<=3; i++)
242 d_confo_auto_bh.set(i) = 0. ;
243
245
248
249 }
250 else { // Isotropic coordinates with the maximal slicing
251
252 // Sets C/M^2 for each case of the lapse boundary condition
253 // --------------------------------------------------------
254 double cc ;
255
256 if (bc_lapconf_nd) { // Neumann boundary condition
257 if (bc_lapconf_fs) { // First condition
258 // d(\alpha \psi)/dr = 0
259 // ---------------------
260 cc = 2. * (sqrt(13.) - 1.) / 3. ;
261 }
262 else { // Second condition
263 // d(\alpha \psi)/dr = (\alpha \psi)/(2 rah)
264 // -----------------------------------------
265 cc = 4. / 3. ;
266 }
267 }
268 else { // Dirichlet boundary condition
269 if (bc_lapconf_fs) { // First condition
270 // (\alpha \psi) = 1/2
271 // -------------------
272 cout << "!!!!! WARNING: Not yet prepared !!!!!" << endl ;
273 abort() ;
274 }
275 else { // Second condition
276 // (\alpha \psi) = 1/sqrt(2.) \psi_KS
277 // -----------------------------------
278 cout << "!!!!! WARNING: Not yet prepared !!!!!" << endl ;
279 abort() ;
280 // cc = 2. * sqrt(2.) ;
281 }
282 }
283
284 //----------------------------------------------
285 // Metric quantities from the analytic solution
286 //----------------------------------------------
287
288 Scalar r_are(mp) ;
290 r_are.std_spectral_base() ;
291
292 lapconf_auto_bh = sqrt(1. - 2.*mass/r_are/rr
293 + cc*cc*pow(mass/r_are/rr, 4.)) * sqrt(r_are) ;
297
298 confo_auto_bh = sqrt(r_are) ;
302
303 shift_auto_bh = mass * mass * cc * ll / rr / rr / pow(r_are, 3.) ;
306 for (int i=1; i<=3; i++)
308
309 //---------------------------------
310 // Derivative of metric quantities
311 //---------------------------------
312
314 for (int i=1; i<=3; i++)
316
318
320 for (int i=1; i<=3; i++) {
321 d_lapconf_auto_bh.set(i) = sqrt(r_are)
322 * (mass/r_are/rr - 2.*cc*cc*pow(mass/r_are/rr,4.))
323 * ll(i) / rr
324 + 0.5 * sqrt(r_are)
325 * (sqrt(1. - 2.*mass/r_are/rr + cc*cc*pow(mass/r_are/rr,4.))-1.)
326 * sqrt(1. - 2.*mass/r_are/rr + cc*cc*pow(mass/r_are/rr,4.))
327 * ll(i) / rr ;
328 }
332
336 for (int i=1; i<=3; i++)
338
340 for (int i=1; i<=3; i++) {
341 for (int j=1; j<=3; j++) {
342 d_shift_auto_rs.set(i,j) = shift_auto_rs(j).deriv(i) ;
343 }
344 }
345
347
349 for (int i=1; i<=3; i++) {
350 for (int j=1; j<=3; j++) {
351 d_shift_auto_bh.set(i,j) =
352 mass*mass*cc*(flat.con()(i,j)
353 -3.*sqrt(1. - 2.*mass/r_are/rr
354 +cc*cc*pow(mass/r_are/rr,4.))
355 *ll(i)*ll(j))
356 / pow(r_are*rr,3.) ;
357 }
358 }
362
366 for (int i=1; i<=3; i++) {
367 for (int j=1; j<=3; j++) {
368 d_shift_auto.set(i,j).raccord(1) ;
369 }
370 }
371
373 for (int i=1; i<=3; i++)
375
377
379 for (int i=1; i<=3; i++) {
380 d_confo_auto_bh.set(i) = 0.5*sqrt(r_are)
381 *(sqrt(1. - 2.*mass/r_are/rr +cc*cc*pow(mass/r_are/rr,4.)) - 1.)
382 * ll(i) / rr ;
383 }
387
391 for (int i=1; i<=3; i++)
392 d_confo_auto.set(i).raccord(1) ;
393
394 }
395
396 //-------------------------
397 // Total metric quantities
398 //-------------------------
399
404
409
413 for (int i=1; i<=3; i++) {
414 shift_auto.set(i).raccord(1) ;
415 }
416
420 for (int i=1; i<=3; i++) {
421 shift_tot.set(i).raccord(1) ;
422 }
423
428
432 confo_tot.raccord(1) ;
433
436
439
440 // The derived quantities are obsolete
441 // -----------------------------------
442
443 del_deriv() ;
444
445}
446}
const Scalar r_coord(bool neumann, bool first) const
Expresses the areal radial coordinate by that in spatially isotropic coordinates.
Map & mp
Mapping associated with the black hole.
Definition blackhole.h:80
Metric_flat flat
Flat metric defined on the mapping (Spherical components with respect to the mapping of the black hol...
Definition blackhole.h:143
bool kerrschild
true for a Kerr-Schild background, false for a conformally flat background
Definition blackhole.h:85
double mass_bh
Gravitational mass of BH.
Definition blackhole.h:88
Class for black holes in black hole-neutron star binaries.
Definition hole_bhns.h:65
Scalar confo_auto
Conformal factor generated by the black hole.
Definition hole_bhns.h:163
Scalar lapconf_auto
Lapconf function generated by the black hole.
Definition hole_bhns.h:95
Vector shift_tot
Total shift vector ;.
Definition hole_bhns.h:138
Vector d_lapconf_auto
Derivative of the lapconf function generated by the black hole.
Definition hole_bhns.h:120
Tensor d_shift_auto
Derivative of the shift vector generated by the black hole.
Definition hole_bhns.h:151
Vector d_lapconf_auto_rs
Derivative of the part of the lapconf function from the numerical computation.
Definition hole_bhns.h:112
Scalar confo_auto_bh
Part of the conformal factor from the analytic background.
Definition hole_bhns.h:160
Tensor d_shift_auto_rs
Derivative of the part of the shift vector from the numerical computation.
Definition hole_bhns.h:143
Vector shift_auto_rs
Part of the shift vector from the numerical computation.
Definition hole_bhns.h:126
virtual void del_deriv() const
Deletes all the derived quantities.
Definition hole_bhns.C:392
Tensor d_shift_auto_bh
Derivative of the part of the shift vector from the analytic background.
Definition hole_bhns.h:148
bool bc_lapconf_fs
true for the first type BC for the lapconf function, false for the second type BH
Definition hole_bhns.h:78
void update_metric_bhns(const Star_bhns &star, const Hole_bhns &hole_prev, double relax)
Computes metric coefficients from known potentials with relaxation when the companion is a black hole...
Vector d_lapconf_auto_bh
Derivative of the part of the lapconf function from the analytic background.
Definition hole_bhns.h:117
Vector d_confo_auto_rs
Derivative of the part of the conformal factor from the numerical computation.
Definition hole_bhns.h:174
Scalar confo_auto_rs
Part of the conformal factor from the numerical computation.
Definition hole_bhns.h:157
Scalar lapse_tot
Total lapse function.
Definition hole_bhns.h:107
Scalar lapse_auto
Lapse function of the "black hole" part.
Definition hole_bhns.h:104
Vector shift_auto_bh
Part of the shift vector from the analytic background.
Definition hole_bhns.h:129
Scalar confo_comp
Conformal factor generated by the companion star.
Definition hole_bhns.h:166
Scalar lapconf_comp
Lapconf function generated by the companion star.
Definition hole_bhns.h:98
Scalar lapconf_auto_bh
Part of the lapconf function from the analytic background.
Definition hole_bhns.h:92
Vector d_confo_auto
Derivative of the conformal factor generated by the black hole.
Definition hole_bhns.h:182
Vector shift_auto
Shift vector generated by the black hole.
Definition hole_bhns.h:132
Scalar lapconf_auto_rs
Part of the lapconf function from the numerical computation.
Definition hole_bhns.h:89
bool bc_lapconf_nd
true for the Neumann type BC for the lapconf function, false for the Dirichlet type BH
Definition hole_bhns.h:73
Scalar lapconf_tot
Total lapconf function.
Definition hole_bhns.h:101
Vector d_confo_auto_bh
Derivative of the part of the conformal factor from the analytic background.
Definition hole_bhns.h:179
Vector shift_comp
Shift vector generated by the companion star.
Definition hole_bhns.h:135
Scalar confo_tot
Total conformal factor.
Definition hole_bhns.h:169
Base class for coordinate mappings.
Definition map.h:670
Coord cosp
Definition map.h:724
const Base_vect_cart & get_bvect_cart() const
Returns the Cartesian basis associated with the coordinates (x,y,z) of the mapping,...
Definition map.h:791
Coord sint
Definition map.h:721
Coord r
r coordinate centered on the grid
Definition map.h:718
Coord sinp
Definition map.h:723
Coord cost
Definition map.h:722
virtual const Sym_tensor & con() const
Read-only access to the contravariant representation.
Tensor field of valence 0 (or component of a tensorial field).
Definition scalar.h:387
void raccord(int n)
Performs the matching of the nucleus with respect to the first shell.
const Scalar & deriv(int i) const
Returns of *this , where .
virtual void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition scalar.C:353
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
virtual void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition scalar.C:324
void import(const Scalar &ci)
Assignment to another Scalar defined on a different mapping.
Class for stars in black hole-neutron star binaries.
Definition star_bhns.h:59
const Scalar & get_confo_auto() const
Returns the part of the conformal factor generated by the star.
Definition star_bhns.h:383
const Vector & get_shift_auto() const
Returns the part of the shift vector generated by the star.
Definition star_bhns.h:364
const Scalar & get_lapconf_auto() const
Returns the part of the lapconf function generated by the star.
Definition star_bhns.h:339
const Map & get_mp() const
Returns the mapping.
Definition star.h:355
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
virtual void change_triad(const Base_vect &)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
Scalar & set(int)
Read/write access to a component.
Definition vector.C:296
Cmp sqrt(const Cmp &)
Square root.
Definition cmp_math.C:220
Cmp pow(const Cmp &, int)
Power .
Definition cmp_math.C:348
void set_triad(const Base_vect &new_triad)
Assigns a new vectorial basis (triad) of decomposition.
Definition tensor.C:519
virtual void inc_dzpuis(int inc=1)
Increases by inc units the value of dzpuis and changes accordingly the values in the compactified ext...
Definition tensor.C:816
const Base_vect * get_triad() const
Returns the vectorial basis (triad) on which the components are defined.
Definition tensor.h:866
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
virtual void set_etat_zero()
Sets the logical state of all components to ETATZERO (zero state).
Definition tensor.C:497
Scalar & set(const Itbl &ind)
Returns the value of a component (read/write version).
Definition tensor.C:654
virtual void std_spectral_base()
Sets the standard spectal bases of decomposition for each component.
Definition tensor.C:926
Lorene prototypes.
Definition app_hor.h:64
Standard units of space, time and mass.