28char blackhole_extr_curv_C[] =
"$Header: /cvsroot/Lorene/C++/Source/Black_hole/blackhole_extr_curv.C,v 1.4 2014/10/13 08:52:46 j_novak Exp $" ;
58#include "utilitaires.h"
107 for (
int i=1; i<=3; i++) {
108 for (
int j=1; j<=3; j++) {
111 - 2. * divshift *
flat.
con()(i,j) / 3. ;
118 lapse_bh = 1. /
sqrt(1. + 2. * mass / rr) ;
126 for (
int i=1; i<=3; i++) {
127 for (
int j=1; j<=3; j++) {
128 curv_taij.
set(i,j) = -2. * lapse_bh * lapse_bh * mass
130 - 2. * ll(i) * ll(j) * divshift / 3.) / rr ;
139 for (
int i=1; i<=3; i++) {
140 for (
int j=1; j<=3; j++) {
141 resi_taij.
set(i,j) = 2. * lapse_bh * lapse_bh * mass
144 - lapse_bh*lapse_bh*(9.+14.*mass/rr)*ll(i)*ll(j) )
155 * (flat_taij + curv_taij + resi_taij) /
lapconf ;
163 for (
int i=1; i<=3; i++) {
164 for (
int j=1; j<=3; j++) {
165 taij_bh.
set(i,j) = 2.*
pow(lapse_bh,6.)*mass*(2.+3.*mass/rr)
166 *( (1.+2.*mass/rr) *
flat.
con()(i,j)
167 - (3.+2.*mass/rr) * ll(i) * ll(j) )
215 for (
int i=1; i<=3; i++) {
216 for (
int j=1; j<=3; j++) {
223 - 2. * divshift *
flat.
cov()(i,j) / 3. ;
232 for (
int i=1; i<=3; i++) {
233 for (
int j=1; j<=3; j++) {
234 curv_dshift.
set(i,j) = 2. * mass
235 *( ll(j) *( ll(1)*(
shift_rs(1).deriv(i))
238 + ll(i) *( ll(1)*(
shift_rs(1).deriv(j))
241 - 2. * divshift * ll(i) * ll(j) / 3. ) / rr ;
250 for (
int i=1; i<=3; i++) {
251 for (
int j=1; j<=3; j++) {
252 tmp1.
set(i,j) = 2. * mass
253 *(ll(j)*( (
flat.
cov()(i,1)+2.*mass*ll(i)*ll(1)/rr)
255 + (
flat.
cov()(i,2)+2.*mass*ll(i)*ll(2)/rr)
257 + (
flat.
cov()(i,3)+2.*mass*ll(i)*ll(3)/rr)
260 + ll(i)*( (
flat.
cov()(j,1)+2.*mass*ll(j)*ll(1)/rr)
262 + (
flat.
cov()(j,2)+2.*mass*ll(j)*ll(2)/rr)
264 + (
flat.
cov()(j,3)+2.*mass*ll(j)*ll(3)/rr)
275 for (
int i=1; i<=3; i++) {
276 for (
int j=1; j<=3; j++) {
277 tmp2.
set(i,j) = 2. * mass * lapse_bh * lapse_bh
280 - (9.+28.*mass/rr+24.*mass*mass/rr/rr)*ll(i)*ll(j))
291 * (flat_dshift + curv_dshift + tmp1 + tmp2) /
lapconf ;
299 for (
int i=1; i<=3; i++) {
300 for (
int j=1; j<=3; j++) {
301 taij_bh_down.
set(i,j) = 2.*
pow(lapse_bh,4.)*mass*(2.+3.*mass/rr)
340 taij_quad_rsrs = 0. ;
342 for (
int i=1; i<=3; i++) {
343 for (
int j=1; j<=3; j++) {
344 taij_quad_rsrs += taij_rs_down(i,j) *
taij_rs(i,j) ;
350 taij_quad_rsbh1 = 0. ;
352 for (
int i=1; i<=3; i++) {
353 for (
int j=1; j<=3; j++) {
354 taij_quad_rsbh1 += taij_rs_down(i,j) * taij_bh(i,j) ;
360 taij_quad_rsbh2 = 0. ;
362 for (
int i=1; i<=3; i++) {
363 for (
int j=1; j<=3; j++) {
364 taij_quad_rsbh2 += taij_bh_down(i,j) *
taij_rs(i,j) ;
369 taij_quad_rs = taij_quad_rsrs + taij_quad_rsbh1 + taij_quad_rsbh2 ;
373 taij_quad_bh = 8.*
pow(lapse_bh,10.)*mass*mass*(2.+3.*mass/rr)
407 for (
int i=1; i<=3; i++) {
408 for (
int j=1; j<=3; j++) {
410 - 2. * divshift *
flat.
con()(i,j) / 3. ;
428 for (
int i=1; i<=3; i++) {
429 for (
int j=1; j<=3; j++) {
436 - 2. * divshift *
flat.
cov()(i,j) / 3. ;
443 for (
int i=1; i<=3; i++) {
444 for (
int j=1; j<=3; j++) {
445 taij_down.
set(i,j) = 0.5 *
pow(
confo, 7.) * flat_dshift(i,j)
455 for (
int i=1; i<=3; i++) {
456 for (
int j=1; j<=3; j++) {
Vector shift_rs
Part of the shift vector from the numerical computation.
Scalar taij_quad
Part of the scalar generated by .
Scalar taij_quad_rs
Part of the scalar.
Sym_tensor taij
Trace of the extrinsic curvature.
Sym_tensor taij_rs
Part of the extrinsic curvature tensor.
Map & mp
Mapping associated with the black hole.
Scalar lapconf
A function (lapse function * conformal factor) lapconf generated by the black hole.
Vector shift
Shift vector generated by the black hole.
Metric_flat flat
Flat metric defined on the mapping (Spherical components with respect to the mapping of the black hol...
bool kerrschild
true for a Kerr-Schild background, false for a conformally flat background
Scalar confo
Conformal factor generated by the black hole.
double mass_bh
Gravitational mass of BH.
void extr_curv_bh()
Computes taij , taij_quad from shift , lapse , confo .
const Base_vect_cart & get_bvect_cart() const
Returns the Cartesian basis associated with the coordinates (x,y,z) of the mapping,...
Coord r
r coordinate centered on the grid
virtual const Sym_tensor & con() const
Read-only access to the contravariant representation.
virtual const Sym_tensor & cov() const
Read-only access to the covariant representation.
Tensor field of valence 0 (or component of a tensorial field).
void raccord(int n)
Performs the matching of the nucleus with respect to the first shell.
virtual void std_spectral_base()
Sets the spectral bases of the Valeur va to the standard ones for a scalar field.
virtual void inc_dzpuis(int inc=1)
Increases by inc units the value of dzpuis and changes accordingly the values of the Scalar in the co...
Class intended to describe valence-2 symmetric tensors.
Tensor field of valence 1.
virtual void std_spectral_base()
Sets the standard spectal bases of decomposition for each component.
Scalar & set(int)
Read/write access to a component.
Cmp sqrt(const Cmp &)
Square root.
Cmp pow(const Cmp &, int)
Power .
virtual void inc_dzpuis(int inc=1)
Increases by inc units the value of dzpuis and changes accordingly the values in the compactified ext...
void annule_domain(int l)
Sets the Tensor to zero in a given domain.
virtual void set_etat_qcq()
Sets the logical state of all components to ETATQCQ (ordinary state).
Scalar & set(const Itbl &ind)
Returns the value of a component (read/write version).
virtual void std_spectral_base()
Sets the standard spectal bases of decomposition for each component.
Standard units of space, time and mass.