Lorene::Star_QI Class Reference
[Stationary compact objects (under development)]

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

#include <compobj.h>

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

List of all members.

Public Member Functions

 Star_QI (Map &mp_i)
 Standard constructor.
 Star_QI (const Star_QI &)
 Copy constructor.
 Star_QI (Map &mp_i, FILE *fich)
 Constructor from a file (see sauve(FILE*) ).
virtual ~Star_QI ()
 Destructor.
void operator= (const Star_QI &)
 Assignment to another Star_QI.
const Scalarget_logn () const
 Returns the logarithm of the lapse N.
const Scalarget_tnphi () const
 Returns the component $\tilde N^\varphi = N^\varphi r\sin\theta$ of the shift vector.
const Scalarget_nuf () const
 Returns the part of the Metric potential $\nu = \ln N$ = logn generated by the matter terms.
const Scalarget_nuq () const
 Returns the Part of the Metric potential $\nu = \ln N$ = logn generated by the quadratic terms.
const Scalarget_dzeta () const
 Returns the Metric potential $\zeta = \ln(AN)$.
const Scalarget_tggg () const
 Returns the Metric potential $\tilde G = (NB-1) r\sin\theta$.
const Vectorget_w_shift () const
 Returns the vector $W^i$ used in the decomposition of shift , following Shibata's prescription [Prog.
const Scalarget_khi_shift () const
 Returns the scalar $\chi$ used in the decomposition of shift following Shibata's prescription [Prog.
virtual void sauve (FILE *) const
 Save in a file.
virtual double mass_g () const
 Gravitational mass.
virtual double angu_mom () const
 Angular momentum.
virtual double grv2 () const
 Error on the virial identity GRV2.
virtual double grv3 (ostream *ost=0x0) const
 Error on the virial identity GRV3.
virtual double mom_quad () const
 Quadrupole moment.
void update_metric ()
 Computes metric coefficients from known potentials.
void fait_shift ()
 Computes shift from w_shift and khi_shift according to Shibata's prescription [Prog.
void fait_nphi ()
 Computes tnphi and nphi from the Cartesian components of the shift, stored in shift .
const Scalarget_bbb () const
 Returns the metric factor B.
const Scalarget_a_car () const
 Returns the square of the metric factor A.
const Scalarget_b_car () const
 Returns the square of the metric factor B.
const Scalarget_nphi () const
 Returns the metric coefficient $N^\varphi$.
const Scalarget_ak_car () const
 Returns the scalar $A^2 K_{ij} K^{ij}$.
void gyoto_data (const char *file_name) const
 Save in a file for GYOTO.
virtual double r_isco (int lmin, ostream *ost=0x0) const
 Coordinate r of the innermost stable circular orbit (ISCO).
virtual double f_isco (int lmin) const
 Orbital frequency at the innermost stable circular orbit (ISCO).
virtual double espec_isco (int lmin) const
 Energy of a particle at the ISCO.
virtual double lspec_isco (int lmin) const
 Angular momentum of a particle at the ISCO.
virtual double r_mb (int lmin, ostream *ost=0x0) const
 Coordinate r of the marginally bound circular orbit (R_mb).
virtual void extrinsic_curvature ()
 Computes the extrinsic curvature and ak_car from nphi , nn and b_car .
Mapset_mp ()
 Read/write of the mapping.
const Mapget_mp () const
 Returns the mapping.
const Scalarget_nn () const
 Returns the lapse function N .
const Vectorget_beta () const
 Returns the shift vector $\beta^i$.
const Metricget_gamma () const
 Returns the 3-metric $\gamma_{ij}$.
const Scalarget_ener_euler () const
 Returns the total energy density E in the Eulerian frame.
const Vectorget_mom_euler () const
 Returns the total 3-momentum density $P^i$ in the Eulerian frame.
const Sym_tensorget_stress_euler () const
 Returns the stress tensor $S_{ij}$ with respect to the Eulerian observer.
const Sym_tensorget_kk () const
 Returns the extrinsic curvature tensor $K_{ij}$.
virtual double adm_mass () const
 ADM mass (computed as a surface integral at spatial infinity).

Static Public Member Functions

static double lambda_grv2 (const Scalar &sou_m, const Scalar &sou_q)
 Computes the coefficient $\lambda$ which ensures that the GRV2 virial identity is satisfied.

Protected Member Functions

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

Protected Attributes

Scalar logn
 Logarithm of the lapse N .
Scalar tnphi
 Component $\tilde N^\varphi = N^\varphi r\sin\theta$ of the shift vector.
Scalar nuf
 Part of the Metric potential $\nu = \ln N$ = logn generated by the matter terms.
Scalar nuq
 Part of the Metric potential $\nu = \ln N$ = logn generated by the quadratic terms.
Scalar dzeta
 Metric potential $\zeta = \ln(AN)$.
Scalar tggg
 Metric potential $\tilde G = (NB-1) r\sin\theta$.
Vector w_shift
 Vector $W^i$ used in the decomposition of shift , following Shibata's prescription [Prog.
Scalar khi_shift
 Scalar $\chi$ used in the decomposition of shift , following Shibata's prescription [Prog.
Scalar ssjm1_nuf
 Effective source at the previous step for the resolution of the Poisson equation for nuf by means of Map_et::poisson .
Scalar ssjm1_nuq
 Effective source at the previous step for the resolution of the Poisson equation for nuq by means of Map_et::poisson .
Scalar ssjm1_dzeta
 Effective source at the previous step for the resolution of the Poisson equation for dzeta .
Scalar ssjm1_tggg
 Effective source at the previous step for the resolution of the Poisson equation for tggg .
Scalar ssjm1_khi
 Effective source at the previous step for the resolution of the Poisson equation for the scalar $\chi$ by means of Map_et::poisson .
Vector ssjm1_wshift
 Effective source at the previous step for the resolution of the vector Poisson equation for $W^i$.
double * p_grv2
 Error on the virial identity GRV2.
double * p_grv3
 Error on the virial identity GRV3.
double * p_mom_quad
 Quadrupole moment.
double * p_mass_g
 Gravitational mass (ADM mass as a volume integral).
Scalar a_car
 Square of the metric factor A.
Scalar bbb
 Metric factor B.
Scalar b_car
 Square of the metric factor B.
Scalar nphi
 Metric coefficient $N^\varphi$.
Scalar ak_car
 Scalar $A^2 K_{ij} K^{ij}$.
double * p_angu_mom
 Angular momentum.
double * p_r_isco
 Coordinate r of the ISCO.
double * p_f_isco
 Orbital frequency of the ISCO.
double * p_espec_isco
 Specific energy of a particle at the ISCO.
double * p_lspec_isco
 Specific angular momentum of a particle at the ISCO.
double * p_r_mb
 Coordinate r of the marginally bound orbit.
Mapmp
 Mapping describing the coordinate system (r,theta,phi).
Scalar nn
 Lapse function N .
Vector beta
 Shift vector $\beta^i$.
Metric gamma
 3-metric $\gamma_{ij}$
Scalar ener_euler
 Total energy density E in the Eulerian frame.
Vector mom_euler
 Total 3-momentum density $P^i$ in the Eulerian frame.
Sym_tensor stress_euler
 Stress tensor $S_{ij}$ with respect to the Eulerian observer.
Sym_tensor kk
 Extrinsic curvature tensor $K_{ij}$.
double * p_adm_mass
 ADM mass.

Friends

ostream & operator<< (ostream &, const Compobj &)
 Display.

Detailed Description

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

()

The time slice $t=\mathrm{const}$ has the topology of $R^3$ and 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 490 of file compobj.h.


Constructor & Destructor Documentation

Lorene::Star_QI::Star_QI ( Map mp_i  ) 

Standard constructor.

Parameters:
mp_i Mapping on which the star is contructed

Definition at line 71 of file star_QI.C.

References Lorene::Compobj::beta, dzeta, khi_shift, logn, nuf, nuq, set_der_0x0(), Lorene::Tensor::set_etat_zero(), ssjm1_dzeta, ssjm1_khi, ssjm1_nuf, ssjm1_nuq, ssjm1_tggg, ssjm1_wshift, tggg, tnphi, and w_shift.

Lorene::Star_QI::Star_QI ( const Star_QI st  ) 

Copy constructor.

Definition at line 116 of file star_QI.C.

References set_der_0x0().

Lorene::Star_QI::Star_QI ( Map mp_i,
FILE *  fich 
)

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

Parameters:
mp_i Mapping on which the star is constructed
fich input file (must have been created by the function Star_QI::sauve )

Definition at line 140 of file star_QI.C.

References dzeta, fait_nphi(), fait_shift(), Lorene::Map::get_bvect_cart(), Lorene::Map::get_mg(), khi_shift, logn, Lorene::Compobj::mp, nuf, nuq, set_der_0x0(), ssjm1_dzeta, ssjm1_khi, ssjm1_nuf, ssjm1_nuq, ssjm1_tggg, ssjm1_wshift, tggg, update_metric(), and w_shift.

Lorene::Star_QI::~Star_QI (  )  [virtual]

Destructor.

Definition at line 214 of file star_QI.C.

References del_deriv().


Member Function Documentation

double Lorene::Compobj::adm_mass (  )  const [virtual, inherited]
double Lorene::Star_QI::angu_mom (  )  const [virtual]
void Lorene::Star_QI::del_deriv (  )  const [protected, virtual]

Deletes all the derived quantities.

Reimplemented from Lorene::Compobj_QI.

Reimplemented in Lorene::Boson_star.

Definition at line 225 of file star_QI.C.

References p_grv2, p_grv3, p_mass_g, p_mom_quad, and set_der_0x0().

double Lorene::Compobj_QI::espec_isco ( int  lmin  )  const [virtual, inherited]

Energy of a particle at the ISCO.

Definition at line 323 of file compobj_QI_global.C.

References Lorene::Compobj_QI::p_espec_isco, and Lorene::Compobj_QI::r_isco().

void Lorene::Compobj_QI::extrinsic_curvature (  )  [virtual, inherited]
double Lorene::Compobj_QI::f_isco ( int  lmin  )  const [virtual, inherited]

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

Definition at line 289 of file compobj_QI_global.C.

References Lorene::Compobj_QI::p_f_isco, and Lorene::Compobj_QI::r_isco().

void Lorene::Star_QI::fait_nphi (  ) 

Computes tnphi and nphi from the Cartesian components of the shift, stored in shift .

Definition at line 389 of file star_QI.C.

References Lorene::Compobj::beta, Lorene::Map::comp_p_from_cartesian(), Lorene::Scalar::div_rsint(), Lorene::Compobj::mp, Lorene::Compobj_QI::nphi, and tnphi.

void Lorene::Star_QI::fait_shift (  ) 
const Scalar& Lorene::Compobj_QI::get_a_car (  )  const [inline, inherited]

Returns the square of the metric factor A.

Definition at line 381 of file compobj.h.

References Lorene::Compobj_QI::a_car.

const Scalar& Lorene::Compobj_QI::get_ak_car (  )  const [inline, inherited]

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 Lorene::Compobj_QI::ak_car.

const Scalar& Lorene::Compobj_QI::get_b_car (  )  const [inline, inherited]

Returns the square of the metric factor B.

Definition at line 384 of file compobj.h.

References Lorene::Compobj_QI::b_car.

const Scalar& Lorene::Compobj_QI::get_bbb (  )  const [inline, inherited]

Returns the metric factor B.

Definition at line 378 of file compobj.h.

References Lorene::Compobj_QI::bbb.

const Vector& Lorene::Compobj::get_beta (  )  const [inline, inherited]

Returns the shift vector $\beta^i$.

Definition at line 216 of file compobj.h.

References Lorene::Compobj::beta.

const Scalar& Lorene::Star_QI::get_dzeta (  )  const [inline]

Returns the Metric potential $\zeta = \ln(AN)$.

Definition at line 658 of file compobj.h.

References dzeta.

const Scalar& Lorene::Compobj::get_ener_euler (  )  const [inline, inherited]

Returns the total energy density E in the Eulerian frame.

Definition at line 222 of file compobj.h.

References Lorene::Compobj::ener_euler.

const Metric& Lorene::Compobj::get_gamma (  )  const [inline, inherited]

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

Definition at line 219 of file compobj.h.

References Lorene::Compobj::gamma.

const Scalar& Lorene::Star_QI::get_khi_shift (  )  const [inline]

Returns the scalar $\chi$ used in the decomposition of shift following Shibata's prescription [Prog.

Theor. Phys. 101 , 1199 (1999)] :

\[ N^i = {7\over 8} W^i - {1\over 8} \left(\nabla^i\chi+\nabla^iW^kx_k\right) \]

NB: w_shift contains the components of $W^i$ with respect to the Cartesian triad associated with the mapping mp .

Definition at line 689 of file compobj.h.

References khi_shift.

const Sym_tensor& Lorene::Compobj::get_kk (  )  const [inline, inherited]

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

Definition at line 231 of file compobj.h.

References Lorene::Compobj::kk.

const Scalar& Lorene::Star_QI::get_logn (  )  const [inline]

Returns the logarithm of the lapse N.

Definition at line 639 of file compobj.h.

References logn.

const Vector& Lorene::Compobj::get_mom_euler (  )  const [inline, inherited]

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.

const Map& Lorene::Compobj::get_mp (  )  const [inline, inherited]

Returns the mapping.

Definition at line 210 of file compobj.h.

References Lorene::Compobj::mp.

const Scalar& Lorene::Compobj::get_nn (  )  const [inline, inherited]

Returns the lapse function N .

Definition at line 213 of file compobj.h.

References Lorene::Compobj::nn.

const Scalar& Lorene::Compobj_QI::get_nphi (  )  const [inline, inherited]

Returns the metric coefficient $N^\varphi$.

Definition at line 387 of file compobj.h.

References Lorene::Compobj_QI::nphi.

const Scalar& Lorene::Star_QI::get_nuf (  )  const [inline]

Returns the part of the Metric potential $\nu = \ln N$ = logn generated by the matter terms.

Definition at line 650 of file compobj.h.

References nuf.

const Scalar& Lorene::Star_QI::get_nuq (  )  const [inline]

Returns the Part of the Metric potential $\nu = \ln N$ = logn generated by the quadratic terms.

Definition at line 655 of file compobj.h.

References nuq.

const Sym_tensor& Lorene::Compobj::get_stress_euler (  )  const [inline, inherited]

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.

const Scalar& Lorene::Star_QI::get_tggg (  )  const [inline]

Returns the Metric potential $\tilde G = (NB-1) r\sin\theta$.

Definition at line 661 of file compobj.h.

References tggg.

const Scalar& Lorene::Star_QI::get_tnphi (  )  const [inline]

Returns the component $\tilde N^\varphi = N^\varphi r\sin\theta$ of the shift vector.

Definition at line 645 of file compobj.h.

References tnphi.

const Vector& Lorene::Star_QI::get_w_shift (  )  const [inline]

Returns the vector $W^i$ used in the decomposition of shift , following Shibata's prescription [Prog.

Theor. Phys. 101 , 1199 (1999)] :

\[ N^i = {7\over 8} W^i - {1\over 8} \left(\nabla^i\chi+\nabla^iW^kx_k\right) \]

NB: w_shift contains the components of $W^i$ with respect to the Cartesian triad associated with the mapping mp .

Definition at line 675 of file compobj.h.

References w_shift.

double Lorene::Star_QI::grv2 (  )  const [virtual]
double Lorene::Star_QI::grv3 ( ostream *  ost = 0x0  )  const [virtual]

Error on the virial identity GRV3.

The error is computed as the integral defined by Eq. (43) of [Gourgoulhon and Bonazzola, Class. Quantum Grav. 11, 443 (1994)] divided by the integral of the matter terms.

Parameters:
ost output stream to give details of the computation; if set to 0x0 [default value], no details will be given.

Definition at line 151 of file star_QI_global.C.

References Lorene::Compobj_QI::a_car, Lorene::Compobj_QI::ak_car, Lorene::Compobj_QI::bbb, Lorene::Scalar::derive_cov(), Lorene::Scalar::dsdr(), dzeta, Lorene::Map::flat_met_spher(), Lorene::Compobj::gamma, Lorene::Scalar::get_dzpuis(), Lorene::Scalar::get_etat(), Lorene::Scalar::integrale(), Lorene::log(), logn, Lorene::Compobj::mp, Lorene::Valeur::mult_ct(), p_grv3, Lorene::Scalar::set_dzpuis(), Lorene::Scalar::set_spectral_va(), Lorene::Scalar::srdsdt(), Lorene::Valeur::ssint(), Lorene::Scalar::std_spectral_base(), Lorene::Compobj::stress_euler, Lorene::Valeur::sx(), Lorene::Tensor::trace(), and Lorene::Map_radial::xsr.

void Lorene::Compobj_QI::gyoto_data ( const char *  file_name  )  const [inherited]
double Lorene::Star_QI::lambda_grv2 ( const Scalar sou_m,
const Scalar sou_q 
) [static]

Computes the coefficient $\lambda$ which ensures that the GRV2 virial identity is satisfied.

$\lambda$ is the coefficient by which one must multiply the quadratic source term $\sigma_q$ of the 2-D Poisson equation

\[ \Delta_2 u = \sigma_m + \sigma_q \]

in order that the total source does not contain any monopolar term, i.e. in order that

\[ \int_0^{2\pi} \int_0^{+\infty} \sigma(r, \theta) \, r \, dr \, d\theta = 0 \ , \]

where $\sigma = \sigma_m + \sigma_q$. $\lambda$ is computed according to the formula

\[ \lambda = - { \int_0^{2\pi} \int_0^{+\infty} \sigma_m(r, \theta) \, r \, dr \, d\theta \over \int_0^{2\pi} \int_0^{+\infty} \sigma_q(r, \theta) \, r \, dr \, d\theta } \ . \]

Then, by construction, the new source $\sigma' = \sigma_m + \lambda \sigma_q$ has a vanishing monopolar term.

Parameters:
sou_m [input] matter source term $\sigma_m$
sou_q [input] quadratic source term $\sigma_q$
Returns:
value of $\lambda$

Definition at line 324 of file star_QI_global.C.

References Lorene::Valeur::c, Lorene::Scalar::check_dzpuis(), Lorene::Valeur::coef_i(), Lorene::Map_radial::dxdr, Lorene::Map_af::get_alpha(), Lorene::Map_af::get_beta(), Lorene::Valeur::get_etat(), Lorene::Scalar::get_etat(), Lorene::Tbl::get_etat(), Lorene::Mg3d::get_grille3d(), Lorene::Map::get_mg(), Lorene::Tensor::get_mp(), Lorene::Mg3d::get_np(), Lorene::Mg3d::get_nr(), Lorene::Mg3d::get_nt(), Lorene::Mg3d::get_nzone(), Lorene::Scalar::get_spectral_va(), Lorene::Mg3d::get_type_r(), Lorene::Map_af::set_alpha(), Lorene::Map_af::set_beta(), Lorene::Tbl::t, Lorene::Mtbl::t, Lorene::Map::val_r(), Lorene::Grille3d::x, and Lorene::Map_radial::xsr.

double Lorene::Compobj_QI::lspec_isco ( int  lmin  )  const [virtual, inherited]

Angular momentum of a particle at the ISCO.

Definition at line 306 of file compobj_QI_global.C.

References Lorene::Compobj_QI::p_lspec_isco, and Lorene::Compobj_QI::r_isco().

double Lorene::Star_QI::mass_g (  )  const [virtual]
double Lorene::Star_QI::mom_quad (  )  const [virtual]

Quadrupole moment.

The quadrupole moment Q is defined according to Eq. (7) of [Salgado, Bonazzola, Gourgoulhon and Haensel, Astron. Astrophys. 291 , 155 (1994)]. At the Newtonian limit it is related to the component ${\bar I}_{zz}$ of the MTW (1973) reduced quadrupole moment ${\bar I}_{ij}$ by: $Q = -3/2 {\bar I}_{zz}$. Note that Q is the negative of the quadrupole moment defined by Laarakkers and Poisson, Astrophys. J. 512 , 282 (1999).

Definition at line 253 of file star_QI_global.C.

References Lorene::Compobj_QI::a_car, Lorene::Compobj_QI::ak_car, Lorene::Compobj_QI::bbb, Lorene::Scalar::check_dzpuis(), Lorene::Scalar::derive_cov(), Lorene::Compobj::ener_euler, Lorene::Map::flat_met_spher(), Lorene::Compobj::gamma, Lorene::Scalar::get_etat(), Lorene::Scalar::inc_dzpuis(), Lorene::Scalar::integrale(), Lorene::log(), logn, Lorene::Compobj::mp, Lorene::Valeur::mult_ct(), Lorene::Scalar::mult_r(), p_mom_quad, Lorene::Scalar::set_spectral_va(), Lorene::Scalar::std_spectral_base(), Lorene::Compobj::stress_euler, and Lorene::Tensor::trace().

void Lorene::Star_QI::operator= ( const Star_QI st  ) 

Assignment to another Star_QI.

Reimplemented from Lorene::Compobj_QI.

Reimplemented in Lorene::Boson_star.

Definition at line 253 of file star_QI.C.

References del_deriv(), dzeta, khi_shift, logn, nuf, nuq, ssjm1_dzeta, ssjm1_khi, ssjm1_nuf, ssjm1_nuq, ssjm1_tggg, ssjm1_wshift, tggg, tnphi, and w_shift.

ostream & Lorene::Star_QI::operator>> ( ostream &  ost  )  const [protected, virtual]

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

Reimplemented from Lorene::Compobj_QI.

Reimplemented in Lorene::Boson_star.

Definition at line 303 of file star_QI.C.

References angu_mom(), dzeta, grv2(), grv3(), logn, mass_g(), mom_quad(), nuf, nuq, Lorene::pow(), and Lorene::Scalar::val_grid_point().

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

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

Parameters:
lmin index 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.
ost output 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(), Lorene::Compobj_QI::bbb, Lorene::Scalar::dsdr(), Lorene::Map::get_mg(), Lorene::Mg3d::get_nzone(), Lorene::Scalar::get_spectral_va(), Lorene::Compobj::mp, Lorene::Compobj::nn, Lorene::Compobj_QI::nphi, Lorene::Compobj_QI::p_espec_isco, Lorene::Compobj_QI::p_f_isco, Lorene::Compobj_QI::p_lspec_isco, Lorene::Compobj_QI::p_r_isco, Lorene::Map::r, Lorene::sqrt(), Lorene::Scalar::std_spectral_base(), Lorene::Valeur::val_point(), Lorene::Map::val_r(), and Lorene::zerosec().

double Lorene::Compobj_QI::r_mb ( int  lmin,
ostream *  ost = 0x0 
) const [virtual, inherited]
void Lorene::Star_QI::sauve ( FILE *  fich  )  const [virtual]

Save in a file.

Reimplemented from Lorene::Compobj_QI.

Reimplemented in Lorene::Boson_star.

Definition at line 282 of file star_QI.C.

References dzeta, khi_shift, nuf, nuq, Lorene::Tensor::sauve(), Lorene::Scalar::sauve(), ssjm1_dzeta, ssjm1_khi, ssjm1_nuf, ssjm1_nuq, ssjm1_tggg, ssjm1_wshift, tggg, and w_shift.

void Lorene::Star_QI::set_der_0x0 (  )  const [protected, virtual]

Sets to 0x0 all the pointers on derived quantities.

Reimplemented from Lorene::Compobj_QI.

Reimplemented in Lorene::Boson_star.

Definition at line 238 of file star_QI.C.

References p_grv2, p_grv3, p_mass_g, and p_mom_quad.

Map& Lorene::Compobj::set_mp (  )  [inline, inherited]

Read/write of the mapping.

Definition at line 203 of file compobj.h.

References Lorene::Compobj::mp.

void Lorene::Star_QI::update_metric (  )  [virtual]

Computes metric coefficients from known potentials.

The calculation is performed starting from the quantities logn , dzeta , tggg and shift , which are supposed to be up to date. From these, the following fields are updated: nnn , a_car , bbb and b_car, as well as the 3-metric gamma.

Reimplemented from Lorene::Compobj_QI.

Definition at line 413 of file star_QI.C.

References Lorene::Compobj_QI::a_car, Lorene::Compobj_QI::b_car, Lorene::Compobj_QI::bbb, del_deriv(), Lorene::Scalar::div_rsint(), dzeta, Lorene::exp(), logn, Lorene::Compobj::nn, Lorene::Scalar::std_spectral_base(), and tggg.


Friends And Related Function Documentation

ostream& operator<< ( ostream &  ,
const Compobj  
) [friend, inherited]

Display.


Member Data Documentation

Scalar Lorene::Compobj_QI::a_car [protected, inherited]

Square of the metric factor A.

Definition at line 290 of file compobj.h.

Scalar Lorene::Compobj_QI::ak_car [protected, inherited]

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.

Scalar Lorene::Compobj_QI::b_car [protected, inherited]

Square of the metric factor B.

Definition at line 296 of file compobj.h.

Scalar Lorene::Compobj_QI::bbb [protected, inherited]

Metric factor B.

Definition at line 293 of file compobj.h.

Vector Lorene::Compobj::beta [protected, inherited]

Shift vector $\beta^i$.

Definition at line 141 of file compobj.h.

Metric potential $\zeta = \ln(AN)$.

Definition at line 516 of file compobj.h.

Scalar Lorene::Compobj::ener_euler [protected, inherited]

Total energy density E in the Eulerian frame.

Definition at line 147 of file compobj.h.

Metric Lorene::Compobj::gamma [protected, inherited]

3-metric $\gamma_{ij}$

Definition at line 144 of file compobj.h.

Scalar $\chi$ used in the decomposition of shift , following Shibata's prescription [Prog.

Theor. Phys. 101 , 1199 (1999)] :

\[ N^i = {7\over 8} W^i - {1\over 8} \left(\nabla^i\chi+\nabla^iW^kx_k\right) \]

Definition at line 542 of file compobj.h.

Sym_tensor Lorene::Compobj::kk [protected, inherited]

Extrinsic curvature tensor $K_{ij}$.

Definition at line 156 of file compobj.h.

Logarithm of the lapse N .

Definition at line 498 of file compobj.h.

Vector Lorene::Compobj::mom_euler [protected, inherited]

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

Definition at line 150 of file compobj.h.

Map& Lorene::Compobj::mp [protected, inherited]

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

Definition at line 135 of file compobj.h.

Scalar Lorene::Compobj::nn [protected, inherited]

Lapse function N .

Definition at line 138 of file compobj.h.

Scalar Lorene::Compobj_QI::nphi [protected, inherited]

Metric coefficient $N^\varphi$.

Definition at line 299 of file compobj.h.

Part of the Metric potential $\nu = \ln N$ = logn generated by the matter terms.

Definition at line 508 of file compobj.h.

Part of the Metric potential $\nu = \ln N$ = logn generated by the quadratic terms.

Definition at line 513 of file compobj.h.

double* Lorene::Compobj::p_adm_mass [mutable, protected, inherited]

ADM mass.

Definition at line 161 of file compobj.h.

double* Lorene::Compobj_QI::p_angu_mom [mutable, protected, inherited]

Angular momentum.

Definition at line 324 of file compobj.h.

double* Lorene::Compobj_QI::p_espec_isco [mutable, protected, inherited]

Specific energy of a particle at the ISCO.

Definition at line 328 of file compobj.h.

double* Lorene::Compobj_QI::p_f_isco [mutable, protected, inherited]

Orbital frequency of the ISCO.

Definition at line 326 of file compobj.h.

double* Lorene::Star_QI::p_grv2 [mutable, protected]

Error on the virial identity GRV2.

Definition at line 588 of file compobj.h.

double* Lorene::Star_QI::p_grv3 [mutable, protected]

Error on the virial identity GRV3.

Definition at line 589 of file compobj.h.

double* Lorene::Compobj_QI::p_lspec_isco [mutable, protected, inherited]

Specific angular momentum of a particle at the ISCO.

Definition at line 330 of file compobj.h.

double* Lorene::Star_QI::p_mass_g [mutable, protected]

Gravitational mass (ADM mass as a volume integral).

Definition at line 591 of file compobj.h.

double* Lorene::Star_QI::p_mom_quad [mutable, protected]

Quadrupole moment.

Definition at line 590 of file compobj.h.

double* Lorene::Compobj_QI::p_r_isco [mutable, protected, inherited]

Coordinate r of the ISCO.

Definition at line 325 of file compobj.h.

double* Lorene::Compobj_QI::p_r_mb [mutable, protected, inherited]

Coordinate r of the marginally bound orbit.

Definition at line 331 of file compobj.h.

Effective source at the previous step for the resolution of the Poisson equation for dzeta .

Definition at line 559 of file compobj.h.

Effective source at the previous step for the resolution of the Poisson equation for the scalar $\chi$ by means of Map_et::poisson .

$\chi$ is an intermediate quantity for the resolution of the elliptic equation for the shift vector $N^i$

Definition at line 572 of file compobj.h.

Effective source at the previous step for the resolution of the Poisson equation for nuf by means of Map_et::poisson .

Definition at line 548 of file compobj.h.

Effective source at the previous step for the resolution of the Poisson equation for nuq by means of Map_et::poisson .

Definition at line 554 of file compobj.h.

Effective source at the previous step for the resolution of the Poisson equation for tggg .

Definition at line 564 of file compobj.h.

Effective source at the previous step for the resolution of the vector Poisson equation for $W^i$.

$W^i$ is an intermediate quantity for the resolution of the elliptic equation for the shift vector $N^i$ (Components with respect to the Cartesian triad associated with the mapping mp )

Definition at line 581 of file compobj.h.

Sym_tensor Lorene::Compobj::stress_euler [protected, inherited]

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

Definition at line 153 of file compobj.h.

Metric potential $\tilde G = (NB-1) r\sin\theta$.

Definition at line 519 of file compobj.h.

Component $\tilde N^\varphi = N^\varphi r\sin\theta$ of the shift vector.

Definition at line 503 of file compobj.h.

Vector $W^i$ used in the decomposition of shift , following Shibata's prescription [Prog.

Theor. Phys. 101 , 1199 (1999)] :

\[ N^i = {7\over 8} W^i - {1\over 8} \left(\nabla^i\chi+\nabla^iW^kx_k\right) \]

NB: w_shift contains the components of $W^i$ with respect to the Cartesian triad associated with the mapping mp .

Definition at line 532 of file compobj.h.


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

Generated on 7 Dec 2019 for LORENE by  doxygen 1.6.1