LORENE
Lorene::Compobj_QI Class Reference

Base class for axisymmetric stationary compact objects in Quasi-Isotropic coordinates (under development). More...

#include <compobj.h>

Inheritance diagram for Lorene::Compobj_QI:
Lorene::Compobj Lorene::AltBH_QI Lorene::Kerr_QI Lorene::Star_QI Lorene::Boson_star

Public Member Functions

 Compobj_QI (Map &map_i)
 Standard constructor. More...
 
 Compobj_QI (const Compobj_QI &)
 Copy constructor. More...
 
 Compobj_QI (Map &map_i, FILE *)
 Constructor from a file (see sauve(FILE* )). More...
 
virtual ~Compobj_QI ()
 Destructor. More...
 
void operator= (const Compobj_QI &)
 Assignment to another Compobj_QI. More...
 
const Scalarget_bbb () const
 Returns the metric factor B. More...
 
const Scalarget_a_car () const
 Returns the square of the metric factor A. More...
 
const Scalarget_b_car () const
 Returns the square of the metric factor B. More...
 
const Scalarget_nphi () const
 Returns the metric coefficient $N^\varphi$. More...
 
const Scalarget_ak_car () const
 Returns the scalar $A^2 K_{ij} K^{ij}$. More...
 
virtual void sauve (FILE *) const
 Save in a file. More...
 
void gyoto_data (const char *file_name) const
 Save in a file for GYOTO. More...
 
virtual double angu_mom () const
 Angular momentum. More...
 
virtual double r_isco (int lmin, ostream *ost=0x0) const
 Coordinate r of the innermost stable circular orbit (ISCO). More...
 
virtual double f_isco (int lmin) const
 Orbital frequency at the innermost stable circular orbit (ISCO). More...
 
virtual double espec_isco (int lmin) const
 Energy of a particle at the ISCO. More...
 
virtual double lspec_isco (int lmin) const
 Angular momentum of a particle at the ISCO. More...
 
virtual double r_mb (int lmin, ostream *ost=0x0) const
 Coordinate r of the marginally bound circular orbit (R_mb). More...
 
virtual void update_metric ()
 Updates the 3-metric $\gamma_{ij}$ from A and B and the shift vector $\beta^i$ from $N^\phi$. More...
 
virtual void extrinsic_curvature ()
 Computes the extrinsic curvature and ak_car from nphi , nn and b_car . More...
 
Mapset_mp ()
 Read/write of the mapping. More...
 
const Mapget_mp () const
 Returns the mapping. More...
 
const Scalarget_nn () const
 Returns the lapse function N . More...
 
const Vectorget_beta () const
 Returns the shift vector $\beta^i$. More...
 
const Metricget_gamma () const
 Returns the 3-metric $\gamma_{ij}$. More...
 
const Scalarget_ener_euler () const
 Returns the total energy density E in the Eulerian frame. More...
 
const Vectorget_mom_euler () const
 Returns the total 3-momentum density $P^i$ in the Eulerian frame. More...
 
const Sym_tensorget_stress_euler () const
 Returns the stress tensor $S_{ij}$ with respect to the Eulerian observer. More...
 
const Sym_tensorget_kk () const
 Returns the extrinsic curvature tensor $K_{ij}$. More...
 
virtual double adm_mass () const
 ADM mass (computed as a surface integral at spatial infinity) More...
 

Protected Member Functions

virtual void del_deriv () const
 Deletes all the derived quantities. More...
 
void set_der_0x0 () const
 Sets to 0x0 all the pointers on derived quantities. More...
 
virtual ostream & operator>> (ostream &) const
 Operator >> (virtual function called by the operator <<). More...
 

Protected Attributes

Scalar a_car
 Square of the metric factor A. More...
 
Scalar bbb
 Metric factor B. More...
 
Scalar b_car
 Square of the metric factor B. More...
 
Scalar nphi
 Metric coefficient $N^\varphi$. More...
 
Scalar ak_car
 Scalar $A^2 K_{ij} K^{ij}$. More...
 
double * p_angu_mom
 Angular momentum. More...
 
double * p_r_isco
 Coordinate r of the ISCO. More...
 
double * p_f_isco
 Orbital frequency of the ISCO. More...
 
double * p_espec_isco
 Specific energy of a particle at the ISCO. More...
 
double * p_lspec_isco
 Specific angular momentum of a particle at the ISCO. More...
 
double * p_r_mb
 Coordinate r of the marginally bound orbit. More...
 
Mapmp
 Mapping describing the coordinate system (r,theta,phi) More...
 
Scalar nn
 Lapse function N . More...
 
Vector beta
 Shift vector $\beta^i$. More...
 
Metric gamma
 3-metric $\gamma_{ij}$ More...
 
Scalar ener_euler
 Total energy density E in the Eulerian frame. More...
 
Vector mom_euler
 Total 3-momentum density $P^i$ in the Eulerian frame. More...
 
Sym_tensor stress_euler
 Stress tensor $S_{ij}$ with respect to the Eulerian observer. More...
 
Sym_tensor kk
 Extrinsic curvature tensor $K_{ij}$. More...
 
double * p_adm_mass
 ADM mass. More...
 

Detailed Description

Base class for axisymmetric stationary compact objects in Quasi-Isotropic coordinates (under development).

()

The metric is expressed in Quasi-Isotropic (QI) coordinates :

\[ ds^2 = - N^2 dt^2 + A^2 (dr^2 + r^2 d\theta^2) + B^2 r^2 \sin^2\theta (d\varphi - N^\varphi dt)^2 \]

Definition at line 283 of file compobj.h.

Constructor & Destructor Documentation

◆ Compobj_QI() [1/3]

Lorene::Compobj_QI::Compobj_QI ( Map map_i)

Standard constructor.

Parameters
mp_iMapping on which the object is defined

Definition at line 89 of file compobj_QI.C.

References a_car, ak_car, b_car, bbb, nphi, set_der_0x0(), and Lorene::Scalar::std_spectral_base().

◆ Compobj_QI() [2/3]

Lorene::Compobj_QI::Compobj_QI ( const Compobj_QI co)

Copy constructor.

Definition at line 113 of file compobj_QI.C.

References set_der_0x0().

◆ Compobj_QI() [3/3]

Lorene::Compobj_QI::Compobj_QI ( Map map_i,
FILE *  fich 
)

Constructor from a file (see sauve(FILE* )).

Parameters
mp_iMapping on which the object is defined
fichinput file (must have been created by the function sauve)

Definition at line 128 of file compobj_QI.C.

References b_car, bbb, extrinsic_curvature(), Lorene::Map::get_mg(), Lorene::Compobj::mp, Lorene::Compobj::nn, set_der_0x0(), and update_metric().

◆ ~Compobj_QI()

Lorene::Compobj_QI::~Compobj_QI ( )
virtual

Destructor.

Definition at line 155 of file compobj_QI.C.

References del_deriv().

Member Function Documentation

◆ adm_mass()

double Lorene::Compobj::adm_mass ( ) const
virtualinherited

◆ angu_mom()

double Lorene::Compobj_QI::angu_mom ( ) const
virtual

Angular momentum.

Reimplemented in Lorene::Star_QI.

Definition at line 93 of file compobj_QI_global.C.

References p_angu_mom.

◆ del_deriv()

void Lorene::Compobj_QI::del_deriv ( ) const
protectedvirtual

Deletes all the derived quantities.

Reimplemented from Lorene::Compobj.

Reimplemented in Lorene::AltBH_QI, Lorene::Kerr_QI, Lorene::Star_QI, and Lorene::Boson_star.

Definition at line 166 of file compobj_QI.C.

References Lorene::Compobj::del_deriv(), p_angu_mom, p_espec_isco, p_f_isco, p_lspec_isco, p_r_isco, p_r_mb, and set_der_0x0().

◆ espec_isco()

double Lorene::Compobj_QI::espec_isco ( int  lmin) const
virtual

Energy of a particle at the ISCO.

Definition at line 323 of file compobj_QI_global.C.

References p_espec_isco, and r_isco().

◆ extrinsic_curvature()

void Lorene::Compobj_QI::extrinsic_curvature ( )
virtual

◆ f_isco()

double Lorene::Compobj_QI::f_isco ( int  lmin) const
virtual

Orbital frequency at the innermost stable circular orbit (ISCO).

Definition at line 289 of file compobj_QI_global.C.

References p_f_isco, and r_isco().

◆ get_a_car()

const Scalar& Lorene::Compobj_QI::get_a_car ( ) const
inline

Returns the square of the metric factor A.

Definition at line 381 of file compobj.h.

References a_car.

◆ get_ak_car()

const Scalar& Lorene::Compobj_QI::get_ak_car ( ) const
inline

Returns the scalar $A^2 K_{ij} K^{ij}$.

For axisymmetric stars, this quantity is related to the derivatives of $N^\varphi$ by

\[ A^2 K_{ij} K^{ij} = {B^2 \over 2 N^2} \, r^2\sin^2\theta \, \left[ \left( {\partial N^\varphi \over \partial r} \right) ^2 + {1\over r^2} \left( {\partial N^\varphi \over \partial \theta} \right) ^2 \right] \ . \]

In particular it is related to the quantities $k_1$ and $k_2$ introduced by Eqs. (3.7) and (3.8) of Bonazzola et al. Astron. Astrophys. 278 , 421 (1993) by

\[ A^2 K_{ij} K^{ij} = 2 A^2 (k_1^2 + k_2^2) \ . \]

Definition at line 407 of file compobj.h.

References ak_car.

◆ get_b_car()

const Scalar& Lorene::Compobj_QI::get_b_car ( ) const
inline

Returns the square of the metric factor B.

Definition at line 384 of file compobj.h.

References b_car.

◆ get_bbb()

const Scalar& Lorene::Compobj_QI::get_bbb ( ) const
inline

Returns the metric factor B.

Definition at line 378 of file compobj.h.

References bbb.

◆ get_beta()

const Vector& Lorene::Compobj::get_beta ( ) const
inlineinherited

Returns the shift vector $\beta^i$.

Definition at line 216 of file compobj.h.

References Lorene::Compobj::beta.

◆ get_ener_euler()

const Scalar& Lorene::Compobj::get_ener_euler ( ) const
inlineinherited

Returns the total energy density E in the Eulerian frame.

Definition at line 222 of file compobj.h.

References Lorene::Compobj::ener_euler.

◆ get_gamma()

const Metric& Lorene::Compobj::get_gamma ( ) const
inlineinherited

Returns the 3-metric $\gamma_{ij}$.

Definition at line 219 of file compobj.h.

References Lorene::Compobj::gamma.

◆ get_kk()

const Sym_tensor& Lorene::Compobj::get_kk ( ) const
inlineinherited

Returns the extrinsic curvature tensor $K_{ij}$.

Definition at line 231 of file compobj.h.

References Lorene::Compobj::kk.

◆ get_mom_euler()

const Vector& Lorene::Compobj::get_mom_euler ( ) const
inlineinherited

Returns the total 3-momentum density $P^i$ in the Eulerian frame.

Definition at line 225 of file compobj.h.

References Lorene::Compobj::mom_euler.

◆ get_mp()

const Map& Lorene::Compobj::get_mp ( ) const
inlineinherited

Returns the mapping.

Definition at line 210 of file compobj.h.

References Lorene::Compobj::mp.

◆ get_nn()

const Scalar& Lorene::Compobj::get_nn ( ) const
inlineinherited

Returns the lapse function N .

Definition at line 213 of file compobj.h.

References Lorene::Compobj::nn.

◆ get_nphi()

const Scalar& Lorene::Compobj_QI::get_nphi ( ) const
inline

Returns the metric coefficient $N^\varphi$.

Definition at line 387 of file compobj.h.

References nphi.

◆ get_stress_euler()

const Sym_tensor& Lorene::Compobj::get_stress_euler ( ) const
inlineinherited

Returns the stress tensor $S_{ij}$ with respect to the Eulerian observer.

Definition at line 228 of file compobj.h.

References Lorene::Compobj::stress_euler.

◆ gyoto_data()

◆ lspec_isco()

double Lorene::Compobj_QI::lspec_isco ( int  lmin) const
virtual

Angular momentum of a particle at the ISCO.

Definition at line 306 of file compobj_QI_global.C.

References p_lspec_isco, and r_isco().

◆ operator=()

void Lorene::Compobj_QI::operator= ( const Compobj_QI co)

Assignment to another Compobj_QI.

Definition at line 198 of file compobj_QI.C.

References a_car, ak_car, b_car, bbb, del_deriv(), nphi, and Lorene::Compobj::operator=().

◆ operator>>()

ostream & Lorene::Compobj_QI::operator>> ( ostream &  ost) const
protectedvirtual

Operator >> (virtual function called by the operator <<).

Reimplemented from Lorene::Compobj.

Reimplemented in Lorene::AltBH_QI, Lorene::Kerr_QI, Lorene::Star_QI, and Lorene::Boson_star.

Definition at line 267 of file compobj_QI.C.

References a_car, ak_car, b_car, nphi, Lorene::Compobj::operator>>(), r_isco(), r_mb(), and Lorene::Scalar::val_grid_point().

◆ r_isco()

double Lorene::Compobj_QI::r_isco ( int  lmin,
ostream *  ost = 0x0 
) const
virtual

Coordinate r of the innermost stable circular orbit (ISCO).

Parameters
lminindex of the innermost domain in which the ISCO is searched: the ISCO is searched inwards from the last but one domain to the domain of index lmin.
ostoutput stream to give details of the computation; if set to 0x0 [default value], no details will be given.

Definition at line 111 of file compobj_QI_global.C.

References Lorene::Param::add_int(), Lorene::Param::add_scalar(), Lorene::Tensor::annule_domain(), bbb, Lorene::Scalar::dsdr(), Lorene::Map::get_mg(), Lorene::Mg3d::get_nzone(), Lorene::Scalar::get_spectral_va(), Lorene::Compobj::mp, Lorene::Compobj::nn, nphi, p_espec_isco, p_f_isco, p_lspec_isco, p_r_isco, Lorene::Map::r, Lorene::sqrt(), Lorene::Scalar::std_spectral_base(), Lorene::Valeur::val_point(), Lorene::Map::val_r(), and Lorene::zerosec().

◆ r_mb()

◆ sauve()

void Lorene::Compobj_QI::sauve ( FILE *  fich) const
virtual

Save in a file.

Reimplemented from Lorene::Compobj.

Reimplemented in Lorene::AltBH_QI, Lorene::Kerr_QI, Lorene::Star_QI, and Lorene::Boson_star.

Definition at line 219 of file compobj_QI.C.

References a_car, bbb, Lorene::Compobj::nn, nphi, and Lorene::Scalar::sauve().

◆ set_der_0x0()

void Lorene::Compobj_QI::set_der_0x0 ( ) const
protected

Sets to 0x0 all the pointers on derived quantities.

Definition at line 181 of file compobj_QI.C.

References p_angu_mom, p_espec_isco, p_f_isco, p_lspec_isco, p_r_isco, and p_r_mb.

◆ set_mp()

Map& Lorene::Compobj::set_mp ( )
inlineinherited

Read/write of the mapping.

Definition at line 203 of file compobj.h.

References Lorene::Compobj::mp.

◆ update_metric()

void Lorene::Compobj_QI::update_metric ( )
virtual

Member Data Documentation

◆ a_car

Scalar Lorene::Compobj_QI::a_car
protected

Square of the metric factor A.

Definition at line 290 of file compobj.h.

◆ ak_car

Scalar Lorene::Compobj_QI::ak_car
protected

Scalar $A^2 K_{ij} K^{ij}$.

For axisymmetric stars, this quantity is related to the derivatives of $N^\varphi$ by

\[ A^2 K_{ij} K^{ij} = {B^2 \over 2 N^2} \, r^2\sin^2\theta \, \left[ \left( {\partial N^\varphi \over \partial r} \right) ^2 + {1\over r^2} \left( {\partial N^\varphi \over \partial \theta} \right) ^2 \right] \ . \]

In particular it is related to the quantities $k_1$ and $k_2$ introduced by Eqs.~(3.7) and (3.8) of Bonazzola et al. Astron. Astrophys. 278 , 421 (1993) by

\[ A^2 K_{ij} K^{ij} = 2 A^2 (k_1^2 + k_2^2) \ . \]

Definition at line 318 of file compobj.h.

◆ b_car

Scalar Lorene::Compobj_QI::b_car
protected

Square of the metric factor B.

Definition at line 296 of file compobj.h.

◆ bbb

Scalar Lorene::Compobj_QI::bbb
protected

Metric factor B.

Definition at line 293 of file compobj.h.

◆ beta

Vector Lorene::Compobj::beta
protectedinherited

Shift vector $\beta^i$.

Definition at line 141 of file compobj.h.

◆ ener_euler

Scalar Lorene::Compobj::ener_euler
protectedinherited

Total energy density E in the Eulerian frame.

Definition at line 147 of file compobj.h.

◆ gamma

Metric Lorene::Compobj::gamma
protectedinherited

3-metric $\gamma_{ij}$

Definition at line 144 of file compobj.h.

◆ kk

Sym_tensor Lorene::Compobj::kk
protectedinherited

Extrinsic curvature tensor $K_{ij}$.

Definition at line 156 of file compobj.h.

◆ mom_euler

Vector Lorene::Compobj::mom_euler
protectedinherited

Total 3-momentum density $P^i$ in the Eulerian frame.

Definition at line 150 of file compobj.h.

◆ mp

Map& Lorene::Compobj::mp
protectedinherited

Mapping describing the coordinate system (r,theta,phi)

Definition at line 135 of file compobj.h.

◆ nn

Scalar Lorene::Compobj::nn
protectedinherited

Lapse function N .

Definition at line 138 of file compobj.h.

◆ nphi

Scalar Lorene::Compobj_QI::nphi
protected

Metric coefficient $N^\varphi$.

Definition at line 299 of file compobj.h.

◆ p_adm_mass

double* Lorene::Compobj::p_adm_mass
mutableprotectedinherited

ADM mass.

Definition at line 161 of file compobj.h.

◆ p_angu_mom

double* Lorene::Compobj_QI::p_angu_mom
mutableprotected

Angular momentum.

Definition at line 324 of file compobj.h.

◆ p_espec_isco

double* Lorene::Compobj_QI::p_espec_isco
mutableprotected

Specific energy of a particle at the ISCO.

Definition at line 328 of file compobj.h.

◆ p_f_isco

double* Lorene::Compobj_QI::p_f_isco
mutableprotected

Orbital frequency of the ISCO.

Definition at line 326 of file compobj.h.

◆ p_lspec_isco

double* Lorene::Compobj_QI::p_lspec_isco
mutableprotected

Specific angular momentum of a particle at the ISCO.

Definition at line 330 of file compobj.h.

◆ p_r_isco

double* Lorene::Compobj_QI::p_r_isco
mutableprotected

Coordinate r of the ISCO.

Definition at line 325 of file compobj.h.

◆ p_r_mb

double* Lorene::Compobj_QI::p_r_mb
mutableprotected

Coordinate r of the marginally bound orbit.

Definition at line 331 of file compobj.h.

◆ stress_euler

Sym_tensor Lorene::Compobj::stress_euler
protectedinherited

Stress tensor $S_{ij}$ with respect to the Eulerian observer.

Definition at line 153 of file compobj.h.


The documentation for this class was generated from the following files: