28char star_rot_C[] =
"$Header: /cvsroot/Lorene/C++/Source/Star/star_rot.C,v 1.9 2015/05/19 09:30:56 j_novak Exp $" ;
74#include "utilitaires.h"
99 w_shift(
mpi, CON, mp.get_bvect_cart()),
101 tkij(
mpi, COV, mp.get_bvect_cart()),
108 ssjm1_wshift(
mpi, CON, mp.get_bvect_cart())
158 relativistic(et.relativistic),
172 khi_shift(et.khi_shift),
175 ssjm1_nuf(et.ssjm1_nuf),
176 ssjm1_nuq(et.ssjm1_nuq),
177 ssjm1_dzeta(et.ssjm1_dzeta),
178 ssjm1_tggg(et.ssjm1_tggg),
179 ssjm1_khi(et.ssjm1_khi),
180 ssjm1_wshift(et.ssjm1_wshift)
201 w_shift(
mpi, CON, mp.get_bvect_cart()),
203 tkij(
mpi, COV, mp.get_bvect_cart()),
210 ssjm1_wshift(
mpi, CON, mp.get_bvect_cart())
439 ost <<
"Uniformly rotating star" <<
endl ;
440 ost <<
"-----------------------" <<
endl ;
444 <<
" rad/s f : " <<
freq * f_unit <<
" Hz" <<
endl ;
445 ost <<
"Rotation period : " << 1000. / (
freq * f_unit) <<
" ms"
450 ost <<
"Differentially rotating star" <<
endl ;
451 ost <<
"----------------------------" <<
endl ;
454 ost <<
"Central value of Omega : " <<
omega_c * f_unit
455 <<
" rad/s f : " <<
freq * f_unit <<
" Hz" <<
endl ;
456 ost <<
"Central rotation period : " << 1000. / (
freq * f_unit) <<
" ms"
460 ost <<
"Relativistic star" <<
endl ;
463 ost <<
"Newtonian star" <<
endl ;
476 ost <<
"Error on the virial identity GRV2 : " <<
grv2() <<
endl ;
478 ost <<
"Error on the virial identity GRV3 : " <<
xgrv3 <<
endl ;
484 ost <<
"Q / (M R_circ^2) : "
486 ost <<
"c^4 Q / (G^2 M^3) : "
490 ost <<
"Angular momentum J : "
491 <<
angu_mom()/( qpig / (4*
M_PI) * msol*msol) <<
" G M_sol^2 / c"
493 ost <<
"c J / (G M^2) : "
505 ost <<
"Circumferential equatorial radius R_circ : "
508 ost <<
"Surface area : " <<
area()/(km*km) <<
" km^2" <<
endl ;
511 ost <<
"Coordinate equatorial radius r_eq : " <<
ray_eq()/km <<
" km"
519 ost <<
"Equatorial value of the velocity U: "
521 ost <<
"Redshift at the equator (forward) : " <<
z_eqf() <<
endl ;
522 ost <<
"Redshift at the equator (backward): " <<
z_eqb() <<
endl ;
526 ost <<
"Central value of log(N) : "
529 ost <<
"Central value of dzeta=log(AN) : "
539 ost <<
" ... w_shift (NB: components in the star Cartesian frame) [c] : "
541 <<
w_shift(1).val_grid_point(0, 0, 0, 0) <<
" "
542 <<
w_shift(2).val_grid_point(0, 0, 0, 0) <<
" "
545 ost <<
"Central value of khi_shift [km c] : "
552 "Relative discrepancy in each domain between the metric coef. A^2 and B^2 : "
571 ost <<
endl <<
"Innermost stable circular orbit (ISCO) : " <<
endl ;
573 ost <<
" circumferential radius r_isco = "
575 ost <<
" (approx. 6M + 1st order in j : "
577 ost <<
" (approx. 6M : "
579 ost <<
" orbital frequency f_isco = "
581 ost <<
" (approx. 1st order in j : "
597 <<
" rad/s f : " <<
freq * f_unit <<
" Hz" <<
endl ;
598 ost <<
"Rotation period : " << 1000. / (
freq * f_unit) <<
" ms"
602 <<
" x 0.1 fm^-3" <<
endl ;
604 <<
" rho_nuc c^2" <<
endl ;
606 <<
" rho_nuc c^2" <<
endl ;
608 ost <<
"Central value of log(N) : "
611 ost <<
"Central value of dzeta=log(AN) : "
629 ost <<
"Equatorial value of the velocity U: "
633 <<
"Coordinate equatorial radius r_eq = "
678 ost <<
endl <<
"Quantities in polytropic units : " <<
endl ;
679 ost <<
"==============================" <<
endl ;
707 d_khi.dec_dzpuis(2) ;
726 xk.set(1).set_domain(nz-1) =
sintcosp.domain(nz-1) ;
727 xk.set(1).set_dzpuis(-1) ;
729 xk.set(2).set_domain(nz-1) =
sintsinp.domain(nz-1) ;
730 xk.set(2).set_dzpuis(-1) ;
732 xk.set(3).set_domain(nz-1) = cost.
domain(nz-1) ;
733 xk.set(3).set_dzpuis(-1) ;
734 xk.std_spectral_base() ;
750 if ( (
beta(1).get_etat() == ETATZERO) && (
beta(2).get_etat() == ETATZERO) ) {
Polytropic equation of state (relativistic case).
Equation of state base class.
Time evolution with partial storage (*** under development ***).
Base class for coordinate mappings.
const Base_vect_cart & get_bvect_cart() const
Returns the Cartesian basis associated with the coordinates (x,y,z) of the mapping,...
Coord y
y coordinate centered on the grid
virtual void comp_p_from_cartesian(const Scalar &v_x, const Scalar &v_y, Scalar &v_p) const =0
Computes the Spherical component (with respect to bvect_spher ) of a vector given by its cartesian c...
const Metric_flat & flat_met_cart() const
Returns the flat metric associated with the Cartesian coordinates and with components expressed in th...
Coord z
z coordinate centered on the grid
Coord x
x coordinate centered on the grid
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Tensor field of valence 0 (or component of a tensorial field).
virtual void sauve(FILE *) const
Save in a file.
virtual void std_spectral_base()
Sets the spectral bases of the Valeur va to the standard ones for a scalar field.
double val_grid_point(int l, int k, int j, int i) const
Returns the value of the field at a specified grid point.
void div_rsint()
Division by everywhere; dzpuis is not changed.
const Tbl & domain(int l) const
Read-only of the value in a given domain.
const Vector & derive_con(const Metric &gam) const
Returns the "contravariant" derivative of *this with respect to some metric , by raising the index of...
Class for isolated rotating stars.
double * p_angu_mom
Angular momentum.
virtual void display_poly(ostream &) const
Display in polytropic units.
virtual double mean_radius() const
Mean star radius from the area .
double * p_r_isco
Circumferential radius of the ISCO.
double * p_aplat
Flatening r_pole/r_eq.
Star_rot(Map &mp_i, int nzet_i, bool relat, const Eos &eos_i)
Standard constructor.
double * p_mom_quad
Quadrupole moment.
Scalar ssjm1_tggg
Effective source at the previous step for the resolution of the Poisson equation for tggg .
Sym_tensor tkij
Tensor related to the extrinsic curvature tensor by .
virtual double mom_quad() const
Quadrupole moment.
double * p_z_eqf
Forward redshift factor at equator.
Scalar tggg
Metric potential .
Scalar tnphi
Component of the shift vector.
void fait_shift()
Computes shift from w_shift and khi_shift according to Shibata's prescription [Prog.
virtual ~Star_rot()
Destructor.
virtual double z_pole() const
Redshift factor at North pole.
virtual double r_isco(ostream *ost=0x0) const
Circumferential radius of the innermost stable circular orbit (ISCO).
virtual double z_eqb() const
Backward redshift factor at equator.
double unsurc2
: unsurc2=1 for a relativistic star, 0 for a Newtonian one.
Vector w_shift
Vector used in the decomposition of shift , following Shibata's prescription [Prog.
Scalar b_car
Square of the metric factor B.
void operator=(const Star_rot &)
Assignment to another Star_rot.
double * p_lspec_isco
Specific angular momentum of a particle on the ISCO.
Scalar bbb
Metric factor B.
virtual double mass_b() const
Baryon mass.
virtual double get_omega_c() const
Returns the central value of the rotation angular velocity ([f_unit] )
double omega
Rotation angular velocity ([f_unit] )
virtual double r_circ() const
Circumferential radius.
Scalar uuu
Norm of u_euler.
virtual double z_eqf() const
Forward redshift factor at equator.
Scalar nphi
Metric coefficient .
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
double * p_f_isco
Orbital frequency of the ISCO.
virtual double f_isco() const
Orbital frequency at the innermost stable circular orbit (ISCO).
virtual double mass_g() const
Gravitational mass.
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
virtual void sauve(FILE *) const
Save in a file.
double * p_area
Integrated surface area.
virtual double aplat() const
Flatening r_pole/r_eq.
double * p_grv3
Error on the virial identity GRV3.
Scalar ssjm1_nuq
Effective source at the previous step for the resolution of the Poisson equation for nuq by means of ...
double * p_espec_isco
Specific energy of a particle on the ISCO.
Scalar nuf
Part of the Metric potential = logn generated by the matter terms.
bool relativistic
Indicator of relativity: true for a relativistic star, false for a Newtonian one.
double * p_f_eq
Orbital frequency at the equator.
Scalar ssjm1_dzeta
Effective source at the previous step for the resolution of the Poisson equation for dzeta .
virtual double tsw() const
Ratio T/W.
virtual void del_hydro_euler()
Sets to ETATNONDEF (undefined state) the hydrodynamical quantities relative to the Eulerian observer.
virtual double grv3(ostream *ost=0x0) const
Error on the virial identity GRV3.
Scalar ssjm1_khi
Effective source at the previous step for the resolution of the Poisson equation for the scalar by m...
Scalar dzeta
Metric potential .
virtual double area() const
Integrated surface area in .
Scalar a_car
Square of the metric factor A.
Scalar ssjm1_nuf
Effective source at the previous step for the resolution of the Poisson equation for nuf by means of ...
Scalar nuq
Part of the Metric potential = logn generated by the quadratic terms.
double * p_grv2
Error on the virial identity GRV2.
virtual void del_deriv() const
Deletes all the derived quantities.
Vector ssjm1_wshift
Effective source at the previous step for the resolution of the vector Poisson equation for .
double * p_r_circ
Circumferential radius.
virtual double grv2() const
Error on the virial identity GRV2.
virtual double angu_mom() const
Angular momentum.
double * p_z_eqb
Backward redshift factor at equator.
virtual void partial_display(ostream &) const
Printing of some informations, excluding all global quantities.
double * p_z_pole
Redshift factor at North pole.
Scalar khi_shift
Scalar used in the decomposition of shift , following Shibata's prescription [Prog.
void fait_nphi()
Computes tnphi and nphi from the Cartesian components of the shift, stored in shift .
Scalar ener
Total energy density in the fluid frame.
Scalar logn
Logarithm of the lapse N .
Scalar nn
Lapse function N .
virtual void del_hydro_euler()
Sets to ETATNONDEF (undefined state) the hydrodynamical quantities relative to the Eulerian observer.
const Eos & eos
Equation of state of the stellar matter.
Scalar nbar
Baryon density in the fluid frame.
virtual void del_deriv() const
Deletes all the derived quantities.
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
double ray_eq() const
Coordinate radius at , [r_unit].
virtual void sauve(FILE *) const
Save in a file.
Scalar press
Fluid pressure.
Map & mp
Mapping associated with the star.
void operator=(const Star &)
Assignment to another Star.
int nzet
Number of domains of *mp occupied by the star.
Tensor field of valence 1.
Cmp sqrt(const Cmp &)
Square root.
Tbl diffrel(const Cmp &a, const Cmp &b)
Relative difference between two Cmp (norme version).
Cmp pow(const Cmp &, int)
Power .
int fread_be(int *aa, int size, int nb, FILE *fich)
Reads integer(s) from a binary file according to the big endian convention.
int fwrite_be(const int *aa, int size, int nb, FILE *fich)
Writes integer(s) into a binary file according to the big endian convention.
virtual void sauve(FILE *) const
Save in a binary file.
void set_triad(const Base_vect &new_triad)
Assigns a new vectorial basis (triad) of decomposition.
const Tensor & derive_con(const Metric &gam) const
Returns the "contravariant" derivative of this with respect to some metric , by raising the last inde...
virtual void set_etat_zero()
Sets the logical state of all components to ETATZERO (zero state).
Tenseur contract(const Tenseur &, int id1, int id2)
Self contraction of two indices of a Tenseur .
Standard units of space, time and mass.