|
LORENE
|
Spheroidal 2-surfaces embedded in a time-slice of the 3+1 formalism. More...
#include <spheroid.h>
Public Member Functions | |
| Spheroid (const Map_af &map, double radius) | |
| The delta tensorial fields linked to Christoffel symbols. | |
| Spheroid (const Scalar &h_in, const Metric &gamij, const Sym_tensor &Kij) | |
Constructor of a spheroid embedded in a 3-slice (Time_slice ) of 3+1 formalism. | |
| Spheroid (const Spheroid &) | |
| Copy constructor. | |
| Spheroid (FILE *) | |
Constructor from a file (see sauve(FILE*) ). | |
| virtual | ~Spheroid () |
| Destructor. | |
| void | operator= (const Spheroid &) |
| Assignment to another Spheroid. | |
| void | set_ephi (const Scalar &) |
| Assigns the conformal Killing vector field to phi. | |
| const Scalar & | get_hsurf () const |
Returns the field h_surf. | |
| const Metric & | get_qab () const |
Returns the metric ![]() | |
| const Scalar & | get_ricci () const |
Returns the 2-ricci scalar ![]() | |
| const Sym_tensor & | get_hh () const |
Returns the symmetric tensor ![]() | |
| const Sym_tensor & | get_qq () const |
returns the 3-d degenerate 2-metric ![]() | |
| const Tensor & | get_proj () const |
returns the 3-d projector on 2-surface ![]() | |
| const Tensor & | get_jac2d () const |
returns the 2-d jacobian of coordinate transformation ![]() | |
| const Scalar & | get_trk () const |
Returns the trace K on the 2-surface. | |
| const Vector & | get_ll () const |
Returns the vector ![]() | |
| const Vector & | get_ss () const |
Returns the vector ![]() | |
| const Vector & | get_ephi () const |
| Returns the conformal Killing symmetry vector on the 2-surface. | |
| const Sym_tensor & | get_jj () const |
Returns the symmetric tensor ![]() | |
| const Scalar & | get_fff () const |
Returns the normalization scalar F. | |
| const Scalar & | get_ggg () const |
Returns the normalization scalar G. | |
| bool | get_issphere () const |
| Returns the flag saying whether or not the horizon is geometrically round. | |
| Scalar & | set_hsurf () |
Sets the field h_surf. | |
| Metric & | set_qab () |
Sets the modified metric (non degenerated) ![]() | |
| Scalar & | set_ricci () |
Sets the 2-Ricci scalar ![]() | |
| Sym_tensor & | set_qq () |
Sets the degenerated metric ![]() | |
| Tensor & | set_proj () |
Sets the projector ![]() | |
| Sym_tensor & | set_hh () |
Sets the symmetric tensor ![]() | |
| Scalar & | set_trk () |
Sets the trace K on the 2-surface. | |
| Vector & | set_ll () |
Sets the vector ![]() | |
| Vector & | set_ss () |
Sets the vector ![]() | |
| Sym_tensor & | set_jj () |
Sets the symmetric tensor ![]() | |
| Scalar & | set_fff () |
Sets the normalization factor F. | |
| Scalar & | set_ggg () |
Sets the normalization factor G. | |
| bool | set_issphere () |
| Sets the boolean linked to geometrical shape of the horizon. | |
| void | update_from_tslice (const Metric &gamij, const Sym_tensor &Kij) |
| Updates from the 3-slice data. | |
| const Scalar & | sqrt_q () const |
| Computes the normal vector field to the 2-surface. | |
| double | area () const |
| Computes the area of the 2-surface. | |
| double | angu_mom () const |
| Computes the angular momentum with respect to a divergence-free vector field tangent to the 2-surface. | |
| double | mass () const |
Computes the mass as defined from the calculus of angular momentum, done with respect to a divergence free tangent vector field ![]() | |
| double | multipole_mass (const int order) const |
| Computes the mass multipole of a given order for the spheroid, assumed to be spherical. | |
| double | multipole_angu (const int order) const |
| Computes the angular multipole of a given order for the spheroid, assumed to be spherical. | |
| double | epsilon_A_minus_one () const |
| Computation of the refined Penrose parameter for axisymmetric spacetimes, and its difference wrt one. | |
| double | epsilon_P_minus_one () const |
| Computation of the classical Penrose parameter, and its difference wrt one. | |
| const Scalar & | theta_plus () const |
Computes the outgoing null expansion ![]() | |
| const Scalar & | theta_minus () const |
Computes the ingoing null expansion ![]() | |
| const Sym_tensor & | shear () const |
Computes the shear of the 2-surface ![]() | |
| Tensor | derive_cov2dflat (const Tensor &uu) const |
| Computes the round covariant derivative on the spheroid. | |
| const Tensor & | delta () const |
| Computes the delta coefficients for covariant derivative. | |
| Tensor | derive_cov2d (const Tensor &uu) const |
| Computes the total covariant derivative on the spheroid. | |
| virtual void | sauve (FILE *) const |
| Save in a file. | |
Protected Member Functions | |
| virtual void | del_deriv () const |
| Deletes all the derived quantities. | |
| void | set_der_0x0 () const |
Sets to 0x0 all the pointers on derived quantities. | |
Protected Attributes | |
| Scalar | h_surf |
The location of the 2-surface as r = h_surf ![]() | |
| Tensor | jac2d |
| The jacobian of the adaptation of coordinates (contravariant/covariant representation). | |
| Tensor | proj |
| The 3-d projector on the 2-surface (contravariant-covariant form). | |
| Sym_tensor | |
| The 3-d covariant degenerated 2-metric on the surface. | |
| Vector | ss |
| The adapted normal vector field to spheroid in the 3-slice. | |
| Vector | ephi |
The conformal Killing vector field on the 2-surface (set to by default to the axial vector associated with coordinate ![]() | |
| Metric | qab |
| Scalar | ricci |
Induced metric on the 2-surface ![]() | |
| Sym_tensor | hh |
| The ricci scalar on the 2-surface. | |
| Scalar | trk |
| Trace K of the extrinsic curvature of the 3-slice. | |
| Vector | ll |
| Normal-tangent component of the extrinsic curvature of the 3-slice. | |
| Sym_tensor | jj |
| Tangent components of the extrinsic curvature of the 3-slice. | |
| Scalar | fff |
| Normalization function for the outgoing null vector l. | |
| Scalar | ggg |
| Normalization function for the ingoing null vector k. | |
| Scalar | zeta |
| bool | issphere |
| Flag to know whether the horizon is geometrically round or distorted. | |
| Scalar * | p_sqrt_q |
Surface element ![]() | |
| double * | p_area |
| The area of the 2-surface. | |
| double * | p_angu_mom |
| The angular momentum. | |
| double * | p_mass |
| Mass defined from angular momentum. | |
| double * | p_multipole_mass |
| Mass multipole for the spheroid. | |
| double * | p_multipole_angu |
| Angular momentum multipole for the spheroid. | |
| double * | p_epsilon_A_minus_one |
| double * | p_epsilon_P_minus_one |
| Refined Penrose parameter, difference wrt one. | |
| Scalar * | p_theta_plus |
| Classical Penrose parameter, difference wrt one. | |
| Scalar * | p_theta_minus |
| Null ingoing expansion. | |
| Sym_tensor * | p_shear |
| The shear tensor. | |
| Tensor * | p_delta |
Friends | |
| ostream & | operator<< (ostream &, const Spheroid &) |
| Display. | |
Spheroidal 2-surfaces embedded in a time-slice of the 3+1 formalism.
()
Definition at line 84 of file spheroid.h.
| Lorene::Spheroid::Spheroid | ( | const Map_af & | map, |
| double | radius ) |
The delta tensorial fields linked to Christoffel symbols.
Standard constructor. The input mapping must be defined on a mono-domain angular grid (see Mg3d::get_angu_mono_domain() for details).
Definition at line 94 of file spheroid.C.
References ephi, fff, Lorene::flat_met_spher(), Lorene::get_bvect_spher(), ggg, h_surf, hh, issphere, jac2d, jj, ll, proj, qq, ricci, set_der_0x0(), ss, and trk.
| Lorene::Spheroid::Spheroid | ( | const Scalar & | h_in, |
| const Metric & | gamij, | ||
| const Sym_tensor & | Kij ) |
Constructor of a spheroid embedded in a 3-slice (Time_slice ) of 3+1 formalism.
This is done from the Time_slice data.
| h_in | : the location of the surface r = h_in (WARNING:must be defined on a mono-domain angular grid) |
| gamij | : the 3-metric on the 3-slice |
| Kij | : the extrinsic curvature of the 3-slice (covariant representation) |
Definition at line 150 of file spheroid.C.
References Lorene::Scalar::allocate_all(), Lorene::Tensor::allocate_all(), Lorene::Tensor::annule_domain(), Lorene::Scalar::annule_hard(), area(), Lorene::Valeur::c_cf, Lorene::Valeur::coef(), Lorene::Metric::con(), Lorene::contract(), Lorene::Metric::cov(), Lorene::Scalar::dec_dzpuis(), del_deriv(), Lorene::Scalar::derive_con(), Lorene::Scalar::derive_cov(), Lorene::Tensor::derive_cov(), Lorene::Tensor::down(), ephi, fff, Lorene::flat_met_spher(), Lorene::get_bvect_spher(), Lorene::Metric::get_mp(), Lorene::Tensor::get_mp(), Lorene::Mg3d::get_np(), Lorene::Mg3d::get_nr(), Lorene::Mg3d::get_nt(), Lorene::Mg3d::get_nzone(), Lorene::Scalar::get_spectral_base(), Lorene::Scalar::get_spectral_va(), ggg, Lorene::Base_val::give_quant_numbers(), h_surf, hh, issphere, jac2d, jj, ll, Lorene::Map(), Lorene::Scalar::mult_cost(), Lorene::Scalar::mult_r(), Lorene::Scalar::mult_sint(), proj, qq, Lorene::Metric::ricci(), ricci, Lorene::Itbl::set(), Lorene::Tensor::set(), Lorene::Vector::set(), set_der_0x0(), Lorene::Tensor::set_etat_zero(), Lorene::Scalar::set_grid_point(), Lorene::Scalar::set_spectral_va(), Lorene::sqrt(), sqrt_q(), Lorene::Scalar::srdsdt(), Lorene::Scalar::srstdsdp(), ss, Lorene::Scalar::std_spectral_base(), Lorene::Tensor::std_spectral_base(), Lorene::Vector::std_spectral_base(), Lorene::Tensor::trace(), trk, Lorene::Tensor::up(), Lorene::Map_radial::val_lx_jk(), Lorene::Valeur::val_point_jk(), and Lorene::Valeur::ylm().
| Lorene::Spheroid::Spheroid | ( | const Spheroid & | sph_in | ) |
| Lorene::Spheroid::Spheroid | ( | FILE * | ) |
Constructor from a file (see sauve(FILE*) ).
References Spheroid().
|
virtual |
| double Lorene::Spheroid::angu_mom | ( | ) | const |
Computes the angular momentum with respect to a divergence-free vector field tangent to the 2-surface.
This is defined as
![\[{\cal J} = \int h^2 L_i \phi^i \sqrt{\det q_{ab}} \sin \theta {\rm d}\theta
{\rm d}\varphi \]](form_748.png)
| phi | : the divergence-free vector field ![]() |
Definition at line 721 of file spheroid.C.
References Lorene::contract(), get_ephi(), h_surf, Lorene::Map_af::integrale_surface(), jac2d, ll, p_angu_mom, Lorene::phi, and sqrt_q().
| double Lorene::Spheroid::area | ( | ) | const |
Computes the area of the 2-surface.
This is defined as
![\[{\cal A} = \int h^2 \sqrt{\det q_{ab}} \sin \theta {\rm d}\theta
{\rm d}\varphi \]](form_747.png)
Definition at line 708 of file spheroid.C.
References h_surf, Lorene::Map_af::integrale_surface(), p_area, and sqrt_q().
|
protectedvirtual |
Deletes all the derived quantities.
Definition at line 650 of file spheroid.C.
References p_angu_mom, p_area, p_epsilon_P_minus_one, p_mass, p_multipole_angu, p_multipole_mass, p_shear, p_sqrt_q, p_theta_minus, p_theta_plus, and set_der_0x0().
| const Tensor & Lorene::Spheroid::delta | ( | ) | const |
Computes the delta coefficients for covariant derivative.
Definition at line 1203 of file spheroid.C.
References derive_cov2dflat(), Lorene::Tensor::set(), Lorene::Tensor::set_etat_zero(), and Lorene::Tensor::set_index_type().
Computes the total covariant derivative on the spheroid.
Definition at line 1240 of file spheroid.C.
References Lorene::contract(), delta(), derive_cov2dflat(), Lorene::Tensor::get_index_type(), Lorene::Tensor::get_valence(), and Lorene::y.
Computes the round covariant derivative on the spheroid.
Definition at line 954 of file spheroid.C.
References Lorene::Scalar::div_r_dzpuis(), Lorene::Scalar::div_tant(), Lorene::Tensor::get_index_type(), Lorene::Tensor::get_mp(), Lorene::Tensor::get_n_comp(), Lorene::Tensor::get_triad(), Lorene::Tensor::get_valence(), Lorene::Tensor::indices(), Lorene::Itbl::set(), Lorene::Tensor::set(), Lorene::srdsdt(), Lorene::srstdsdp(), Lorene::Tensor_sym::sym_index1(), and Lorene::Tensor_sym::sym_index2().
| double Lorene::Spheroid::epsilon_A_minus_one | ( | ) | const |
Computation of the refined Penrose parameter for axisymmetric spacetimes, and its difference wrt one.
Definition at line 864 of file spheroid.C.
References angu_mom(), area(), mass(), Lorene::pow(), and Lorene::sqrt().
| double Lorene::Spheroid::epsilon_P_minus_one | ( | ) | const |
Computation of the classical Penrose parameter, and its difference wrt one.
To use in replacement of epsilon_A_minus_one when the computed spacetime is not axisymmetric.
Definition at line 876 of file spheroid.C.
References angu_mom(), area(), mass(), p_epsilon_P_minus_one, and Lorene::pow().
|
inline |
Returns the conformal Killing symmetry vector on the 2-surface.
Definition at line 253 of file spheroid.h.
References ephi.
|
inline |
|
inline |
|
inline |
|
inline |
|
inline |
Returns the flag saying whether or not the horizon is geometrically round.
Definition at line 265 of file spheroid.h.
References issphere.
|
inline |
returns the 2-d jacobian of coordinate transformation 
Definition at line 241 of file spheroid.h.
References jac2d.
|
inline |
|
inline |
|
inline |
|
inline |
Returns the metric 
Definition at line 226 of file spheroid.h.
|
inline |
|
inline |
|
inline |
|
inline |
| double Lorene::Spheroid::mass | ( | ) | const |
Computes the mass as defined from the calculus of angular momentum, done with respect to a divergence free tangent vector field 
Spheroid has to be a real sphere (flag issphere true), of constant radius
![\[R_{s} \]](form_751.png)
. defined as
![\[ M = \frac{1}{2 R_{s}} \sqrt{R_{s}^{4} + 4{\cal J}^{2}} \]](form_752.png)
Definition at line 736 of file spheroid.C.
References angu_mom(), area(), p_mass, and Lorene::sqrt().
| double Lorene::Spheroid::multipole_angu | ( | const int | order | ) | const |
Computes the angular multipole of a given order for the spheroid, assumed to be spherical.

Definition at line 801 of file spheroid.C.
References area(), Lorene::contract(), get_ephi(), h_surf, Lorene::Map_af::integrale_surface(), jac2d, ll, p_multipole_angu, Lorene::phi, Lorene::Scalar::set_spectral_va(), Lorene::sqrt(), sqrt_q(), Lorene::Scalar::std_spectral_base(), and Lorene::Valeur::ylm().
| double Lorene::Spheroid::multipole_mass | ( | const int | order | ) | const |
Computes the mass multipole of a given order for the spheroid, assumed to be spherical.
WARNING: For technical reasons, only even orders are supported by the code.
Definition at line 748 of file spheroid.C.
References area(), get_ricci(), h_surf, Lorene::Map_af::integrale_surface(), mass(), p_multipole_mass, Lorene::Scalar::set_spectral_va(), Lorene::sqrt(), sqrt_q(), Lorene::Scalar::std_spectral_base(), and Lorene::Valeur::ylm().
| void Lorene::Spheroid::operator= | ( | const Spheroid & | sph_in | ) |
|
virtual |
Save in a file.
Definition at line 1186 of file spheroid.C.
|
protected |
Sets to 0x0 all the pointers on derived quantities.
Definition at line 666 of file spheroid.C.
References p_angu_mom, p_area, p_epsilon_P_minus_one, p_mass, p_multipole_angu, p_multipole_mass, p_shear, p_sqrt_q, p_theta_minus, and p_theta_plus.
| void Lorene::Spheroid::set_ephi | ( | const Scalar & | ) |
Assigns the conformal Killing vector field to phi.
|
inline |
Sets the normalization factor F.
Definition at line 298 of file spheroid.h.
References del_deriv(), and fff.
|
inline |
Sets the normalization factor G.
Definition at line 301 of file spheroid.h.
References del_deriv(), and ggg.
|
inline |
Sets the symmetric tensor 
Definition at line 283 of file spheroid.h.
References del_deriv(), and hh.
|
inline |
Sets the field h_surf.
Definition at line 268 of file spheroid.h.
References del_deriv(), and h_surf.
|
inline |
Sets the boolean linked to geometrical shape of the horizon.
Definition at line 304 of file spheroid.h.
References del_deriv(), and issphere.
|
inline |
Sets the symmetric tensor 
Definition at line 295 of file spheroid.h.
References del_deriv(), and jj.
|
inline |
|
inline |
|
inline |
Sets the modified metric (non degenerated) 
Definition at line 271 of file spheroid.h.
References del_deriv().
|
inline |
Sets the degenerated metric 
Definition at line 277 of file spheroid.h.
References del_deriv(), and qq.
|
inline |
Sets the 2-Ricci scalar 
Definition at line 274 of file spheroid.h.
References del_deriv(), and ricci.
|
inline |
|
inline |
Sets the trace K on the 2-surface.
Definition at line 286 of file spheroid.h.
References del_deriv(), and trk.
| const Sym_tensor & Lorene::Spheroid::shear | ( | ) | const |
| const Scalar & Lorene::Spheroid::sqrt_q | ( | ) | const |
Computes the normal vector field to the 2-surface.
Computes the square root of the determinant of 
Definition at line 693 of file spheroid.C.
References get_qq(), p_sqrt_q, and Lorene::sqrt().
| const Scalar & Lorene::Spheroid::theta_minus | ( | ) | const |
Computes the ingoing null expansion 
Definition at line 909 of file spheroid.C.
References ggg, hh, jj, and p_theta_minus.
| const Scalar & Lorene::Spheroid::theta_plus | ( | ) | const |
Computes the outgoing null expansion 
Definition at line 889 of file spheroid.C.
References fff, hh, jj, and p_theta_plus.
| void Lorene::Spheroid::update_from_tslice | ( | const Metric & | gamij, |
| const Sym_tensor & | Kij ) |
Updates from the 3-slice data.
|
friend |
Display.
References Spheroid().
|
protected |
The conformal Killing vector field on the 2-surface (set to by default to the axial vector associated with coordinate 
Definition at line 113 of file spheroid.h.
|
protected |
Normalization function for the outgoing null vector l.
Definition at line 138 of file spheroid.h.
|
protected |
Normalization function for the ingoing null vector k.
Definition at line 143 of file spheroid.h.
|
protected |
The location of the 2-surface as r = h_surf 
Definition at line 91 of file spheroid.h.
|
protected |
The ricci scalar on the 2-surface.
Extrinsic curvature of the 2-surface in the 3-slice. 
Definition at line 122 of file spheroid.h.
|
protected |
Flag to know whether the horizon is geometrically round or distorted.
Definition at line 151 of file spheroid.h.
|
protected |
The jacobian of the adaptation of coordinates (contravariant/covariant representation).
Definition at line 96 of file spheroid.h.
|
protected |
Tangent components of the extrinsic curvature of the 3-slice.

Definition at line 134 of file spheroid.h.
|
protected |
Normal-tangent component of the extrinsic curvature of the 3-slice.

Definition at line 129 of file spheroid.h.
|
mutableprotected |
The angular momentum.
Definition at line 159 of file spheroid.h.
|
mutableprotected |
The area of the 2-surface.
Definition at line 158 of file spheroid.h.
|
mutableprotected |
Definition at line 168 of file spheroid.h.
|
mutableprotected |
Definition at line 163 of file spheroid.h.
|
mutableprotected |
Refined Penrose parameter, difference wrt one.
Definition at line 164 of file spheroid.h.
|
mutableprotected |
Mass defined from angular momentum.
Definition at line 160 of file spheroid.h.
|
mutableprotected |
Angular momentum multipole for the spheroid.
Definition at line 162 of file spheroid.h.
|
mutableprotected |
Mass multipole for the spheroid.
Definition at line 161 of file spheroid.h.
|
mutableprotected |
The shear tensor.
Definition at line 167 of file spheroid.h.
|
mutableprotected |
Surface element 
Definition at line 157 of file spheroid.h.
|
mutableprotected |
Null ingoing expansion.
Definition at line 166 of file spheroid.h.
|
mutableprotected |
Classical Penrose parameter, difference wrt one.
Null outgoing expansion
Definition at line 165 of file spheroid.h.
|
protected |
The 3-d projector on the 2-surface (contravariant-covariant form).
Definition at line 100 of file spheroid.h.
|
protected |
Definition at line 115 of file spheroid.h.
|
protected |
The 3-d covariant degenerated 2-metric on the surface.
Definition at line 104 of file spheroid.h.
|
protected |
Induced metric on the 2-surface 
Definition at line 117 of file spheroid.h.
|
protected |
The adapted normal vector field to spheroid in the 3-slice.
Definition at line 108 of file spheroid.h.
|
protected |
Trace K of the extrinsic curvature of the 3-slice.
Definition at line 124 of file spheroid.h.
|
protected |
Definition at line 147 of file spheroid.h.