LORENE
etoile.h
1/*
2 * Definition of Lorene classes Etoile
3 * Etoile_bin
4 * Etoile_rot
5 *
6 */
7
8/*
9 * Copyright (c) 2000-2001 Eric Gourgoulhon
10 * Copyright (c) 2000-2001 Keisuke Taniguchi
11 *
12 * This file is part of LORENE.
13 *
14 * LORENE is free software; you can redistribute it and/or modify
15 * it under the terms of the GNU General Public License as published by
16 * the Free Software Foundation; either version 2 of the License, or
17 * (at your option) any later version.
18 *
19 * LORENE is distributed in the hope that it will be useful,
20 * but WITHOUT ANY WARRANTY; without even the implied warranty of
21 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
22 * GNU General Public License for more details.
23 *
24 * You should have received a copy of the GNU General Public License
25 * along with LORENE; if not, write to the Free Software
26 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
27 *
28 */
29
30
31#ifndef __ETOILE_H_
32#define __ETOILE_H_
33
34/*
35 * $Id: etoile.h,v 1.35 2015/06/12 12:38:24 j_novak Exp $
36 * $Log: etoile.h,v $
37 * Revision 1.35 2015/06/12 12:38:24 j_novak
38 * Implementation of the corrected formula for the quadrupole momentum.
39 *
40 * Revision 1.34 2015/06/11 13:50:19 j_novak
41 * Minor corrections
42 *
43 * Revision 1.33 2015/06/10 14:36:39 a_sourie
44 * Corrected the formula for the computation of the quadrupole momentum.
45 *
46 * Revision 1.32 2014/10/13 08:52:34 j_novak
47 * Lorene classes and functions now belong to the namespace Lorene.
48 *
49 * Revision 1.31 2011/10/06 14:55:36 j_novak
50 * equation_of_state() is now virtual to be able to call to the magnetized
51 * Eos_mag.
52 *
53 * Revision 1.30 2010/02/02 21:05:49 e_gourgoulhon
54 * Etoile_bin:equilibrium : suppressed the argument method_khi added by
55 * mistake during previous commit.
56 *
57 * Revision 1.29 2010/02/02 13:34:12 e_gourgoulhon
58 * Marked DEPRECATED (in the documentation).
59 *
60 * Revision 1.28 2008/11/14 13:51:08 e_gourgoulhon
61 * Added the parameter ent_limit to Etoile::equilibrium_spher and
62 * Etoile_bin::equilibrium.
63 *
64 * Revision 1.27 2005/10/05 15:14:47 j_novak
65 * Added a Param* as parameter of Etoile_rot::equilibrium
66 *
67 * Revision 1.26 2005/08/29 15:10:12 p_grandclement
68 * Addition of things needed :
69 * 1) For BBH with different masses
70 * 2) Provisory files for the mixted binaries (Bh and NS) : THIS IS NOT
71 * WORKING YET !!!
72 *
73 * Revision 1.25 2005/01/18 22:35:51 k_taniguchi
74 * Delete a pointer for Etoile::ray_eq(int kk).
75 *
76 * Revision 1.24 2005/01/18 20:34:14 k_taniguchi
77 * Addition of Etoile::ray_eq(int kk).
78 *
79 * Revision 1.23 2005/01/17 20:39:32 k_taniguchi
80 * Addition of Etoile::ray_eq_3pis2().
81 *
82 * Revision 1.22 2004/11/30 20:40:06 k_taniguchi
83 * Addition of the method for calculating a spherical star with falloff
84 * condition at the outer boundary.
85 *
86 * Revision 1.21 2004/05/10 10:18:33 f_limousin
87 * Change to avoid warning in the compilation of Lorene
88 *
89 * Revision 1.20 2004/05/07 12:37:12 f_limousin
90 * Add new member ssjm1_psi
91 *
92 * Revision 1.19 2004/04/08 16:42:31 f_limousin
93 * Add a function velocity_potential with argument ssjm1_psi for the
94 * case of strange stars.
95 *
96 * Revision 1.18 2004/03/22 13:12:41 j_novak
97 * Modification of comments to use doxygen instead of doc++
98 *
99 * Revision 1.17 2003/10/24 12:24:41 k_taniguchi
100 * Suppress the methods of update metric for NS-BH
101 *
102 * Revision 1.16 2003/10/21 11:44:43 k_taniguchi
103 * Delete various things for the Bin_ns_bh project.
104 * They are moved to et_bin_nsbh.h.
105 *
106 * Revision 1.15 2003/10/20 14:50:04 k_taniguchi
107 * Addition of various things for the Bin_ns_bh project
108 * which are related with the part of the neutron star.
109 *
110 * Revision 1.14 2003/10/20 13:11:03 k_taniguchi
111 * Back to version 1.12
112 *
113 * Revision 1.13 2003/10/20 12:15:55 k_taniguchi
114 * Addition of various things for the Bin_ns_bh project
115 * which are related with the part of the neutron star.
116 *
117 * Revision 1.12 2003/06/20 14:13:16 f_limousin
118 * Change to virtual the functions equilibrium_spher() and kinematics().
119 *
120 * Revision 1.11 2003/02/13 16:40:24 p_grandclement
121 * Addition of various things for the Bin_ns_bh project, non of them being
122 * completely tested
123 *
124 * Revision 1.10 2003/02/04 16:20:35 f_limousin
125 * Change to virtual the routine extrinsic_curvature
126 *
127 * Revision 1.9 2003/01/31 16:57:12 p_grandclement
128 * addition of the member Cmp decouple used to compute the K_ij auto, once
129 * the K_ij total is known
130 *
131 * Revision 1.8 2002/12/19 14:48:00 e_gourgoulhon
132 *
133 * Class Etoile_bin: added the new functions:
134 * void update_metric(const Bhole& comp)
135 * void update_metric_der_comp(const Bhole& comp)
136 * to treat the case where the companion is a black hole
137 *
138 * Revision 1.7 2002/12/17 21:17:08 e_gourgoulhon
139 * Class Etoile_bin: suppression of the member p_companion
140 * as well as the associated functions set_companion
141 * and get_companion.
142 *
143 * Revision 1.6 2002/09/13 09:17:33 j_novak
144 * Modif. commentaires
145 *
146 * Revision 1.5 2002/06/17 14:05:16 j_novak
147 * friend functions are now also declared outside the class definition
148 *
149 * Revision 1.4 2002/04/05 09:09:36 j_novak
150 * The inversion of the EOS for 2-fluids polytrope has been modified.
151 * Some errors in the determination of the surface were corrected.
152 *
153 * Revision 1.3 2002/01/11 14:09:34 j_novak
154 * Added newtonian version for 2-fluid stars
155 *
156 * Revision 1.2 2001/12/06 15:11:43 jl_zdunik
157 * Introduction of the new function f_eq() in the class Etoile_rot
158 *
159 * Revision 1.1.1.1 2001/11/20 15:19:27 e_gourgoulhon
160 * LORENE
161 *
162 * Revision 2.61 2001/10/25 08:32:57 eric
163 * Etoile_rot::display_poly passe de protected a public.
164 *
165 * Revision 2.60 2001/10/24 15:35:55 eric
166 * Classe Etoile_rot: ajout de la fonction display_poly.
167 *
168 * Revision 2.59 2001/10/16 14:48:00 eric
169 * La fonction Etoile_rot::get_omega() devient
170 * virtual Etoile_rot::get_omega_c()
171 * (retourne la valeur centrale de Omega).
172 *
173 * Revision 2.58 2001/10/11 09:24:00 eric
174 * La fonction Etoile_rot::equilibrium est desormais virtuelle.
175 *
176 * Revision 2.57 2001/08/06 15:39:04 keisuke
177 * Addition of a new argument to Etoile_bin::equilibrium and equil_regular.
178 *
179 * Revision 2.56 2001/06/25 12:52:33 eric
180 * Classe Etoile_bin : ajout du membre p_companion et des fonctions
181 * associees set_companion() et get_companion().
182 *
183 * Revision 2.55 2001/06/13 14:11:55 eric
184 * Modif commentaires (mise en conformite Doc++ 3.4.7)
185 *
186 * Revision 2.54 2001/03/26 09:29:59 jlz
187 * Classe Etoile_rot: new members p_espec_isco and p_lspec_isco.
188 *
189 * Revision 2.53 2001/02/08 15:37:09 eric
190 * *** empty log message ***
191 *
192 * Revision 2.52 2001/02/08 15:12:42 eric
193 * Ajout de la fonction Etoile_rot::f_eccentric.
194 *
195 * Revision 2.51 2000/11/23 15:43:24 eric
196 * Ajout de l'argument ent_limit a Etoile_rot::equilibrium.
197 *
198 * Revision 2.50 2000/11/19 18:51:13 eric
199 * Etoile_rot: ajout de la fonction (static) lambda_grv2
200 *
201 * Revision 2.49 2000/11/18 23:17:32 eric
202 * Ajout de l'argument ost a Etoile_rot::r_isco.
203 *
204 * Revision 2.48 2000/11/18 21:08:33 eric
205 * Classe Etoile_rot: ajout des fonctions r_isco() et f_isco()
206 * ainsi que des membres associes p_r_isco et p_f_isco.
207 *
208 * Revision 2.47 2000/11/10 15:15:38 eric
209 * Modif arguments Etoile_rot::equilibrium.
210 *
211 * Revision 2.46 2000/10/20 13:10:23 eric
212 * Etoile_rot::equilibrium: ajout de l'argument nzadapt.
213 *
214 * Revision 2.45 2000/10/17 15:59:14 eric
215 * Modif commentaires Etoile_rot::equilibrium
216 *
217 * Revision 2.44 2000/10/12 15:22:29 eric
218 * Modif commentaires Etoile_rot.
219 *
220 * Revision 2.43 2000/10/12 10:20:12 eric
221 * Ajout de la fonction Etoile_rot::fait_nphi().
222 *
223 * Revision 2.42 2000/09/22 15:50:13 keisuke
224 * d_logn_auto_div devient desormais un membre de la classe Etoile
225 * et non plus de la classe derivee Etoile_bin.
226 *
227 * Revision 2.41 2000/09/18 16:14:37 eric
228 * Classe Etoile_rot: ajout du membre tkij et de la fonction
229 * extrinsic curvature().
230 *
231 * Revision 2.40 2000/09/07 14:31:09 keisuke
232 * Ajout des membres logn_auto_regu et get_logn_auto_regu() a la classe Etoile.
233 * Ajout des membres d_logn_auto_regu et get_d_logn_auto_regu()
234 * a la classe Etoile_bin.
235 *
236 * Revision 2.39 2000/09/07 10:25:50 keisuke
237 * Ajout du membre get_logn_auto_div() a la classe Etoile.
238 * Ajout du membre get_d_logn_auto_div() a la classe Etoile_bin.
239 *
240 * Revision 2.38 2000/08/31 11:25:24 eric
241 * Classe Etoile_rot: ajout des membres tnphi et ak_car.
242 *
243 * Revision 2.37 2000/08/29 11:37:06 eric
244 * Ajout des membres k_div et logn_auto_div a la classe Etoile.
245 * Ajout du membre d_logn_auto_div a la classe Etoile_bin.
246 *
247 * Revision 2.36 2000/08/25 10:25:57 keisuke
248 * Ajout de Etoile_bin::equil_regular.
249 *
250 * Revision 2.35 2000/08/18 14:01:07 eric
251 * Modif commentaires.
252 *
253 * Revision 2.34 2000/08/17 12:38:30 eric
254 * Modifs classe Etoile_rot : ajout des membres nuf, nuq, ssjm1_nuf et ssjm1_nuq
255 * Modif arguments Etoile_rot::equilibrium.
256 * .\
257 *
258 * Revision 2.33 2000/08/07 12:11:13 keisuke
259 * Ajout de Etoile::equil_spher_regular.
260 *
261 * Revision 2.32 2000/07/21 12:02:55 eric
262 * Suppression de Etoile_rot::relaxation.
263 *
264 * Revision 2.31 2000/07/20 15:32:28 eric
265 * *** empty log message ***
266 *
267 * Revision 2.30 2000/07/06 09:39:12 eric
268 * Ajout du membre p_xa_barycenter a Etoile_bin, ainsi que de la
269 * fonction associee xa_barycenter().
270 *
271 * Revision 2.29 2000/05/25 13:47:38 eric
272 * Modif Etoile_bin::equilibrium: ajout de l'argument thres_adapt
273 *
274 * Revision 2.28 2000/03/22 16:41:45 eric
275 * Ajout de la sortie de nouvelles erreurs dans Etoile_bin::equilibrium.
276 *
277 * Revision 2.27 2000/03/22 12:54:44 eric
278 * Nouveau prototype d'Etoile_bin::velocity_potential : l'erreur est
279 * retournee en double.
280 * Nouveau prototype d'Etoile_bin::equilibrium : diff_ent est remplace
281 * par le Tbl diff.
282 *
283 * Revision 2.26 2000/03/15 11:04:15 eric
284 * Ajout des fonctions Etoile_bin::set_w_shift() et Etoile_bin::set_khi_shift()
285 * Amelioration des commentaires.
286 *
287 * Revision 2.25 2000/03/08 12:12:49 eric
288 * Ajout de la fonction Etoile_bin::is_irrotational().
289 *
290 * Revision 2.24 2000/03/07 14:48:01 eric
291 * Ajout de la fonction Etoile_bin::extrinsic_curvature().
292 *
293 * Revision 2.23 2000/02/21 13:57:57 eric
294 * classe Etoile_bin: suppression du membre psi: remplacement par psi0.
295 *
296 * Revision 2.22 2000/02/17 16:51:22 eric
297 * Ajout de l'argument diff_ent dans Etoile_bin::equilibrium.
298 *
299 * Revision 2.21 2000/02/17 15:29:38 eric
300 * Ajout de la fonction Etoile_bin::relaxation.
301 *
302 * Revision 2.20 2000/02/17 13:54:21 eric
303 * Ajout de la fonction Etoile_bin::velocity_potential.
304 *
305 * Revision 2.19 2000/02/16 15:05:14 eric
306 * Ajout des membres w_shift et khi_shift.
307 * (sauves dans les fichiers a la place de shift_auto).
308 * Ajout de la fonction Etoile_bin::fait_shift_auto.
309 *
310 * Revision 2.18 2000/02/16 13:47:02 eric
311 * Classe Etoile_bin: ajout des membres ssjm1_khi et ssjm1_wshift.
312 *
313 * Revision 2.17 2000/02/16 11:54:13 eric
314 * Classe Etoile_bin : ajout des membres ssjm1_logn et ssjm1_beta.
315 *
316 * Revision 2.16 2000/02/15 15:40:07 eric
317 * Ajout de Etoile_bin::equilibrium.
318 *
319 * Revision 2.15 2000/02/12 18:40:15 eric
320 * Modif commentaires.
321 *
322 * Revision 2.14 2000/02/12 14:44:26 eric
323 * Ajout des fonctions Etoile_bin::set_logn_comp et set_pot_centri.
324 *
325 * Revision 2.13 2000/02/10 20:22:25 eric
326 * Modif commentaires.
327 *
328 * Revision 2.12 2000/02/10 16:11:24 eric
329 * Classe Etoile_bin : ajout des accesseurs get_psi, etc...
330 * ajout de la fonction fait_d_psi
331 *
332 * Revision 2.11 2000/02/08 19:28:29 eric
333 * La fonction Etoile_bin::scal_prod est rebaptisee Etoile_bin::sprod
334 *
335 * Revision 2.10 2000/02/04 17:15:15 eric
336 * Classe Etoile_bin: ajout du membre ref_triad.
337 *
338 * Revision 2.9 2000/02/04 16:36:48 eric
339 * Ajout des fonctions update_metric* et kinematics.
340 *
341 * Revision 2.8 2000/02/02 10:12:37 eric
342 * Ajout des fonctions de lecture/ecriture mp, nzet, eos, etc...
343 *
344 * Revision 2.7 2000/02/01 15:59:43 eric
345 * Ajout de la fonction Etoile_bin::scal_prod.
346 *
347 * Revision 2.6 2000/01/31 15:56:45 eric
348 * Introduction de la classe derivee Etoile_bin.
349 *
350 * Revision 2.5 2000/01/28 17:17:45 eric
351 * Ajout des fonctions de calcul des quantites globales.
352 *
353 * Revision 2.4 2000/01/27 16:46:59 eric
354 * Ajout des fonctions get_ent(), etc....
355 *
356 * Revision 2.3 2000/01/24 17:19:48 eric
357 * Modif commentaires.
358 *
359 * Revision 2.2 2000/01/24 17:13:04 eric
360 * Le mapping mp n'est plus constant.
361 * Ajout de la fonction equilibrium_spher.
362 *
363 * Revision 2.1 2000/01/24 13:37:19 eric
364 * *** empty log message ***
365 *
366 * Revision 2.0 2000/01/20 17:04:33 eric
367 * *** empty log message ***
368 *
369 *
370 * $Header: /cvsroot/Lorene/C++/Include/etoile.h,v 1.35 2015/06/12 12:38:24 j_novak Exp $
371 *
372 */
373
374// Headers Lorene
375#include "tenseur.h"
376namespace Lorene {
377 class Eos ;
378 class Bhole ;
379
380
381 //---------------------------//
382 // base class Etoile //
383 //---------------------------//
384
424 class Etoile {
425
426 // Data :
427 // -----
428 protected:
429 Map& mp ;
430
432 int nzet ;
433
438
442 double unsurc2 ;
443
449 int k_div ;
450
451 const Eos& eos ;
452
453 // Fluid quantities with respect to the fluid frame
454 // ------------------------------------------------
455
458
462
463 // Fluid quantities with respect to the Eulerian frame
464 // ---------------------------------------------------
466
469
472
475
476 // Metric potentials
477 // -----------------
478
485
492
498
502
507
510
513
516
517 // Derived data :
518 // ------------
519 protected:
521 mutable double* p_ray_eq ;
522
524 mutable double* p_ray_eq_pis2 ;
525
527 mutable double* p_ray_eq_pi ;
528
530 mutable double* p_ray_eq_3pis2 ;
531
533 mutable double* p_ray_pole ;
534
539 mutable Itbl* p_l_surf ;
540
545 mutable Tbl* p_xi_surf ;
546
547 mutable double* p_mass_b ;
548 mutable double* p_mass_g ;
549
550
551 // Constructors - Destructor
552 // -------------------------
553 public:
554
564 Etoile(Map& mp_i, int nzet_i, bool relat, const Eos& eos_i) ;
565
566 Etoile(const Etoile& ) ;
567
575 Etoile(Map& mp_i, const Eos& eos_i, FILE* fich) ;
576
577 virtual ~Etoile() ;
578
579 // Memory management
580 // -----------------
581 protected:
583 virtual void del_deriv() const ;
584
586 virtual void set_der_0x0() const ;
587
591 virtual void del_hydro_euler() ;
592
593
594 // Mutators / assignment
595 // ---------------------
596 public:
598 void operator=(const Etoile&) ;
599
601 Map& set_mp() {return mp; } ;
602
604 void set_enthalpy(const Cmp& ) ;
605
609 virtual void equation_of_state() ;
610
615 virtual void hydro_euler() ;
616
629 virtual void equilibrium_spher(double ent_c, double precis = 1.e-14,
630 const Tbl* ent_limit = 0x0 ) ;
631
642 void equil_spher_regular(double ent_c, double precis = 1.e-14) ;
643
652 virtual void equil_spher_falloff(double ent_c,
653 double precis = 1.e-14) ;
654
655 // Accessors
656 // ---------
657 public:
659 const Map& get_mp() const {return mp; } ;
660
662 int get_nzet() const {return nzet; } ;
663
667 bool is_relativistic() const {return relativistic; } ;
668
670 const Eos& get_eos() const {return eos; } ;
671
673 const Tenseur& get_ent() const {return ent;} ;
674
676 const Tenseur& get_nbar() const {return nbar;} ;
677
679 const Tenseur& get_ener() const {return ener;} ;
680
682 const Tenseur& get_press() const {return press;} ;
683
685 const Tenseur& get_ener_euler() const {return ener_euler;} ;
686
688 const Tenseur& get_s_euler() const {return s_euler;} ;
689
691 const Tenseur& get_gam_euler() const {return gam_euler;} ;
692
694 const Tenseur& get_u_euler() const {return u_euler;} ;
695
701 const Tenseur& get_logn_auto() const {return logn_auto;} ;
702
708 const Tenseur& get_logn_auto_regu() const {return logn_auto_regu;} ;
709
715 const Tenseur& get_logn_auto_div() const {return logn_auto_div;} ;
716
719 const Tenseur& get_d_logn_auto_div() const {return d_logn_auto_div;} ;
720
724 const Tenseur& get_beta_auto() const {return beta_auto;} ;
725
727 const Tenseur& get_nnn() const {return nnn;} ;
728
730 const Tenseur& get_shift() const {return shift;} ;
731
733 const Tenseur& get_a_car() const {return a_car;} ;
734
735 // Outputs
736 // -------
737 public:
738 virtual void sauve(FILE* ) const ;
739
741 friend ostream& operator<<(ostream& , const Etoile& ) ;
742
743 protected:
745 virtual ostream& operator>>(ostream& ) const ;
746
747 // Global quantities
748 // -----------------
749 public:
751 double ray_eq() const ;
752
754 double ray_eq_pis2() const ;
755
757 double ray_eq_pi() const ;
758
760 double ray_eq_3pis2() const ;
761
763 double ray_pole() const ;
764
766 double ray_eq(int kk) const ;
767
775 virtual const Itbl& l_surf() const ;
776
784 const Tbl& xi_surf() const ;
785
787 virtual double mass_b() const ;
788
790 virtual double mass_g() const ;
791
792 };
793 ostream& operator<<(ostream& , const Etoile& ) ;
794
795
796 //---------------------------//
797 // class Etoile_bin //
798 //---------------------------//
799
814 class Etoile_bin : public Etoile {
815
816 // Data :
817 // -----
818 protected:
823
829
834
839
845
850
855
860
865
870
875
880
885
890
896
909
919
926
933
939
945
951
954
960
966
974
984
990
1001
1002 // Derived data :
1003 // ------------
1004 protected:
1006 mutable double* p_xa_barycenter ;
1007
1008
1009 // Constructors - Destructor
1010 // -------------------------
1011 public:
1027 Etoile_bin(Map& mp_i, int nzet_i, bool relat, const Eos& eos_i,
1028 bool irrot, const Base_vect& ref_triad_i) ;
1029
1030
1031 Etoile_bin(const Etoile_bin& ) ;
1032
1045 Etoile_bin(Map& mp_i, const Eos& eos_i, const Base_vect& ref_triad_i,
1046 FILE* fich) ;
1047
1048 virtual ~Etoile_bin() ;
1049
1050
1051 // Memory management
1052 // -----------------
1053 protected:
1055 virtual void del_deriv() const ;
1056
1058 virtual void set_der_0x0() const ;
1059
1063 virtual void del_hydro_euler() ;
1064
1065
1066 // Mutators / assignment
1067 // ---------------------
1068 public:
1070 void operator=(const Etoile_bin& ) ;
1071
1076
1079
1081 Tenseur& set_w_shift() ;
1082
1085
1086 // Accessors
1087 // ---------
1088 public:
1092 bool is_irrotational() const {return irrotational; } ;
1093
1095 const Tenseur& get_psi0() const {return psi0;} ;
1096
1100 const Tenseur& get_d_psi() const {return d_psi;} ;
1101
1106 const Tenseur& get_wit_w() const {return wit_w;} ;
1107
1111 const Tenseur& get_loggam() const {return loggam;} ;
1112
1116 const Tenseur& get_logn_comp() const {return logn_comp;} ;
1117
1121 const Tenseur& get_d_logn_auto() const {return d_logn_auto;} ;
1122
1127
1131 const Tenseur& get_d_logn_comp() const {return d_logn_comp;} ;
1132
1136 const Tenseur& get_beta_comp() const {return beta_comp;} ;
1137
1141 const Tenseur& get_d_beta_auto() const {return d_beta_auto;} ;
1142
1146 const Tenseur& get_d_beta_comp() const {return d_beta_comp;} ;
1147
1152 const Tenseur& get_shift_auto() const {return shift_auto;} ;
1153
1158 const Tenseur& get_shift_comp() const {return shift_comp;} ;
1159
1172 const Tenseur& get_w_shift() const {return w_shift;} ;
1173
1186 const Tenseur& get_khi_shift() const {return khi_shift;} ;
1187
1192 const Tenseur_sym& get_tkij_auto() const {return tkij_auto;} ;
1193
1198 const Tenseur_sym& get_tkij_comp() const {return tkij_comp;} ;
1199
1204 const Tenseur& get_akcar_auto() const {return akcar_auto;} ;
1205
1210 const Tenseur& get_akcar_comp() const {return akcar_comp;} ;
1211
1216 const Tenseur& get_bsn() const {return bsn;} ;
1217
1219 const Tenseur& get_pot_centri() const {return pot_centri;} ;
1224 const Cmp get_decouple() const {return decouple ;}
1225
1226 // Outputs
1227 // -------
1228 public:
1229 virtual void sauve(FILE* ) const ;
1230
1231 protected:
1233 virtual ostream& operator>>(ostream& ) const ;
1234
1235 // Global quantities
1236 // -----------------
1237 public:
1239 virtual double mass_b() const ;
1240
1242 virtual double mass_g() const ;
1243
1252 virtual double xa_barycenter() const ;
1253
1254
1255 // Computational routines
1256 // ----------------------
1257 public:
1265 virtual Tenseur sprod(const Tenseur& t1, const Tenseur& t2) const ;
1266
1279 virtual void hydro_euler() ;
1280
1298 void update_metric(const Etoile_bin& comp) ;
1299
1307 void update_metric(const Etoile_bin& comp, const Etoile_bin& star_prev,
1308 double relax) ;
1309
1324 void update_metric_der_comp(const Etoile_bin& comp) ;
1325
1336 virtual void kinematics(double omega, double x_axe) ;
1337
1341 void fait_d_psi() ;
1342
1351 void fait_shift_auto() ;
1352
1356 virtual void extrinsic_curvature() ;
1357
1398 void equilibrium(double ent_c,
1399 int mermax, int mermax_poisson,
1400 double relax_poisson, int mermax_potvit,
1401 double relax_potvit, double thres_adapt,
1402 const Tbl& fact, Tbl& diff, const Tbl* ent_limit = 0x0) ;
1403
1443 void equil_regular(double ent_c, int mermax, int mermax_poisson,
1444 double relax_poisson, int mermax_potvit,
1445 double relax_potvit, double thres_adapt,
1446 const Tbl& fact, Tbl& diff) ;
1447
1461 double velocity_potential(int mermax, double precis, double relax) ;
1462
1473 void relaxation(const Etoile_bin& star_prev, double relax_ent,
1474 double relax_met, int mer, int fmer_met) ;
1475
1476 friend class Bin_ns_bh ;
1477 };
1478
1479
1480
1481 //---------------------------//
1482 // class Etoile_rot //
1483 //---------------------------//
1484
1496class Etoile_rot : public Etoile {
1497
1498 // Data :
1499 // -----
1500 protected:
1501 double omega ;
1502
1505
1508
1511
1516
1519
1522
1527
1532
1535
1538
1551
1561
1568
1587
1593
1599
1604
1609
1617
1626
1627 // Derived data :
1628 // ------------
1629 protected:
1630
1631 mutable double* p_angu_mom ;
1632 mutable double* p_tsw ;
1633 mutable double* p_grv2 ;
1634 mutable double* p_grv3 ;
1635 mutable double* p_r_circ ;
1636 mutable double* p_area ;
1637 mutable double* p_aplat ;
1638 mutable double* p_z_eqf ;
1639 mutable double* p_z_eqb ;
1640 mutable double* p_z_pole ;
1641 mutable double* p_mom_quad ;
1642 mutable double* p_mom_quad_old ;
1643 mutable double* p_mom_quad_Bo ;
1644 mutable double* p_r_isco ;
1645 mutable double* p_f_isco ;
1647 mutable double* p_espec_isco ;
1649 mutable double* p_lspec_isco ;
1650 mutable double* p_f_eq ;
1651
1652
1653
1654 // Constructors - Destructor
1655 // -------------------------
1656 public:
1665 Etoile_rot(Map& mp_i, int nzet_i, bool relat, const Eos& eos_i) ;
1666
1667
1668 Etoile_rot(const Etoile_rot& ) ;
1669
1677 Etoile_rot(Map& mp_i, const Eos& eos_i, FILE* fich) ;
1678
1679 virtual ~Etoile_rot() ;
1680
1681
1682 // Memory management
1683 // -----------------
1684 protected:
1686 virtual void del_deriv() const ;
1687
1689 virtual void set_der_0x0() const ;
1690
1694 virtual void del_hydro_euler() ;
1695
1696
1697 // Mutators / assignment
1698 // ---------------------
1699 public:
1701 void operator=(const Etoile_rot& ) ;
1702
1703 // Accessors
1704 // ---------
1705 public:
1709 virtual double get_omega_c() const ;
1710
1712 const Tenseur& get_bbb() const {return bbb;} ;
1713
1715 const Tenseur& get_b_car() const {return b_car;} ;
1716
1718 const Tenseur& get_nphi() const {return nphi;} ;
1719
1723 const Tenseur& get_tnphi() const {return tnphi;} ;
1724
1726 const Tenseur& get_uuu() const {return uuu;} ;
1727
1729 const Tenseur& get_logn() const {return logn;} ;
1730
1734 const Tenseur& get_nuf() const {return nuf;} ;
1735
1739 const Tenseur& get_nuq() const {return nuq;} ;
1740
1742 const Tenseur& get_dzeta() const {return dzeta;} ;
1743
1745 const Tenseur& get_tggg() const {return tggg;} ;
1746
1759 const Tenseur& get_w_shift() const {return w_shift;} ;
1760
1773 const Tenseur& get_khi_shift() const {return khi_shift;} ;
1774
1780 const Tenseur_sym& get_tkij() const {return tkij;} ;
1781
1799 const Tenseur& get_ak_car() const {return ak_car;} ;
1800
1801 // Outputs
1802 // -------
1803 public:
1804 virtual void sauve(FILE* ) const ;
1805
1807 virtual void display_poly(ostream& ) const ;
1808
1809 protected:
1811 virtual ostream& operator>>(ostream& ) const ;
1812
1814 virtual void partial_display(ostream& ) const ;
1815
1816 // Global quantities
1817 // -----------------
1818 public:
1819
1827 virtual const Itbl& l_surf() const ;
1828
1829 virtual double mass_b() const ;
1830 virtual double mass_g() const ;
1831 virtual double angu_mom() const ;
1832 virtual double tsw() const ;
1833
1837 virtual double grv2() const ;
1838
1850 virtual double grv3(ostream* ost = 0x0) const ;
1851
1852 virtual double r_circ() const ;
1853 virtual double area() const ;
1854 virtual double mean_radius() const ;
1855 virtual double aplat() const ;
1856 virtual double z_eqf() const ;
1857 virtual double z_eqb() const ;
1858 virtual double z_pole() const ;
1859
1874 virtual double mom_quad() const ;
1875
1883 virtual double mom_quad_old() const ;
1884
1891 virtual double mom_quad_Bo() const ;
1892
1899 virtual double r_isco(ostream* ost = 0x0) const ;
1900
1902 virtual double f_isco() const ;
1903
1905 virtual double espec_isco() const ;
1906
1908 virtual double lspec_isco() const ;
1909
1910
1921 virtual double f_eccentric(double ecc, double periast,
1922 ostream* ost = 0x0) const ;
1923
1925 virtual double f_eq() const ;
1926
1927
1928 // Computational routines
1929 // ----------------------
1930 public:
1941 virtual void hydro_euler() ;
1942
1952 void update_metric() ;
1953
1962 void fait_shift() ;
1963
1967 void fait_nphi() ;
1968
1972 void extrinsic_curvature() ;
1973
2003 static double lambda_grv2(const Cmp& sou_m, const Cmp& sou_q) ;
2004
2084 virtual void equilibrium(double ent_c, double omega0, double fact_omega,
2085 int nzadapt, const Tbl& ent_limit,
2086 const Itbl& icontrol, const Tbl& control,
2087 double mbar_wanted, double aexp_mass,
2088 Tbl& diff, Param* = 0x0) ;
2089
2090
2091};
2092
2093
2094
2095
2096}
2097#endif
Vectorial bases (triads) with respect to which the tensorial components are defined.
Definition base_vect.h:105
Neutron star - black hole binary system.
Definition bin_ns_bh.h:117
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition cmp.h:446
Equation of state base class.
Definition eos.h:190
Class for stars in binary system.
Definition etoile.h:814
Tenseur shift_comp
Part of the shift vector generated principaly by the companion star.
Definition etoile.h:895
void update_metric(const Etoile_bin &comp)
Computes metric coefficients from known potentials, when the companion is another star.
const Tenseur & get_psi0() const
Returns the non-translational part of the velocity potential.
Definition etoile.h:1095
Tenseur psi0
Scalar potential of the non-translational part of fluid 4-velocity (in the irrotational case)
Definition etoile.h:833
const Tenseur & get_d_logn_auto_regu() const
Returns the gradient of logn_auto_regu (Cartesian components with respect to ref_triad )
Definition etoile.h:1126
Tenseur d_beta_auto
Gradient of beta_auto (Cartesian components with respect to ref_triad )
Definition etoile.h:879
Tenseur d_logn_auto
Gradient of logn_auto (Cartesian components with respect to ref_triad )
Definition etoile.h:859
virtual double mass_b() const
Baryon mass.
double * p_xa_barycenter
Absolute coordinate X of the barycenter of the baryon density.
Definition etoile.h:1006
const Tenseur & get_d_logn_auto() const
Returns the gradient of logn_auto (Cartesian components with respect to ref_triad )
Definition etoile.h:1121
virtual void sauve(FILE *) const
Save in a file.
Definition etoile_bin.C:559
Tenseur bsn
3-vector shift, divided by N , of the rotating coordinates, .
Definition etoile.h:950
Tenseur akcar_comp
Part of the scalar generated by shift_auto and shift_comp , i.e.
Definition etoile.h:944
virtual void del_hydro_euler()
Sets to ETATNONDEF (undefined state) the hydrodynamical quantities relative to the Eulerian observer.
Definition etoile_bin.C:467
Cmp ssjm1_khi
Effective source at the previous step for the resolution of the Poisson equation for the scalar by m...
Definition etoile.h:973
const Tenseur & get_wit_w() const
Returns the spatial projection of the fluid 3-velocity with respect to the co-orbiting observer.
Definition etoile.h:1106
Tenseur d_psi
Gradient of (in the irrotational case) (Cartesian components with respect to ref_triad )
Definition etoile.h:838
double velocity_potential(int mermax, double precis, double relax)
Computes the non-translational part of the velocity scalar potential by solving the continuity equat...
void equilibrium(double ent_c, int mermax, int mermax_poisson, double relax_poisson, int mermax_potvit, double relax_potvit, double thres_adapt, const Tbl &fact, Tbl &diff, const Tbl *ent_limit=0x0)
Computes an equilibrium configuration.
void relaxation(const Etoile_bin &star_prev, double relax_ent, double relax_met, int mer, int fmer_met)
Performs a relaxation on ent , logn_auto , beta_auto and shift_auto .
Definition etoile_bin.C:863
const Base_vect & ref_triad
Reference triad ("absolute frame"), with respect to which the components of all the member Tenseur 's...
Definition etoile.h:828
void fait_shift_auto()
Computes shift_auto from w_shift and khi_shift according to Shibata's prescription [Prog.
Definition etoile_bin.C:826
void update_metric_der_comp(const Etoile_bin &comp)
Computes the derivative of metric functions related to the companion star.
const Tenseur & get_logn_comp() const
Returns the part of the lapse logarithm (gravitational potential at the Newtonian limit) generated pr...
Definition etoile.h:1116
bool irrotational
true for an irrotational star, false for a corotating one
Definition etoile.h:822
const Tenseur & get_shift_auto() const
Returns the part of the shift vector generated principaly by the star.
Definition etoile.h:1152
Tenseur w_shift
Vector used in the decomposition of shift_auto , following Shibata's prescription [Prog.
Definition etoile.h:908
Tenseur logn_comp
Part of the lapse logarithm (gravitational potential at the Newtonian limit) generated principaly by ...
Definition etoile.h:854
const Cmp get_decouple() const
Returns the function used to construct tkij_auto from tkij_tot .
Definition etoile.h:1224
Tenseur & set_w_shift()
Read/write of w_shift.
Definition etoile_bin.C:538
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition etoile_bin.C:459
const Tenseur_sym & get_tkij_auto() const
Returns the part of the extrinsic curvature tensor generated by shift_auto .
Definition etoile.h:1192
Tenseur d_logn_comp
Gradient of logn_comp (Cartesian components with respect to ref_triad )
Definition etoile.h:869
Tenseur beta_comp
Part of the logarithm of AN generated principaly by the companion star.
Definition etoile.h:874
virtual void hydro_euler()
Computes the hydrodynamical quantities relative to the Eulerian observer from those in the fluid fram...
const Tenseur & get_w_shift() const
Returns the vector used in the decomposition of shift_auto , following Shibata's prescription [Prog.
Definition etoile.h:1172
Tenseur ssjm1_wshift
Effective source at the previous step for the resolution of the vector Poisson equation for by means...
Definition etoile.h:983
Cmp ssjm1_logn
Effective source at the previous step for the resolution of the Poisson equation for logn_auto by mea...
Definition etoile.h:959
virtual void del_deriv() const
Deletes all the derived quantities.
Definition etoile_bin.C:447
const Tenseur & get_shift_comp() const
Returns the part of the shift vector generated principaly by the companion star.
Definition etoile.h:1158
Cmp ssjm1_psi
Effective source at the previous step for the resolution of the Poisson equation for the scalar by m...
Definition etoile.h:989
const Tenseur & get_loggam() const
Returns the logarithm of the Lorentz factor between the fluid and the co-orbiting observer.
Definition etoile.h:1111
Tenseur pot_centri
Centrifugal potential.
Definition etoile.h:953
void equil_regular(double ent_c, int mermax, int mermax_poisson, double relax_poisson, int mermax_potvit, double relax_potvit, double thres_adapt, const Tbl &fact, Tbl &diff)
Computes an equilibrium configuration by regularizing the diverging source.
Tenseur loggam
Logarithm of the Lorentz factor between the fluid and the co-orbiting observer.
Definition etoile.h:849
Tenseur shift_auto
Part of the shift vector generated principaly by the star.
Definition etoile.h:889
virtual double xa_barycenter() const
Absolute coordinate X of the barycenter of the baryon density, defined according to the formula.
const Tenseur & get_d_beta_auto() const
Returns the gradient of beta_auto (Cartesian components with respect to ref_triad )
Definition etoile.h:1141
bool is_irrotational() const
Returns true for an irrotational motion, false for a corotating one.
Definition etoile.h:1092
const Tenseur & get_akcar_comp() const
Returns the part of the scalar generated by shift_auto and shift_comp , i.e.
Definition etoile.h:1210
Tenseur_sym tkij_auto
Part of the extrinsic curvature tensor generated by shift_auto .
Definition etoile.h:925
Tenseur akcar_auto
Part of the scalar generated by shift_auto , i.e.
Definition etoile.h:938
Cmp ssjm1_beta
Effective source at the previous step for the resolution of the Poisson equation for beta_auto by mea...
Definition etoile.h:965
Cmp decouple
Function used to construct the part of generated by the star from the total .
Definition etoile.h:1000
Tenseur wit_w
Spatial projection of the fluid 3-velocity with respect to the co-orbiting observer.
Definition etoile.h:844
const Tenseur & get_akcar_auto() const
Returns the part of the scalar generated by shift_auto , i.e.
Definition etoile.h:1204
virtual double mass_g() const
Gravitational mass.
Tenseur & set_logn_comp()
Read/write the part of the lapse logarithm (gravitational potential at the Newtonian limit) generated...
Definition etoile_bin.C:524
Tenseur_sym tkij_comp
Part of the extrinsic curvature tensor generated by shift_comp .
Definition etoile.h:932
Tenseur & set_khi_shift()
Read/write of khi_shift.
Definition etoile_bin.C:545
virtual Tenseur sprod(const Tenseur &t1, const Tenseur &t2) const
Performs the scalar product of two tensors by contracting the last index of t1 with the first index o...
Definition etoile_bin.C:748
const Tenseur & get_pot_centri() const
Returns the centrifugal potential.
Definition etoile.h:1219
virtual void extrinsic_curvature()
Computes tkij_auto and akcar_auto from shift_auto , nnn and a_car .
Tenseur d_beta_comp
Gradient of beta_comp (Cartesian components with respect to ref_triad )
Definition etoile.h:884
Tenseur d_logn_auto_regu
Gradient of logn_auto_regu (Cartesian components with respect to ref_triad )
Definition etoile.h:864
virtual ~Etoile_bin()
Destructor.
Definition etoile_bin.C:437
const Tenseur & get_bsn() const
Returns the shift vector, divided by N , of the rotating coordinates, .
Definition etoile.h:1216
Tenseur & set_pot_centri()
Read/write the centrifugal potential.
Definition etoile_bin.C:531
void fait_d_psi()
Computes the gradient of the total velocity potential .
Definition etoile_bin.C:764
const Tenseur & get_d_logn_comp() const
Returns the gradient of logn_comp (Cartesian components with respect to ref_triad )
Definition etoile.h:1131
const Tenseur & get_d_beta_comp() const
Returns the gradient of beta_comp (Cartesian components with respect to ref_triad )
Definition etoile.h:1146
const Tenseur & get_beta_comp() const
Returns the part of the logarithm of AN generated principaly by the companion star.
Definition etoile.h:1136
virtual void kinematics(double omega, double x_axe)
Computes the quantities bsn and pot_centri .
const Tenseur_sym & get_tkij_comp() const
Returns the part of the extrinsic curvature tensor generated by shift_comp .
Definition etoile.h:1198
void operator=(const Etoile_bin &)
Assignment to another Etoile_bin.
Definition etoile_bin.C:482
const Tenseur & get_khi_shift() const
Returns the scalar used in the decomposition of shift_auto following Shibata's prescription [Prog.
Definition etoile.h:1186
const Tenseur & get_d_psi() const
Returns the gradient of the velocity potential (Cartesian components with respect to ref_triad )
Definition etoile.h:1100
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition etoile_bin.C:582
Tenseur khi_shift
Scalar used in the decomposition of shift_auto , following Shibata's prescription [Prog.
Definition etoile.h:918
Class for isolated rotating stars *** DEPRECATED : use class Star_rot instead ***.
Definition etoile.h:1496
Tenseur ssjm1_wshift
Effective source at the previous step for the resolution of the vector Poisson equation for .
Definition etoile.h:1625
const Tenseur & get_b_car() const
Returns the square of the metric factor B.
Definition etoile.h:1715
Tenseur uuu
Norm of u_euler.
Definition etoile.h:1518
const Tenseur & get_tnphi() const
Returns the component of the shift vector.
Definition etoile.h:1723
const Tenseur & get_w_shift() const
Returns the vector used in the decomposition of shift , following Shibata's prescription [Prog.
Definition etoile.h:1759
virtual void del_deriv() const
Deletes all the derived quantities.
Definition etoile_rot.C:335
double omega
Rotation angular velocity ([f_unit] )
Definition etoile.h:1501
double * p_z_pole
Redshift factor at North pole.
Definition etoile.h:1640
virtual double espec_isco() const
Energy of a particle on the ISCO.
virtual double mom_quad() const
Quadrupole moment.
double * p_z_eqf
Forward redshift factor at equator.
Definition etoile.h:1638
void operator=(const Etoile_rot &)
Assignment to another Etoile_rot.
Definition etoile_rot.C:404
virtual double r_circ() const
Circumferential radius.
virtual double f_eccentric(double ecc, double periast, ostream *ost=0x0) const
Computation of frequency of eccentric orbits.
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition etoile_rot.C:364
Tenseur & logn
Metric potential = logn_auto.
Definition etoile.h:1521
virtual void hydro_euler()
Computes the hydrodynamical quantities relative to the Eulerian observer from those in the fluid fram...
const Tenseur & get_nphi() const
Returns the metric coefficient .
Definition etoile.h:1718
virtual const Itbl & l_surf() const
Description of the stellar surface: returns a 2-D Itbl containing the values of the domain index l on...
const Tenseur & get_khi_shift() const
Returns the scalar used in the decomposition of shift following Shibata's prescription [Prog.
Definition etoile.h:1773
double * p_mom_quad_old
Part of the quadrupole moment.
Definition etoile.h:1642
const Tenseur & get_tggg() const
Returns the Metric potential .
Definition etoile.h:1745
Cmp ssjm1_dzeta
Effective source at the previous step for the resolution of the Poisson equation for dzeta .
Definition etoile.h:1603
Tenseur nuq
Part of the Metric potential = logn generated by the quadratic terms.
Definition etoile.h:1531
const Tenseur_sym & get_tkij() const
Returns the tensor related to the extrinsic curvature tensor by .
Definition etoile.h:1780
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition etoile_rot.C:466
double * p_z_eqb
Backward redshift factor at equator.
Definition etoile.h:1639
virtual double mass_g() const
Gravitational mass.
virtual double f_isco() const
Orbital frequency at the innermost stable circular orbit (ISCO).
double * p_r_circ
Circumferential radius.
Definition etoile.h:1635
Tenseur khi_shift
Scalar used in the decomposition of shift , following Shibata's prescription [Prog.
Definition etoile.h:1560
virtual double mom_quad_Bo() const
Part of the quadrupole moment.
virtual double aplat() const
Flatening r_pole/r_eq.
virtual double grv3(ostream *ost=0x0) const
Error on the virial identity GRV3.
const Tenseur & get_bbb() const
Returns the metric factor B.
Definition etoile.h:1712
virtual double tsw() const
Ratio T/W.
Tenseur tggg
Metric potential .
Definition etoile.h:1537
Tenseur nuf
Part of the Metric potential = logn generated by the matter terms.
Definition etoile.h:1526
Cmp ssjm1_tggg
Effective source at the previous step for the resolution of the Poisson equation for tggg .
Definition etoile.h:1608
Tenseur nphi
Metric coefficient .
Definition etoile.h:1510
virtual double mass_b() const
Baryon mass.
virtual double r_isco(ostream *ost=0x0) const
Circumferential radius of the innermost stable circular orbit (ISCO).
Definition et_rot_isco.C:84
static double lambda_grv2(const Cmp &sou_m, const Cmp &sou_q)
Computes the coefficient which ensures that the GRV2 virial identity is satisfied.
const Tenseur & get_uuu() const
Returns the norm of u_euler.
Definition etoile.h:1726
const Tenseur & get_logn() const
Returns the metric potential = logn_auto.
Definition etoile.h:1729
Tenseur bbb
Metric factor B.
Definition etoile.h:1504
double * p_area
Surface area.
Definition etoile.h:1636
void update_metric()
Computes metric coefficients from known potentials.
void extrinsic_curvature()
Computes tkij and ak_car from shift , nnn and b_car .
Tenseur ak_car
Scalar .
Definition etoile.h:1586
virtual ~Etoile_rot()
Destructor.
Definition etoile_rot.C:325
Tenseur & dzeta
Metric potential = beta_auto.
Definition etoile.h:1534
virtual double mean_radius() const
Mean radius.
double * p_grv3
Error on the virial identity GRV3.
Definition etoile.h:1634
double * p_grv2
Error on the virial identity GRV2.
Definition etoile.h:1633
double * p_aplat
Flatening r_pole/r_eq.
Definition etoile.h:1637
Cmp ssjm1_nuf
Effective source at the previous step for the resolution of the Poisson equation for nuf by means of ...
Definition etoile.h:1592
double * p_mom_quad_Bo
Part of the quadrupole moment.
Definition etoile.h:1643
virtual double lspec_isco() const
Angular momentum of a particle on the ISCO.
const Tenseur & get_nuf() const
Returns the part of the Metric potential = logn generated by the matter terms.
Definition etoile.h:1734
virtual double grv2() const
Error on the virial identity GRV2.
void fait_nphi()
Computes tnphi and nphi from the Cartesian components of the shift, stored in shift .
Definition etoile_rot.C:781
double * p_mom_quad
Quadrupole moment.
Definition etoile.h:1641
double * p_f_eq
Orbital frequency at the equator.
Definition etoile.h:1650
virtual void del_hydro_euler()
Sets to ETATNONDEF (undefined state) the hydrodynamical quantities relative to the Eulerian observer.
Definition etoile_rot.C:389
virtual double area() const
Surface area.
double * p_angu_mom
Angular momentum.
Definition etoile.h:1631
Cmp ssjm1_khi
Effective source at the previous step for the resolution of the Poisson equation for the scalar by m...
Definition etoile.h:1616
virtual void equilibrium(double ent_c, double omega0, double fact_omega, int nzadapt, const Tbl &ent_limit, const Itbl &icontrol, const Tbl &control, double mbar_wanted, double aexp_mass, Tbl &diff, Param *=0x0)
Computes an equilibrium configuration.
virtual double angu_mom() const
Angular momentum.
Tenseur_sym tkij
Tensor related to the extrinsic curvature tensor by .
Definition etoile.h:1567
virtual double z_eqf() const
Forward redshift factor at equator.
const Tenseur & get_dzeta() const
Returns the Metric potential = beta_auto.
Definition etoile.h:1742
virtual void display_poly(ostream &) const
Display in polytropic units.
Definition etoile_rot.C:690
double * p_f_isco
Orbital frequency of the ISCO.
Definition etoile.h:1645
virtual double mom_quad_old() const
Part of the quadrupole moment.
const Tenseur & get_nuq() const
Returns the Part of the Metric potential = logn generated by the quadratic terms.
Definition etoile.h:1739
double * p_lspec_isco
Specific angular momentum of a particle on the ISCO.
Definition etoile.h:1649
double * p_tsw
Ratio T/W.
Definition etoile.h:1632
Tenseur b_car
Square of the metric factor B.
Definition etoile.h:1507
virtual double z_eqb() const
Backward redshift factor at equator.
Tenseur tnphi
Component of the shift vector.
Definition etoile.h:1515
const Tenseur & get_ak_car() const
Returns the scalar .
Definition etoile.h:1799
double * p_espec_isco
Specific energy of a particle on the ISCO.
Definition etoile.h:1647
double * p_r_isco
Circumferential radius of the ISCO.
Definition etoile.h:1644
virtual void sauve(FILE *) const
Save in a file.
Definition etoile_rot.C:441
virtual double f_eq() const
Orbital frequency at the equator.
virtual double get_omega_c() const
Returns the central value of the rotation angular velocity ([f_unit] )
Definition etoile_rot.C:680
Cmp ssjm1_nuq
Effective source at the previous step for the resolution of the Poisson equation for nuq by means of ...
Definition etoile.h:1598
virtual double z_pole() const
Redshift factor at North pole.
Tenseur w_shift
Vector used in the decomposition of shift , following Shibata's prescription [Prog.
Definition etoile.h:1550
virtual void partial_display(ostream &) const
Printing of some informations, excluding all global quantities.
Definition etoile_rot.C:630
void fait_shift()
Computes shift from w_shift and khi_shift according to Shibata's prescription [Prog.
Definition etoile_rot.C:748
Base class for stars *** DEPRECATED : use class Star instead ***.
Definition etoile.h:424
virtual void equilibrium_spher(double ent_c, double precis=1.e-14, const Tbl *ent_limit=0x0)
Computes a spherical static configuration.
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition etoile.C:396
Tenseur logn_auto_div
Divergent part (if k_div!=0 ) of the logarithm of the part of the lapse N generated principaly by t...
Definition etoile.h:497
int nzet
Number of domains of *mp occupied by the star.
Definition etoile.h:432
const Tenseur & get_nnn() const
Returns the total lapse function N.
Definition etoile.h:727
double * p_mass_b
Baryon mass.
Definition etoile.h:547
int get_nzet() const
Returns the number of domains occupied by the star.
Definition etoile.h:662
double * p_mass_g
Gravitational mass.
Definition etoile.h:548
Tbl * p_xi_surf
Description of the stellar surface: 2-D Tbl containing the values of the radial coordinate on the su...
Definition etoile.h:545
void operator=(const Etoile &)
Assignment to another Etoile.
Definition etoile.C:430
Tenseur logn_auto_regu
Regular part of the logarithm of the part of the lapse N generated principaly by the star.
Definition etoile.h:491
double ray_eq_pi() const
Coordinate radius at , [r_unit].
const Tenseur & get_shift() const
Returns the total shift vector .
Definition etoile.h:730
int k_div
Index of regularity of the gravitational potential logn_auto .
Definition etoile.h:449
double ray_eq() const
Coordinate radius at , [r_unit].
Tenseur nnn
Total lapse function.
Definition etoile.h:509
const Tenseur & get_logn_auto_regu() const
Returns the regular part of the logarithm of the part of the lapse N generated principaly by the star...
Definition etoile.h:708
Tenseur logn_auto
Total of the logarithm of the part of the lapse N generated principaly by the star.
Definition etoile.h:484
const Map & get_mp() const
Returns the mapping.
Definition etoile.h:659
Tenseur nbar
Baryon density in the fluid frame.
Definition etoile.h:459
const Eos & eos
Equation of state of the stellar matter.
Definition etoile.h:451
bool is_relativistic() const
Returns true for a relativistic star, false for a Newtonian one.
Definition etoile.h:667
virtual void equation_of_state()
Computes the proper baryon and energy density, as well as pressure from the enthalpy.
Definition etoile.C:566
Tenseur u_euler
Fluid 3-velocity with respect to the Eulerian observer.
Definition etoile.h:474
Itbl * p_l_surf
Description of the stellar surface: 2-D Itbl containing the values of the domain index l on the surfa...
Definition etoile.h:539
Tenseur gam_euler
Lorentz factor between the fluid and Eulerian observers.
Definition etoile.h:471
Map & mp
Mapping associated with the star.
Definition etoile.h:429
double * p_ray_eq_3pis2
Coordinate radius at , .
Definition etoile.h:530
virtual void del_deriv() const
Deletes all the derived quantities.
Definition etoile.C:378
virtual void equil_spher_falloff(double ent_c, double precis=1.e-14)
Computes a spherical static configuration with the outer boundary condition at a finite radius.
double * p_ray_eq_pi
Coordinate radius at , .
Definition etoile.h:527
const Tenseur & get_ent() const
Returns the enthalpy field.
Definition etoile.h:673
const Eos & get_eos() const
Returns the equation of state.
Definition etoile.h:670
virtual double mass_b() const
Baryon mass.
const Tenseur & get_beta_auto() const
Returns the logarithm of the part of the product AN generated principaly by the star.
Definition etoile.h:724
const Tenseur & get_d_logn_auto_div() const
Returns the gradient of logn_auto_div.
Definition etoile.h:719
Tenseur ener
Total energy density in the fluid frame.
Definition etoile.h:460
friend ostream & operator<<(ostream &, const Etoile &)
Display.
Definition etoile.C:506
void equil_spher_regular(double ent_c, double precis=1.e-14)
Computes a spherical static configuration.
double ray_eq_pis2() const
Coordinate radius at , [r_unit].
const Tenseur & get_a_car() const
Returns the total conformal factor .
Definition etoile.h:733
virtual ~Etoile()
Destructor.
Definition etoile.C:367
double * p_ray_pole
Coordinate radius at .
Definition etoile.h:533
const Tenseur & get_logn_auto_div() const
Returns the divergent part of the logarithm of the part of the lapse N generated principaly by the st...
Definition etoile.h:715
Tenseur d_logn_auto_div
Gradient of logn_auto_div (if k_div!=0 )
Definition etoile.h:501
virtual const Itbl & l_surf() const
Description of the stellar surface: returns a 2-D Itbl containing the values of the domain index l ...
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition etoile.C:511
double ray_eq_3pis2() const
Coordinate radius at , [r_unit].
Map & set_mp()
Read/write of the mapping.
Definition etoile.h:601
Tenseur press
Fluid pressure.
Definition etoile.h:461
virtual double mass_g() const
Gravitational mass.
bool relativistic
Indicator of relativity: true for a relativistic star, false for a Newtonian one.
Definition etoile.h:437
const Tenseur & get_logn_auto() const
Returns the logarithm of the part of the lapse N generated principaly by the star.
Definition etoile.h:701
const Tenseur & get_u_euler() const
Returns the fluid 3-velocity with respect to the Eulerian observer.
Definition etoile.h:694
double * p_ray_eq_pis2
Coordinate radius at , .
Definition etoile.h:524
const Tenseur & get_s_euler() const
Returns the trace of the stress tensor in the Eulerian frame.
Definition etoile.h:688
Tenseur ener_euler
Total energy density in the Eulerian frame.
Definition etoile.h:465
virtual void sauve(FILE *) const
Save in a file.
Definition etoile.C:483
Tenseur shift
Total shift vector.
Definition etoile.h:512
virtual void hydro_euler()
Computes the hydrodynamical quantities relative to the Eulerian observer from those in the fluid fram...
Definition etoile.C:684
Tenseur s_euler
Trace of the stress tensor in the Eulerian frame.
Definition etoile.h:468
Tenseur ent
Log-enthalpy (relativistic case) or specific enthalpy (Newtonian case)
Definition etoile.h:457
virtual void del_hydro_euler()
Sets to ETATNONDEF (undefined state) the hydrodynamical quantities relative to the Eulerian observer.
Definition etoile.C:410
Tenseur beta_auto
Logarithm of the part of the product AN generated principaly by by the star.
Definition etoile.h:506
const Tenseur & get_ener() const
Returns the proper total energy density.
Definition etoile.h:679
const Tbl & xi_surf() const
Description of the stellar surface: returns a 2-D Tbl containing the values of the radial coordinat...
const Tenseur & get_press() const
Returns the fluid pressure.
Definition etoile.h:682
Tenseur a_car
Total conformal factor .
Definition etoile.h:515
void set_enthalpy(const Cmp &)
Assignment of the enthalpy field.
Definition etoile.C:465
const Tenseur & get_ener_euler() const
Returns the total energy density with respect to the Eulerian observer.
Definition etoile.h:685
double unsurc2
: unsurc2=1 for a relativistic star, 0 for a Newtonian one.
Definition etoile.h:442
double * p_ray_eq
Coordinate radius at , .
Definition etoile.h:521
const Tenseur & get_nbar() const
Returns the proper baryon density.
Definition etoile.h:676
const Tenseur & get_gam_euler() const
Returns the Lorentz factor between the fluid and Eulerian observers.
Definition etoile.h:691
double ray_pole() const
Coordinate radius at [r_unit].
Basic integer array class.
Definition itbl.h:122
Base class for coordinate mappings.
Definition map.h:670
Parameter storage.
Definition param.h:125
Basic array class.
Definition tbl.h:161
Class intended to describe tensors with a symmetry on the two last indices *** DEPRECATED : use class...
Definition tenseur.h:1253
Tensor handling *** DEPRECATED : use class Tensor instead ***.
Definition tenseur.h:301
Lorene prototypes.
Definition app_hor.h:64