LORENE
et_rot_bifluid.h
1/*
2 * Definition of Lorene class Et_rot_bifluid
3 *
4 */
5
6/*
7 * Copyright (c) 2001 Jerome Novak
8 * (c) 2015 Aurelien Sourie
9 *
10 * This file is part of LORENE.
11 *
12 * LORENE is free software; you can redistribute it and/or modify
13 * it under the terms of the GNU General Public License as published by
14 * the Free Software Foundation; either version 2 of the License, or
15 * (at your option) any later version.
16 *
17 * LORENE is distributed in the hope that it will be useful,
18 * but WITHOUT ANY WARRANTY; without even the implied warranty of
19 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
20 * GNU General Public License for more details.
21 *
22 * You should have received a copy of the GNU General Public License
23 * along with LORENE; if not, write to the Free Software
24 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
25 *
26 */
27
28
29#ifndef __ET_ROT_BIFLUID_H_
30#define __ET_ROT_BIFLUID_H_
31
32/*
33 * $Id: et_rot_bifluid.h,v 1.20 2015/06/26 14:10:08 j_novak Exp $
34 * $Log: et_rot_bifluid.h,v $
35 * Revision 1.20 2015/06/26 14:10:08 j_novak
36 * Modified comments.
37 *
38 * Revision 1.19 2015/06/11 13:50:18 j_novak
39 * Minor corrections
40 *
41 * Revision 1.18 2015/06/10 14:39:17 a_sourie
42 * New class Eos_bf_tabul for tabulated 2-fluid EoSs and associated functions for the computation of rotating stars with such EoSs.
43 *
44 * Revision 1.17 2014/10/13 08:52:34 j_novak
45 * Lorene classes and functions now belong to the namespace Lorene.
46 *
47 * Revision 1.16 2013/11/25 13:50:55 j_novak
48 * The inheritance from Etoile_rot is no longer virtual.
49 *
50 * Revision 1.15 2011/10/06 14:55:36 j_novak
51 * equation_of_state() is now virtual to be able to call to the magnetized
52 * Eos_mag.
53 *
54 * Revision 1.14 2004/09/01 10:56:05 r_prix
55 * added option of converging baryon-mass to equilibrium_bi()
56 *
57 * Revision 1.13 2004/03/22 13:12:41 j_novak
58 * Modification of comments to use doxygen instead of doc++
59 *
60 * Revision 1.12 2003/12/04 14:13:32 r_prix
61 * added method get_typeos {return typeos}; and fixed some comments.
62 *
63 * Revision 1.11 2003/11/20 14:01:45 r_prix
64 * changed member names to better conform to Lorene coding standards:
65 * J_euler -> j_euler, EpS_euler -> enerps_euler, Delta_car -> delta_car
66 *
67 * Revision 1.10 2003/11/18 18:32:36 r_prix
68 * added new class-member: EpS_euler := ener_euler + s_euler
69 * has the advantage of a nice Newtonian limit -> rho
70 * (ener_euler is no longer used in this class!)
71 *
72 * Revision 1.9 2003/11/13 12:02:03 r_prix
73 * - adapted/extended some of the documentation
74 * - changed xxx2 -> Delta_car
75 * - added members J_euler, sphph_euler, representing 3+1 components of Tmunu
76 * (NOTE: these are not 2-fluid specific, and should ideally be moved into Class Etoile!)
77 *
78 * Revision 1.8 2003/09/17 08:27:50 j_novak
79 * New methods: mass_b1() and mass_b2().
80 *
81 * Revision 1.7 2002/10/09 07:54:29 j_novak
82 * Et_rot_bifluid and Et_rot_mag inheritate virtually from Etoile_rot
83 *
84 * Revision 1.6 2002/09/13 09:17:33 j_novak
85 * Modif. commentaires
86 *
87 * Revision 1.5 2002/04/05 09:09:36 j_novak
88 * The inversion of the EOS for 2-fluids polytrope has been modified.
89 * Some errors in the determination of the surface were corrected.
90 *
91 * Revision 1.4 2002/01/16 15:03:28 j_novak
92 * *** empty log message ***
93 *
94 * Revision 1.3 2002/01/08 14:43:53 j_novak
95 * better determination of surfaces for 2-fluid stars
96 *
97 * Revision 1.2 2002/01/03 15:30:27 j_novak
98 * Some comments modified.
99 *
100 * Revision 1.1.1.1 2001/11/20 15:19:27 e_gourgoulhon
101 * LORENE
102 *
103 * Revision 1.3 2001/10/03 09:49:06 novak
104 * *** empty log message ***
105 *
106 * Revision 1.2 2001/08/28 14:14:10 novak
107 * overrided l_surf function
108 *
109 * Revision 1.1 2001/06/22 15:38:52 novak
110 * Initial revision
111 *
112 *
113 * $Header: /cvsroot/Lorene/C++/Include/et_rot_bifluid.h,v 1.20 2015/06/26 14:10:08 j_novak Exp $
114 *
115 */
116
117// Headers Lorene
118#include "eos_bifluid.h"
119#include "etoile.h"
120
121namespace Lorene {
122
123 // Local prototype (for determining the surface)
124 Cmp prolonge_c1(const Cmp& uu, const int nzet) ;
125
147 class Et_rot_bifluid : public Etoile_rot {
148
149 // Data :
150 // -----
151 protected:
152 const Eos_bifluid& eos ;
153
154 double omega2 ;
155
156 // Fluid quantities with respect to the fluid frame
157 // ------------------------------------------------
158
161
163
167
168 // Fluid quantities with respect to the Eulerian frame
169 // ---------------------------------------------------
170
171 // FIXME: the following three variables are not specific to 2-fluid stars
172 // and should ideally be moved to class Etoile!
173
176
183
186
189
190 Tenseur j_euler21_1;// To compute epsN (1st version)
191 Tenseur j_euler22_1;// To compute epsP (1st version)
192
193 Tenseur j_euler11_2;//To compute In (2nd version)
194 Tenseur j_euler12_2;//To compute Ip (2nd version)
195
196 Tenseur j_euler21_2;// To compute epsN (2nd version)
197 Tenseur j_euler22_2;// To compute epsP (2nd version)
198
199
202
205
208
214
215 // Derived data :
216 // ------------
217 protected:
219 mutable double* p_ray_eq2 ;
220
222 mutable double* p_ray_eq2_pis2 ;
223
225 mutable double* p_ray_eq2_pi ;
226
228 mutable double* p_ray_pole2 ;
229
234 mutable Itbl* p_l_surf2 ;
235
240 mutable Tbl* p_xi_surf2 ;
241
242 mutable double* p_r_circ2 ;
243 mutable double* p_area2 ;
244 mutable double* p_aplat2 ;
245
246 mutable double* p_mass_b1 ;
247 mutable double* p_mass_b2 ;
248
249 mutable double* p_angu_mom_1;
250 mutable double* p_angu_mom_2;
251
252 mutable double* p_angu_mom_1_part1_1;
253 mutable double* p_angu_mom_2_part1_1;
254
255 mutable double* p_angu_mom_1_part2_1;
256 mutable double* p_angu_mom_2_part2_1;
257
258 mutable double* p_angu_mom_1_part1_2;
259 mutable double* p_angu_mom_2_part1_2;
260
261 mutable double* p_angu_mom_1_part2_2;
262 mutable double* p_angu_mom_2_part2_2;
263
264 // Constructors - Destructor
265 // -------------------------
266 public:
267
268 Et_rot_bifluid(Map& mp_i, int nzet_i, bool relat,
269 const Eos_bifluid& eos_i) ;
270
272
277 Et_rot_bifluid(Map& mp_i, const Eos_bifluid& eos_i, FILE* fich) ;
278
279 virtual ~Et_rot_bifluid() ;
280
281 // Memory management
282 // -----------------
283 protected:
284
286 virtual void del_deriv() const ;
287
289 virtual void set_der_0x0() const ;
290
294 virtual void del_hydro_euler() ;
295
296 // Mutators / assignment
297 // ---------------------
298 public:
299
301 void operator=(const Et_rot_bifluid&) ;
302
304 void set_enthalpies(const Cmp&, const Cmp&) ;
305
314 void equilibrium_spher_bi(double ent_c, double ent_c2,
315 double precis = 1.e-14) ;
316
327 void equil_spher_regular(double ent_c, double ent_c2,
328 double precis = 1.e-14) ;
329
330 // Accessors
331 // ---------
332 public:
333
335 const Eos_bifluid& get_eos() const {return eos; } ;
336
338 const Tenseur& get_ent2() const {return ent2 ; } ;
339
341 const Tenseur& get_nbar2() const {return nbar2 ; } ;
342
344 const Tenseur& get_K_nn() const {return K_nn ; } ;
346 const Tenseur& get_K_np() const {return K_np ; } ;
348 const Tenseur& get_K_pp() const {return K_pp ; } ;
349
351 const Tenseur& get_delta_car() const {return delta_car ; } ;
352
354 const Tenseur& get_gam_euler2() const {return gam_euler2 ; } ;
355
357 double get_omega2() const {return omega2 ; } ;
358
360 const Tenseur& get_uuu2() const {return uuu2 ; } ;
361
362 // Outputs
363 // -------
364 public:
365 virtual void sauve(FILE *) const ;
366
368 virtual ostream& operator>>(ostream& ) const ;
369
371 virtual void partial_display(ostream& ) const ;
372
373 // Global quantities
374 // -----------------
375 public:
376
384 virtual const Itbl& l_surf() const ;
385
393 const Itbl& l_surf2() const ;
394
402 const Tbl& xi_surf2() const ;
403
405 double ray_eq2() const ;
406
408 double ray_eq2_pis2() const ;
409
411 double ray_eq2_pi() const ;
412
414 double ray_pole2() const ;
415
417 double mass_b1() const ;
418
420 double mass_b2() const ;
421
422 virtual double mass_b() const ;
423 virtual double mass_g() const ;
424 virtual double angu_mom() const ;
425
430 virtual double grv2() const ;
431
443 virtual double grv3(ostream* ost = 0x0) const ;
444
445 virtual double r_circ2() const ;
446 virtual double area2() const ;
447 virtual double mean_radius2() const ;
448 virtual double aplat2() const ;
449
464 virtual double mom_quad() const ;
465
473 virtual double mom_quad_old() const ;
474
481 virtual double mom_quad_Bo() const ;
482
483 virtual double angu_mom_1() const ;
484 virtual double angu_mom_2() const ;
485 virtual double angu_mom_1_part1_1() const ;
486 virtual double angu_mom_2_part1_1() const ;
487 virtual double angu_mom_1_part2_1() const ;
488 virtual double angu_mom_2_part2_1() const ;
489 virtual double angu_mom_1_part1_2() const ;
490 virtual double angu_mom_2_part1_2() const ;
491 virtual double angu_mom_1_part2_2() const ;
492 virtual double angu_mom_2_part2_2() const ;
493
494 // Computational routines
495 // ----------------------
496 public:
507 virtual void hydro_euler() ;
508
513 virtual void equation_of_state() ;
514
573 void equilibrium_bi(double ent_c, double ent_c2, double omega0,
574 double omega20, const Tbl& ent_limit,
575 const Tbl& ent2_limit, const Itbl& icontrol,
576 const Tbl& control, Tbl& diff,
577 int mer_mass, double mbar1_wanted, double mbar2_wanted,
578 double aexp_mass);
579 };
580
581}
582#endif
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition cmp.h:446
2-fluids equation of state base class.
Class for two-fluid rotating relativistic stars.
const Tenseur & get_K_pp() const
Returns the coefficient Kpp.
virtual void del_hydro_euler()
Sets to ETATNONDEF (undefined state) the hydrodynamical quantities relative to the Eulerian observer.
double * p_angu_mom_1_part2_2
To compute Xn (2nd version)
void equilibrium_bi(double ent_c, double ent_c2, double omega0, double omega20, const Tbl &ent_limit, const Tbl &ent2_limit, const Itbl &icontrol, const Tbl &control, Tbl &diff, int mer_mass, double mbar1_wanted, double mbar2_wanted, double aexp_mass)
Computes an equilibrium configuration.
Tenseur gam_euler2
Lorentz factor between the fluid 2 and Eulerian observers
Tenseur enerps_euler
the combination : useful because in the Newtonian limit .
Tenseur K_np
Coefficient Knp.
virtual ~Et_rot_bifluid()
Destructor.
double * p_mass_b1
Baryon mass of fluid 1.
virtual double angu_mom_1_part1_1() const
To compute In (1st version)
double * p_mass_b2
Baryon mass of fluid 2.
virtual double grv3(ostream *ost=0x0) const
Error on the virial identity GRV3.
virtual const Itbl & l_surf() const
Description of the surface of fluid 1: returns a 2-D Itbl containing the values of the domain index...
Tenseur j_euler12_1
To compute Ip (1st version)
const Tenseur & get_ent2() const
Returns the enthalpy field for fluid 2.
virtual double mom_quad_Bo() const
Part of the quadrupole moment.
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
void set_enthalpies(const Cmp &, const Cmp &)
Sets both enthalpy profiles.
double * p_angu_mom_1_part1_2
To compute In (2nd version)
Tenseur j_euler2
To compute Jp.
Tenseur sphph_euler
The component of the stress tensor .
const Tenseur & get_nbar2() const
Returns the proper baryon density for fluid 2.
double * p_angu_mom_2_part2_1
To compute Xp (1st version)
virtual double area2() const
Surface area for fluid 2.
Tenseur K_pp
Coefficient Kpp.
void operator=(const Et_rot_bifluid &)
Assignment to another Et_rot_bifluid.
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
const Eos_bifluid & get_eos() const
Returns the equation of state.
virtual double angu_mom_2_part1_2() const
To compute Ip (2nd version)
Tbl * p_xi_surf2
Description of the surface of fluid 2: 2-D Tbl containing the values of the radial coordinate on the...
virtual double angu_mom_2_part1_1() const
To compute Ip (1st version)
double * p_angu_mom_1
Angular momentum of fluid 1.
virtual double angu_mom_2() const
Angular momentum of fluid 2.
Tenseur j_euler11_1
To compute In (1st version)
virtual void del_deriv() const
Deletes all the derived quantities.
virtual double angu_mom_2_part2_2() const
To compute Xp (2nd version)
double * p_ray_eq2
Coordinate radius at , .
virtual void sauve(FILE *) const
Save in a file.
const Tenseur & get_K_np() const
Returns the coefficient Knp.
double * p_angu_mom_1_part2_1
To compute Xn (1st version)
double ray_pole2() const
Coordinate radius for fluid 2 at [r_unit].
void equilibrium_spher_bi(double ent_c, double ent_c2, double precis=1.e-14)
Computes a spherical static configuration.
const Tenseur & get_K_nn() const
Returns the coefficient Knn.
virtual double mass_b() const
Total Baryon mass.
Tenseur j_euler1
To compute Jn.
Tenseur j_euler
Total angular momentum (flat-space!) 3-vector , which is related to of the "3+1" decomposition,...
virtual void equation_of_state()
Computes the proper baryon and energy densities, as well as pressure and the coefficients Knn,...
const Tenseur & get_gam_euler2() const
Returns the Lorentz factor between the fluid 2 and Eulerian observers.
double get_omega2() const
Returns the rotation angular velocity of fluid 2([f_unit] )
const Tbl & xi_surf2() const
Description of the surface of fluid 2: returns a 2-D Tbl containing the values of the radial coordi...
double * p_ray_pole2
Coordinate radius at .
Tenseur K_nn
Coefficient Knn.
virtual double mass_g() const
Gravitational mass.
virtual double angu_mom_2_part2_1() const
To compute Xp (1st version)
double * p_aplat2
Flatening r_pole/r_eq of fluid no.2.
double * p_angu_mom_2_part2_2
To compute Xp (2nd version)
void equil_spher_regular(double ent_c, double ent_c2, double precis=1.e-14)
Computes a spherical static configuration.
double * p_area2
Surface area of fluid no.2.
double * p_angu_mom_2
Angular momentum of fluid 2.
Tenseur delta_car
The "relative velocity" (squared) of the two fluids.
const Itbl & l_surf2() const
Description of the surface of fluid 2: returns a 2-D Itbl containing the values of the domain index...
double omega2
Rotation angular velocity for fluid 2 ([f_unit] )
const Eos_bifluid & eos
Equation of state for two-fluids model.
double * p_angu_mom_2_part1_2
To compute Ip (2nd version)
double mass_b1() const
Baryon mass of fluid 1.
virtual double angu_mom_1_part1_2() const
To compute In (2nd version)
virtual double angu_mom_1_part2_1() const
To compute Xn (1st version)
Tenseur uuu2
Norm of the (fluid no.2) 3-velocity with respect to the eulerian observer.
double ray_eq2() const
Coordinate radius for fluid 2 at , [r_unit].
virtual double r_circ2() const
Circumferential radius for fluid 2.
double mass_b2() const
Baryon mass of fluid 2.
virtual double mom_quad() const
Quadrupole moment.
virtual double aplat2() const
Flatening r_pole/r_eq for fluid 2.
virtual double grv2() const
Error on the virial identity GRV2.
double * p_angu_mom_1_part1_1
To compute In (1st version)
double * p_r_circ2
Circumferential radius of fluid no.2.
Tenseur ent2
Log-enthalpy for the second fluid.
double * p_ray_eq2_pi
Coordinate radius at , .
virtual double angu_mom_1_part2_2() const
To compute Xn (2nd version)
const Tenseur & get_uuu2() const
Returns the norm of the fluid 2 3-velocity with respect to the eulerian frame.
virtual void partial_display(ostream &) const
Printing of some informations, excluding all global quantities.
double * p_angu_mom_2_part1_1
To compute Ip (1st version)
Itbl * p_l_surf2
Description of the surface of fluid 2: 2-D Itbl containing the values of the domain index l on the su...
Tenseur nbar2
Baryon density in the fluid frame, for fluid 2.
virtual double angu_mom_1() const
Angular momentum of fluid 1.
virtual double angu_mom() const
Angular momentum.
double * p_ray_eq2_pis2
Coordinate radius at , .
double ray_eq2_pis2() const
Coordinate radius for fluid 2 at , [r_unit].
double ray_eq2_pi() const
Coordinate radius for fluid 2 at , [r_unit].
const Tenseur & get_delta_car() const
Returns the "relative velocity" (squared) of the two fluids.
virtual void hydro_euler()
Computes the hydrodynamical quantities relative to the Eulerian observer from those in the fluid fram...
virtual double mean_radius2() const
Mean radius for fluid 2.
virtual double mom_quad_old() const
Part of the quadrupole moment.
Class for isolated rotating stars *** DEPRECATED : use class Star_rot instead ***.
Definition etoile.h:1496
Basic integer array class.
Definition itbl.h:122
Base class for coordinate mappings.
Definition map.h:670
Basic array class.
Definition tbl.h:161
Tensor handling *** DEPRECATED : use class Tensor instead ***.
Definition tenseur.h:301
Lorene prototypes.
Definition app_hor.h:64