LORENE
Lorene::Boson_star Class Reference

Class for stationary axisymmetric boson stars (under development). More...

#include <boson_star.h>

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

Public Member Functions

 Boson_star (Map &mp_i, double m, int k)
 Standard constructor. More...
 
 Boson_star (const Boson_star &)
 Copy constructor. More...
 
 Boson_star (Map &mp_i, FILE *fich)
 Constructor from a file (see sauve(FILE*) ). More...
 
virtual ~Boson_star ()
 Destructor. More...
 
void operator= (const Boson_star &)
 Assignment to another Boson_star. More...
 
const Scalarget_rphi () const
 Returns the real part of the scalar field. More...
 
const Scalarget_iphi () const
 Returns the imaginary part of the scalar field. More...
 
Scalarset_rphi ()
 Sets a value to the real part of the scalar field. More...
 
Scalarset_iphi ()
 Sets a value to the imaginary part of the scalar field. More...
 
virtual void sauve (FILE *) const
 Save in a file. More...
 
void update_ener_mom ()
 Computes the 3+1 components of the energy-momentum tensor (E, P_i and S_{ij}) from the values of the scalar field and the metric. More...
 
virtual void equilibrium (double rphi_c, double iphi_c, int nzadapt, const Tbl &ent_limit, const Itbl &icontrol, const Tbl &control, Tbl &diff, Param *=0x0)
 Solves the equation satisfied by the scalar field. More...
 
const Scalarget_logn () const
 Returns the logarithm of the lapse N. More...
 
const Scalarget_tnphi () const
 Returns the component $\tilde N^\varphi = N^\varphi r\sin\theta$ of the shift vector. More...
 
const Scalarget_nuf () const
 Returns the part of the Metric potential $\nu = \ln N$ = logn generated by the matter terms. More...
 
const Scalarget_nuq () const
 Returns the Part of the Metric potential $\nu = \ln N$ = logn generated by the quadratic terms. More...
 
const Scalarget_dzeta () const
 Returns the Metric potential $\zeta = \ln(AN)$. More...
 
const Scalarget_tggg () const
 Returns the Metric potential $\tilde G = (NB-1) r\sin\theta$. More...
 
const Vectorget_w_shift () const
 Returns the vector $W^i$ used in the decomposition of shift , following Shibata's prescription [Prog. More...
 
const Scalarget_khi_shift () const
 Returns the scalar $\chi$ used in the decomposition of shift
following Shibata's prescription [Prog. More...
 
virtual double mass_g () const
 Gravitational mass. More...
 
virtual double angu_mom () const
 Angular momentum. More...
 
virtual double grv2 () const
 Error on the virial identity GRV2. More...
 
virtual double grv3 (ostream *ost=0x0) const
 Error on the virial identity GRV3. More...
 
virtual double mom_quad () const
 Quadrupole moment. More...
 
void update_metric ()
 Computes metric coefficients from known potentials. More...
 
void fait_shift ()
 Computes shift from w_shift and khi_shift according to Shibata's prescription [Prog. More...
 
void fait_nphi ()
 Computes tnphi and nphi from the Cartesian components of the shift, stored in shift . 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...
 
void gyoto_data (const char *file_name) const
 Save in a file for GYOTO. 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 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...
 

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. More...
 

Protected Member Functions

virtual void del_deriv () const
 Deletes all the derived quantities. More...
 
virtual 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 rphi
 Real part of the scalar field Phi. More...
 
Scalar iphi
 Imaginary part of the scalar field Phi. More...
 
double omega
 Coefficient omega in the time dependence of Phi. More...
 
int kkk
 Coefficient kkk in the azimuthal dependence of Phi. More...
 
double mmm
 Boson mass. More...
 
double m2
 Boson mass squared. More...
 
Scalar logn
 Logarithm of the lapse N . More...
 
Scalar tnphi
 Component $\tilde N^\varphi = N^\varphi r\sin\theta$ of the shift vector. More...
 
Scalar nuf
 Part of the Metric potential $\nu = \ln N$ = logn generated by the matter terms. More...
 
Scalar nuq
 Part of the Metric potential $\nu = \ln N$ = logn generated by the quadratic terms. More...
 
Scalar dzeta
 Metric potential $\zeta = \ln(AN)$. More...
 
Scalar tggg
 Metric potential $\tilde G = (NB-1) r\sin\theta$. More...
 
Vector w_shift
 Vector $W^i$ used in the decomposition of shift , following Shibata's prescription [Prog. More...
 
Scalar khi_shift
 Scalar $\chi$ used in the decomposition of shift , following Shibata's prescription [Prog. More...
 
Scalar ssjm1_nuf
 Effective source at the previous step for the resolution of the Poisson equation for nuf by means of Map_et::poisson . More...
 
Scalar ssjm1_nuq
 Effective source at the previous step for the resolution of the Poisson equation for nuq by means of Map_et::poisson . More...
 
Scalar ssjm1_dzeta
 Effective source at the previous step for the resolution of the Poisson equation for dzeta . More...
 
Scalar ssjm1_tggg
 Effective source at the previous step for the resolution of the Poisson equation for tggg . More...
 
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 . More...
 
Vector ssjm1_wshift
 Effective source at the previous step for the resolution of the vector Poisson equation for $W^i$. More...
 
double * p_grv2
 Error on the virial identity GRV2. More...
 
double * p_grv3
 Error on the virial identity GRV3. More...
 
double * p_mom_quad
 Quadrupole moment. More...
 
double * p_mass_g
 Gravitational mass (ADM mass as a volume integral) More...
 
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

Class for stationary axisymmetric boson stars (under development).

()

Definition at line 64 of file boson_star.h.

Constructor & Destructor Documentation

◆ Boson_star() [1/3]

Lorene::Boson_star::Boson_star ( Map mp_i,
double  m,
int  k 
)

Standard constructor.

Parameters
mp_iMapping on which the star is contructed
mBoson mass //## which unit ?
kCoefficient k in the azimuthal dependence of Phi

Definition at line 67 of file boson_star.C.

References iphi, rphi, and set_der_0x0().

◆ Boson_star() [2/3]

Lorene::Boson_star::Boson_star ( const Boson_star st)

Copy constructor.

Definition at line 87 of file boson_star.C.

References set_der_0x0().

◆ Boson_star() [3/3]

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

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

Parameters
mp_iMapping on which the star is constructed
fichinput file (must have been created by the function Boson_star::sauve )

Definition at line 103 of file boson_star.C.

References Lorene::fread_be(), kkk, m2, mmm, omega, and set_der_0x0().

◆ ~Boson_star()

Lorene::Boson_star::~Boson_star ( )
virtual

Destructor.

Definition at line 122 of file boson_star.C.

References del_deriv().

Member Function Documentation

◆ adm_mass()

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

◆ angu_mom()

double Lorene::Star_QI::angu_mom ( ) const
virtualinherited

◆ del_deriv()

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

Deletes all the derived quantities.

Reimplemented from Lorene::Star_QI.

Definition at line 133 of file boson_star.C.

References Lorene::Star_QI::del_deriv(), and set_der_0x0().

◆ equilibrium()

void Lorene::Boson_star::equilibrium ( double  rphi_c,
double  iphi_c,
int  nzadapt,
const Tbl ent_limit,
const Itbl icontrol,
const Tbl control,
Tbl diff,
Param = 0x0 
)
virtual

Solves the equation satisfied by the scalar field.

Computes an equilibrium configuration.

Parameters
rphi_c[input] Central value of the real part of the scalar field
iphi_c[input] Central value of the imaginary part of the scalar field
nzadapt[input] Number of (inner) domains where the mapping adaptation to an iso-enthalpy surface should be performed
ent_limit[input] 1-D Tbl of dimension nzet which defines the enthalpy at the outer boundary of each domain
icontrol[input] Set of integer parameters (stored as a 1-D Itbl of size 8) to control the iteration:
  • icontrol(0) = mer_max : maximum number of steps
  • icontrol(1) = mer_rot : step at which the rotation is switched on
  • icontrol(2) = mer_change_omega : step at which the rotation velocity is changed to reach the final one
  • icontrol(3) = mer_fix_omega : step at which the final rotation velocity must have been reached
  • icontrol(4) = mer_mass : the absolute value of mer_mass is the step from which the baryon mass is forced to converge, by varying the central enthalpy (mer_mass>0 ) or the angular velocity (mer_mass<0 )
  • icontrol(5) = mermax_poisson : maximum number of steps in Map_et::poisson
  • icontrol(6) = mer_triax : step at which the 3-D perturbation is switched on
  • icontrol(7) = delta_mer_kep : number of steps after mer_fix_omega when omega starts to be increased by fact_omega to search for the Keplerian velocity
control[input] Set of parameters (stored as a 1-D Tbl of size 7) to control the iteration:
  • control(0) = precis : threshold on the enthalpy relative change for ending the computation
  • control(1) = omega_ini : initial angular velocity, switched on only if mer_rot<0 , otherwise 0 is used
  • control(2) = relax : relaxation factor in the main iteration
  • control(3) = relax_poisson : relaxation factor in Map_et::poisson
  • control(4) = thres_adapt : threshold on dH/dr for freezing the adaptation of the mapping
  • control(5) = ampli_triax : relative amplitude of the 3-D perturbation
  • control(6) = precis_adapt : precision for Map_et::adapt
diff[output] 1-D Tbl of size 7 for the storage of some error indicators :
  • diff(0) : Relative change in the enthalpy field between two successive steps
  • diff(1) : Relative error in the resolution of the Poisson equation for nuf
  • diff(2) : Relative error in the resolution of the Poisson equation for nuq
  • diff(3) : Relative error in the resolution of the Poisson equation for dzeta
  • diff(4) : Relative error in the resolution of the Poisson equation for tggg
  • diff(5) : Relative error in the resolution of the equation for shift (x comp.)
  • diff(6) : Relative error in the resolution of the equation for shift (y comp.)

Definition at line 73 of file boson_star_equil.C.

References Lorene::Map::get_mg(), Lorene::Mg3d::get_nzone(), Lorene::Mg3d::get_type_t(), and Lorene::Compobj::mp.

◆ espec_isco()

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

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().

◆ extrinsic_curvature()

◆ f_isco()

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

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().

◆ fait_nphi()

void Lorene::Star_QI::fait_nphi ( )
inherited

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.

◆ fait_shift()

◆ get_a_car()

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

Returns the square of the metric factor A.

Definition at line 381 of file compobj.h.

References Lorene::Compobj_QI::a_car.

◆ get_ak_car()

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

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.

◆ get_b_car()

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

Returns the square of the metric factor B.

Definition at line 384 of file compobj.h.

References Lorene::Compobj_QI::b_car.

◆ get_bbb()

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

Returns the metric factor B.

Definition at line 378 of file compobj.h.

References Lorene::Compobj_QI::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_dzeta()

const Scalar& Lorene::Star_QI::get_dzeta ( ) const
inlineinherited

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

Definition at line 658 of file compobj.h.

References Lorene::Star_QI::dzeta.

◆ 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_iphi()

const Scalar& Lorene::Boson_star::get_iphi ( ) const
inline

Returns the imaginary part of the scalar field.

Definition at line 150 of file boson_star.h.

References iphi.

◆ get_khi_shift()

const Scalar& Lorene::Star_QI::get_khi_shift ( ) const
inlineinherited

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 Lorene::Star_QI::khi_shift.

◆ 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_logn()

const Scalar& Lorene::Star_QI::get_logn ( ) const
inlineinherited

Returns the logarithm of the lapse N.

Definition at line 639 of file compobj.h.

References Lorene::Star_QI::logn.

◆ 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
inlineinherited

Returns the metric coefficient $N^\varphi$.

Definition at line 387 of file compobj.h.

References Lorene::Compobj_QI::nphi.

◆ get_nuf()

const Scalar& Lorene::Star_QI::get_nuf ( ) const
inlineinherited

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 Lorene::Star_QI::nuf.

◆ get_nuq()

const Scalar& Lorene::Star_QI::get_nuq ( ) const
inlineinherited

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 Lorene::Star_QI::nuq.

◆ get_rphi()

const Scalar& Lorene::Boson_star::get_rphi ( ) const
inline

Returns the real part of the scalar field.

Definition at line 146 of file boson_star.h.

References rphi.

◆ 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.

◆ get_tggg()

const Scalar& Lorene::Star_QI::get_tggg ( ) const
inlineinherited

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

Definition at line 661 of file compobj.h.

References Lorene::Star_QI::tggg.

◆ get_tnphi()

const Scalar& Lorene::Star_QI::get_tnphi ( ) const
inlineinherited

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

Definition at line 645 of file compobj.h.

References Lorene::Star_QI::tnphi.

◆ get_w_shift()

const Vector& Lorene::Star_QI::get_w_shift ( ) const
inlineinherited

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 Lorene::Star_QI::w_shift.

◆ grv2()

double Lorene::Star_QI::grv2 ( ) const
virtualinherited

◆ grv3()

double Lorene::Star_QI::grv3 ( ostream *  ost = 0x0) const
virtualinherited

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
ostoutput 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::ak_car, Lorene::Compobj_QI::bbb, Lorene::Scalar::derive_cov(), Lorene::Star_QI::dzeta, Lorene::Map::flat_met_spher(), Lorene::Scalar::get_etat(), Lorene::log(), Lorene::Star_QI::logn, Lorene::Compobj::mp, Lorene::Star_QI::p_grv3, Lorene::Scalar::srdsdt(), and Lorene::Scalar::std_spectral_base().

◆ gyoto_data()

◆ lambda_grv2()

double Lorene::Star_QI::lambda_grv2 ( const Scalar sou_m,
const Scalar sou_q 
)
staticinherited

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::Scalar::check_dzpuis(), Lorene::Map::get_mg(), Lorene::Tensor::get_mp(), Lorene::Mg3d::get_nzone(), Lorene::Mg3d::get_type_r(), and Lorene::Map::val_r().

◆ lspec_isco()

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

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().

◆ mass_g()

◆ mom_quad()

double Lorene::Star_QI::mom_quad ( ) const
virtualinherited

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::log(), Lorene::Star_QI::logn, Lorene::Compobj::mp, Lorene::Scalar::mult_r(), Lorene::Star_QI::p_mom_quad, Lorene::Scalar::std_spectral_base(), Lorene::Compobj::stress_euler, and Lorene::Tensor::trace().

◆ operator=()

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

Assignment to another Boson_star.

Definition at line 151 of file boson_star.C.

References del_deriv(), iphi, kkk, m2, mmm, omega, Lorene::Star_QI::operator=(), and rphi.

◆ operator>>()

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

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

Reimplemented from Lorene::Star_QI.

Definition at line 186 of file boson_star.C.

References iphi, kkk, Lorene::Star_QI::mass_g(), mmm, omega, Lorene::Star_QI::operator>>(), rphi, and Lorene::Scalar::val_grid_point().

◆ r_isco()

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

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(), 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().

◆ r_mb()

◆ sauve()

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

Save in a file.

Reimplemented from Lorene::Star_QI.

Definition at line 172 of file boson_star.C.

References Lorene::fwrite_be(), iphi, kkk, mmm, omega, rphi, Lorene::Star_QI::sauve(), and Lorene::Scalar::sauve().

◆ set_der_0x0()

void Lorene::Boson_star::set_der_0x0 ( ) const
protectedvirtual

Sets to 0x0 all the pointers on derived quantities.

Reimplemented from Lorene::Star_QI.

Definition at line 141 of file boson_star.C.

◆ set_iphi()

Scalar& Lorene::Boson_star::set_iphi ( )
inline

Sets a value to the imaginary part of the scalar field.

Definition at line 158 of file boson_star.h.

References del_deriv(), and iphi.

◆ 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.

◆ set_rphi()

Scalar& Lorene::Boson_star::set_rphi ( )
inline

Sets a value to the real part of the scalar field.

Definition at line 154 of file boson_star.h.

References del_deriv(), and rphi.

◆ update_ener_mom()

void Lorene::Boson_star::update_ener_mom ( )

◆ update_metric()

void Lorene::Star_QI::update_metric ( )
virtualinherited

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, Lorene::Star_QI::del_deriv(), Lorene::Scalar::div_rsint(), Lorene::Star_QI::dzeta, Lorene::exp(), Lorene::Star_QI::logn, Lorene::Compobj::nn, Lorene::Scalar::std_spectral_base(), Lorene::Star_QI::tggg, and Lorene::Compobj_QI::update_metric().

Member Data Documentation

◆ a_car

Scalar Lorene::Compobj_QI::a_car
protectedinherited

Square of the metric factor A.

Definition at line 290 of file compobj.h.

◆ ak_car

Scalar Lorene::Compobj_QI::ak_car
protectedinherited

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
protectedinherited

Square of the metric factor B.

Definition at line 296 of file compobj.h.

◆ bbb

Scalar Lorene::Compobj_QI::bbb
protectedinherited

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.

◆ dzeta

Scalar Lorene::Star_QI::dzeta
protectedinherited

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

Definition at line 516 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.

◆ iphi

Scalar Lorene::Boson_star::iphi
protected

Imaginary part of the scalar field Phi.

Definition at line 76 of file boson_star.h.

◆ khi_shift

Scalar Lorene::Star_QI::khi_shift
protectedinherited

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.

◆ kk

Sym_tensor Lorene::Compobj::kk
protectedinherited

Extrinsic curvature tensor $K_{ij}$.

Definition at line 156 of file compobj.h.

◆ kkk

int Lorene::Boson_star::kkk
protected

Coefficient kkk in the azimuthal dependence of Phi.

Definition at line 84 of file boson_star.h.

◆ logn

Scalar Lorene::Star_QI::logn
protectedinherited

Logarithm of the lapse N .

Definition at line 498 of file compobj.h.

◆ m2

double Lorene::Boson_star::m2
protected

Boson mass squared.

Definition at line 92 of file boson_star.h.

◆ mmm

double Lorene::Boson_star::mmm
protected

Boson mass.

Definition at line 88 of file boson_star.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
protectedinherited

Metric coefficient $N^\varphi$.

Definition at line 299 of file compobj.h.

◆ nuf

Scalar Lorene::Star_QI::nuf
protectedinherited

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

Definition at line 508 of file compobj.h.

◆ nuq

Scalar Lorene::Star_QI::nuq
protectedinherited

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

Definition at line 513 of file compobj.h.

◆ omega

double Lorene::Boson_star::omega
protected

Coefficient omega in the time dependence of Phi.

Definition at line 80 of file boson_star.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
mutableprotectedinherited

Angular momentum.

Definition at line 324 of file compobj.h.

◆ p_espec_isco

double* Lorene::Compobj_QI::p_espec_isco
mutableprotectedinherited

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
mutableprotectedinherited

Orbital frequency of the ISCO.

Definition at line 326 of file compobj.h.

◆ p_grv2

double* Lorene::Star_QI::p_grv2
mutableprotectedinherited

Error on the virial identity GRV2.

Definition at line 588 of file compobj.h.

◆ p_grv3

double* Lorene::Star_QI::p_grv3
mutableprotectedinherited

Error on the virial identity GRV3.

Definition at line 589 of file compobj.h.

◆ p_lspec_isco

double* Lorene::Compobj_QI::p_lspec_isco
mutableprotectedinherited

Specific angular momentum of a particle at the ISCO.

Definition at line 330 of file compobj.h.

◆ p_mass_g

double* Lorene::Star_QI::p_mass_g
mutableprotectedinherited

Gravitational mass (ADM mass as a volume integral)

Definition at line 591 of file compobj.h.

◆ p_mom_quad

double* Lorene::Star_QI::p_mom_quad
mutableprotectedinherited

Quadrupole moment.

Definition at line 590 of file compobj.h.

◆ p_r_isco

double* Lorene::Compobj_QI::p_r_isco
mutableprotectedinherited

Coordinate r of the ISCO.

Definition at line 325 of file compobj.h.

◆ p_r_mb

double* Lorene::Compobj_QI::p_r_mb
mutableprotectedinherited

Coordinate r of the marginally bound orbit.

Definition at line 331 of file compobj.h.

◆ rphi

Scalar Lorene::Boson_star::rphi
protected

Real part of the scalar field Phi.

Definition at line 72 of file boson_star.h.

◆ ssjm1_dzeta

Scalar Lorene::Star_QI::ssjm1_dzeta
protectedinherited

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

Definition at line 559 of file compobj.h.

◆ ssjm1_khi

Scalar Lorene::Star_QI::ssjm1_khi
protectedinherited

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.

◆ ssjm1_nuf

Scalar Lorene::Star_QI::ssjm1_nuf
protectedinherited

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.

◆ ssjm1_nuq

Scalar Lorene::Star_QI::ssjm1_nuq
protectedinherited

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.

◆ ssjm1_tggg

Scalar Lorene::Star_QI::ssjm1_tggg
protectedinherited

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

Definition at line 564 of file compobj.h.

◆ ssjm1_wshift

Vector Lorene::Star_QI::ssjm1_wshift
protectedinherited

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.

◆ 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.

◆ tggg

Scalar Lorene::Star_QI::tggg
protectedinherited

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

Definition at line 519 of file compobj.h.

◆ tnphi

Scalar Lorene::Star_QI::tnphi
protectedinherited

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

Definition at line 503 of file compobj.h.

◆ w_shift

Vector Lorene::Star_QI::w_shift
protectedinherited

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: