LORENE
star.h
1/*
2 * Definition of Lorene classes Star, Star_bin and Star_bin_xcts
3 *
4 */
5
6/*
7 * Copyright (c) 2010 Michal Bejger
8 * Copyright (c) 2004 Francois Limousin
9 *
10 * Copyright (c) 2000-2001 Eric Gourgoulhon (for preceding class Etoile)
11 * Copyright (c) 2000-2001 Keisuke Taniguchi (for preceding class Etoile)
12 *
13 * This file is part of LORENE.
14 *
15 * LORENE is free software; you can redistribute it and/or modify
16 * it under the terms of the GNU General Public License as published by
17 * the Free Software Foundation; either version 2 of the License, or
18 * (at your option) any later version.
19 *
20 * LORENE is distributed in the hope that it will be useful,
21 * but WITHOUT ANY WARRANTY; without even the implied warranty of
22 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
23 * GNU General Public License for more details.
24 *
25 * You should have received a copy of the GNU General Public License
26 * along with LORENE; if not, write to the Free Software
27 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
28 *
29 */
30
31
32#ifndef __STAR_H_
33#define __STAR_H_
34
35/*
36 * $Id: star.h,v 1.32 2014/10/13 08:52:36 j_novak Exp $
37 * $Log: star.h,v $
38 * Revision 1.32 2014/10/13 08:52:36 j_novak
39 * Lorene classes and functions now belong to the namespace Lorene.
40 *
41 * Revision 1.31 2010/12/09 10:36:42 m_bejger
42 * Decouple method removed from declaration in Star_bin_xcts
43 *
44 * Revision 1.30 2010/10/26 18:47:36 m_bejger
45 * Modification of Star_bin_xcts::equilibrium, added table fact_resize
46 *
47 * Revision 1.29 2010/10/18 19:11:53 m_bejger
48 * Changes to Star::equilibrium_spher and Star_bin_xcts::equilibrium to allow for calculations with more than one domain in the star
49 *
50 * Revision 1.28 2010/06/15 08:17:43 m_bejger
51 * Method Star_bin_xcts::set_chi_comp() declared
52 *
53 * Revision 1.27 2010/05/04 07:53:32 m_bejger
54 * Class Star_bin_xcts added (initial version)
55 *
56 * Revision 1.26 2010/03/29 12:00:25 e_gourgoulhon
57 * Removed lnq from documentation
58 * (lnq has to be removed from base class Star; it is meaningfull only
59 * for Star_bin).
60 *
61 * Revision 1.25 2010/01/24 16:07:45 e_gourgoulhon
62 * New class Star_rot.
63 *
64 * Revision 1.24 2007/11/06 16:22:03 j_novak
65 * The data member stress_euler is now a Sym_tensor instead of a Tensor.
66 *
67 * Revision 1.23 2007/06/21 19:48:25 k_taniguchi
68 * Introduction of a method to compute ray_eq_3pis2.
69 *
70 * Revision 1.22 2006/05/31 09:25:47 f_limousin
71 * Modif. of the size of the different domains
72 *
73 * Revision 1.21 2006/04/11 14:26:12 f_limousin
74 * New version of the code : improvement of the computation of some
75 * critical sources, estimation of the dirac gauge, helical symmetry...
76 *
77 * Revision 1.20 2005/09/15 15:56:28 e_gourgoulhon
78 * Made the documentation compliant with Doxygen.
79 *
80 * Revision 1.19 2005/09/13 19:38:32 f_limousin
81 * Reintroduction of the resolution of the equations in cartesian coordinates.
82 *
83 * Revision 1.18 2005/08/13 16:14:11 m_saijo
84 * Corrected the document of the total shift vector
85 *
86 * Revision 1.17 2005/04/08 12:36:45 f_limousin
87 * Just to avoid warnings...
88 *
89 * Revision 1.16 2005/02/24 16:09:29 f_limousin
90 * Change the name of some variables (for instance dcov_logn --> dlogn).
91 * Add also member dlnq but delete dlnpsi_auto and dlogn_auto.
92 *
93 * Revision 1.15 2005/02/17 17:28:18 f_limousin
94 * Change the name of some quantities to be consistent with other classes
95 * (for instance nnn is changed to nn, shift to beta, beta to lnq...)
96 *
97 * Revision 1.14 2005/02/11 18:11:16 f_limousin
98 * Introduction of a member Map_af in derived class Star_bin.
99 *
100 * Revision 1.13 2004/11/11 16:29:48 j_novak
101 * set_der_0x0 is no longer virtual (to be coherent with Tensor/Scalar classes).
102 *
103 * Revision 1.12 2004/11/10 16:31:53 j_novak
104 * Star is no longer an abstract class (mass_b and mass_g are no longer
105 * pure virtual). Modified comments to be readable by doxygen.
106 *
107 * Revision 1.11 2004/07/21 11:48:30 f_limousin
108 * Remove function sprod.
109 *
110 * Revision 1.10 2004/06/22 12:47:01 f_limousin
111 * Change qq, qq_auto and qq_comp to beta, beta_auto and beta_comp.
112 *
113 * Revision 1.9 2004/05/25 14:48:57 f_limousin
114 * Add a parameter for the function equilibrium.
115 *
116 * Revision 1.8 2004/03/23 09:53:50 f_limousin
117 * Minor changes
118 *
119 * Revision 1.7 2004/02/27 09:41:52 f_limousin
120 * Scalars ssjm1_logn, ssjm1_qq ... for all metric coefficients have been
121 * in class Star_bin for the resolution of Poisson equations.
122 * The class Star is now abstract : the computational routines mass_b()
123 * and mass_g() = 0.
124 *
125 * Revision 1.6 2004/01/22 10:06:33 f_limousin
126 * Add methods set_logn_comp() and set_shift_auto().
127 *
128 * Revision 1.5 2004/01/20 15:26:00 f_limousin
129 * New class star and star_bin.
130 *
131 *
132 * $Header: /cvsroot/Lorene/C++/Include/star.h,v 1.32 2014/10/13 08:52:36 j_novak Exp $
133 *
134 */
135
136// Headers Lorene
137#include "tensor.h"
138#include "metric.h"
139
140namespace Lorene {
141class Eos ;
142
143 //---------------------------//
144 // base class Star //
145 //---------------------------//
146
175class Star {
176
177 // Data :
178 // -----
179 protected:
180 Map& mp ;
181
183 int nzet ;
184
185 const Eos& eos ;
186
187 // Fluid quantities with respect to the fluid frame
188 // ------------------------------------------------
189
191
195
196 // Fluid quantities with respect to the Eulerian frame
197 // ---------------------------------------------------
199
202
205
208
213
214 // Metric potentials
215 // -----------------
216
223
226
229
230 // Scalar field \f$\ln Q = \ln(\psi^2 N)\f$
231 //## to be removed from base class Star
232 Scalar lnq ;
233
236
237
238 // Derived data :
239 // ------------
240 protected:
242 mutable double* p_ray_eq ;
243
245 mutable double* p_ray_eq_pis2 ;
246
248 mutable double* p_ray_eq_pi ;
249
251 mutable double* p_ray_eq_3pis2 ;
252
254 mutable double* p_ray_pole ;
255
260 mutable Itbl* p_l_surf ;
261
266 mutable Tbl* p_xi_surf ;
267
268 mutable double* p_mass_b ;
269 mutable double* p_mass_g ;
270
271
272 // Constructors - Destructor
273 // -------------------------
274 public:
275
283 Star(Map& mp_i, int nzet_i, const Eos& eos_i) ;
284
285
286 Star(const Star& ) ;
287
295 Star(Map& mp_i, const Eos& eos_i, FILE* fich) ;
296
297 virtual ~Star() ;
298
299
300 // Memory management
301 // -----------------
302 protected:
304 virtual void del_deriv() const ;
305
307 virtual void set_der_0x0() const ;
308
312 virtual void del_hydro_euler() ;
313
314
315 // Mutators / assignment
316 // ---------------------
317 public:
319 void operator=(const Star&) ;
320
322 Map& set_mp() {return mp; } ;
323
325 void set_enthalpy(const Scalar& ) ;
326
330 void equation_of_state() ;
331
336 virtual void hydro_euler() ;
337
348 virtual void equilibrium_spher(double ent_c, double precis = 1.e-14,
349 const Tbl* pent_limit = 0x0 ) ;
350
351 // Accessors
352 // ---------
353 public:
355 const Map& get_mp() const {return mp; } ;
356
358 int get_nzet() const {return nzet; } ;
359
361 const Eos& get_eos() const {return eos; } ;
362
364 const Scalar& get_ent() const {return ent;} ;
365
367 const Scalar& get_nbar() const {return nbar;} ;
368
370 const Scalar& get_ener() const {return ener;} ;
371
373 const Scalar& get_press() const {return press;} ;
374
376 const Scalar& get_ener_euler() const {return ener_euler;} ;
377
379 const Scalar& get_s_euler() const {return s_euler;} ;
380
382 const Scalar& get_gam_euler() const {return gam_euler;} ;
383
385 const Vector& get_u_euler() const {return u_euler;} ;
386
390 const Tensor& get_stress_euler() const {return stress_euler;} ;
391
396 const Scalar& get_logn() const {return logn;} ;
397
399 const Scalar& get_nn() const {return nn;} ;
400
402 const Vector& get_beta() const {return beta;} ;
403
404 // Returns the scalar field \f$\ln Q\f$.
405 //## to be removed from base class Star
406 const Scalar& get_lnq() const {return lnq;} ;
407
409 const Metric& get_gamma() const {return gamma;} ;
410
411 // Outputs
412 // -------
413 public:
414 virtual void sauve(FILE* ) const ;
415
417 friend ostream& operator<<(ostream& , const Star& ) ;
418
419 protected:
421 virtual ostream& operator>>(ostream& ) const ;
422
423 // Global quantities
424 // -----------------
425 public:
427 double ray_eq() const ;
428
430 double ray_eq_pis2() const ;
431
433 double ray_eq_pi() const ;
434
436 double ray_eq_3pis2() const ;
437
439 double ray_pole() const ;
440
448 virtual const Itbl& l_surf() const ;
449
457 const Tbl& xi_surf() const ;
458
460 virtual double mass_b() const = 0 ;
461
463 virtual double mass_g() const = 0 ;
464
465};
466ostream& operator<<(ostream& , const Star& ) ;
467
468
469 //---------------------------//
470 // class Star_bin //
471 //---------------------------//
472
483class Star_bin : public Star {
484
485 // Data :
486 // -----
487 protected:
492
497
502
508
513
519
522
523
528
533
536
539
544
549
550
553
558
563
566
571
576
577
582
583
589
595
601
607
613
619
624
629
635
636 Vector ssjm1_wbeta ;
637
642
647
652
657
662
667
677
682
683 // Derived data :
684 // ------------
685 protected:
687 mutable double* p_xa_barycenter ;
688
689
690 // Constructors - Destructor
691 // -------------------------
692 public:
703 Star_bin(Map& mp_i, int nzet_i, const Eos& eos_i,
704 bool irrot, bool conf_flat) ;
705
706
707 Star_bin(const Star_bin& ) ;
708
716 Star_bin(Map& mp_i, const Eos& eos_i, FILE* fich) ;
717
718 virtual ~Star_bin() ;
719
720
721 // Memory management
722 // -----------------
723 protected:
725 virtual void del_deriv() const ;
726
728 virtual void set_der_0x0() const ;
729
733 virtual void del_hydro_euler() ;
734
735
736 // Mutators / assignment
737 // ---------------------
738 public:
740 void operator=(const Star_bin& ) ;
741
744
749
751 void set_logn_auto(const Scalar& logn_auto_new) {logn_auto = logn_auto_new ;
752 return ;}
753
755 void set_lnq_auto(const Scalar& lnq_auto_new) {lnq_auto = lnq_auto_new ;
756 return ;}
757
760
762 Vector& set_beta() ;
763
765 void set_conf_flat(bool confflat) {conf_flat = confflat ; return ;}
766
767 // Accessors
768 // ---------
769 public:
773 bool is_irrotational() const {return irrotational; } ;
774
776 const Scalar& get_psi0() const {return psi0;} ;
777
781 const Vector& get_d_psi() const {return d_psi;} ;
782
787 const Vector& get_wit_w() const {return wit_w;} ;
788
792 const Scalar& get_loggam() const {return loggam;} ;
793
798 const Vector& get_bsn() const {return bsn;} ;
799
801 const Scalar& get_pot_centri() const {return pot_centri;} ;
802
806 const Scalar& get_logn_auto() const {return logn_auto;} ;
807
811 const Scalar& get_logn_comp() const {return logn_comp;} ;
812
817 const Vector& get_beta_auto() const {return beta_auto;} ;
818
823 const Vector& get_beta_comp() const {return beta_comp;} ;
824
828 const Scalar& get_lnq_auto() const {return lnq_auto;} ;
829
833 const Scalar& get_lnq_comp() const {return lnq_comp;} ;
834
837 const Vector& get_dcov_logn() const {return dcov_logn;} ;
838
841 const Vector& get_dcon_logn() const {return dcon_logn;} ;
842
846 const Vector& get_dcov_phi() const {return dcov_phi;} ;
847
851 const Vector& get_dcon_phi() const {return dcon_phi;} ;
852
854 const Scalar& get_psi4() const {return psi4;} ;
855
859 const Metric& get_flat() const {return flat;} ;
860
862 const Metric& get_gtilde() const {return gtilde;} ;
863
867 const Sym_tensor& get_hij() const {return hij;} ;
868
873 const Sym_tensor& get_hij_auto() const {return hij_auto;} ;
874
879 const Sym_tensor& get_hij_comp() const {return hij_comp;} ;
880
881
886 const Sym_tensor& get_tkij_auto() const {return tkij_auto;} ;
887
892 const Sym_tensor& get_tkij_comp() const {return tkij_comp;} ;
893
897 const Scalar& get_kcar_auto() const {return kcar_auto;} ;
898
902 const Scalar& get_kcar_comp() const {return kcar_comp;} ;
903
908 const Scalar get_decouple() const {return decouple ;}
909
910 bool is_conf_flat() const {return conf_flat; } ;
911
912 // Outputs
913 // -------
914 public:
915 virtual void sauve(FILE* ) const ;
916
917 protected:
919 virtual ostream& operator>>(ostream& ) const ;
920
921 // Global quantities
922 // -----------------
923 public:
925 virtual double mass_b() const ;
926
928 virtual double mass_g() const ;
929
931 virtual double xa_barycenter() const ;
932
933
934 // Computational routines
935 // ----------------------
936 public:
937
951 virtual void hydro_euler() ;
952
969 void update_metric(const Star_bin& comp, double omega) ;
970
980 void update_metric(const Star_bin& comp, const Star_bin& star_prev,
981 double relax, double omega) ;
982
988 void update_metric_der_comp(const Star_bin& comp, double omega) ;
989
1000 void kinematics(double omega, double x_axe) ;
1001
1005 void fait_d_psi() ;
1006
1012 void extrinsic_curvature(double omega) ;
1013
1014
1031 void equilibrium(double ent_c, int mermax, int mermax_potvit,
1032 int mermax_poisson, double relax_poisson,
1033 double relax_potvit, double thres_adapt, Tbl& diff,
1034 double om) ;
1035
1036
1050 double velocity_potential(int mermax, double precis, double relax) ;
1051
1063 void relaxation(const Star_bin& star_prev, double relax_ent,
1064 double relax_met, int mer, int fmer_met) ;
1065
1066
1068 void test_K_Hi() const ;
1069
1071 void helical(double omega) const ;
1072
1073 friend class Binary ;
1074
1075};
1076
1077 //---------------------------//
1078 // class Star_bin_xcts //
1079 //---------------------------//
1080
1091class Star_bin_xcts : public Star {
1092
1093 // Data :
1094 // -----
1095 protected:
1100
1105
1110
1116
1121
1127
1130
1135
1140
1145
1150
1153
1156
1159
1164
1169
1172
1175
1178
1183
1188
1194
1200
1206
1212
1217
1222
1228
1235
1236 // Derived data :
1237 // ------------
1238 protected:
1240 mutable double* p_xa_barycenter ;
1241
1242 // Constructors - Destructor
1243 // -------------------------
1244 public:
1253 Star_bin_xcts(Map& mp_i, int nzet_i, const Eos& eos_i, bool irrot) ;
1254
1255 Star_bin_xcts(const Star_bin_xcts& ) ;
1256
1264 Star_bin_xcts(Map& mp_i, const Eos& eos_i, FILE* fich) ;
1265
1266 virtual ~Star_bin_xcts() ;
1267
1268
1269 // Memory management
1270 // -----------------
1271 protected:
1273 virtual void del_deriv() const ;
1274
1276 virtual void set_der_0x0() const ;
1277
1281 virtual void del_hydro_euler() ;
1282
1283
1284 // Mutators / assignment
1285 // ---------------------
1286 public:
1288 void operator=(const Star_bin_xcts& ) ;
1289
1292
1294 Scalar& set_Psi_auto() ;
1295
1297 Scalar& set_Psi_comp() ;
1298
1300 Scalar& set_chi_auto() ;
1301
1303 Scalar& set_chi_comp() ;
1304
1306 Vector& set_beta_auto() ;
1307
1309 Vector& set_beta() ;
1310
1311 // Accessors
1312 // ---------
1313 public:
1317 bool is_irrotational() const { return irrotational; } ;
1318
1320 const Scalar& get_psi0() const {return psi0 ; } ;
1321
1325 const Vector& get_d_psi() const {return d_psi ; } ;
1326
1331 const Vector& get_wit_w() const {return wit_w ; } ;
1332
1336 const Scalar& get_loggam() const {return loggam ; } ;
1337
1342 const Vector& get_bsn() const {return bsn ; } ;
1343
1345 const Scalar& get_pot_centri() const {return pot_centri ; } ;
1346
1351 const Vector& get_beta_auto() const {return beta_auto ; } ;
1352
1357 const Vector& get_beta_comp() const {return beta_comp ; } ;
1358
1362 const Scalar& get_Psi_auto() const {return Psi_auto ; } ;
1363
1367 const Scalar& get_Psi_comp() const {return Psi_comp ; } ;
1368
1372 const Scalar& get_chi_auto() const {return chi_auto ; } ;
1373
1377 const Scalar& get_chi_comp() const {return chi_comp ; } ;
1378
1381 const Vector& get_dcov_chi() const {return dcov_chi ; } ;
1382
1386 const Vector& get_dcov_Psi() const {return dcov_Psi ; } ;
1387
1389 const Scalar& get_Psi() const {return Psi ; } ;
1390
1392 const Scalar& get_chi() const {return chi ; } ;
1393
1395 const Scalar& get_psi4() const {return psi4 ; } ;
1396
1400 const Metric& get_flat() const {return flat ; } ;
1401
1406 const Sym_tensor& get_haij_auto() const {return haij_auto ; } ;
1407
1412 const Sym_tensor& get_haij_comp() const {return haij_comp ; } ;
1413
1417 const Scalar& get_hacar_auto() const {return hacar_auto ; } ;
1418
1422 const Scalar& get_hacar_comp() const {return hacar_comp ; } ;
1423
1424 // Outputs
1425 // -------
1426 public:
1427 virtual void sauve(FILE* ) const ;
1428
1429 protected:
1431 virtual ostream& operator>>(ostream& ) const ;
1432
1433 // Global quantities
1434 // -----------------
1435 public:
1437 virtual double mass_b() const ;
1438
1440 virtual double mass_g() const ;
1441
1443 virtual double xa_barycenter() const ;
1444
1445
1446 // Computational routines
1447 // ----------------------
1448 public:
1449
1463 virtual void hydro_euler() ;
1464
1481 void update_metric(const Star_bin_xcts& comp) ;
1482
1492 void update_metric(const Star_bin_xcts& comp,
1493 const Star_bin_xcts& star_prev,
1494 double relax) ;
1495
1499 void update_metric_der_comp(const Star_bin_xcts& comp) ;
1500
1511 void kinematics(double omega, double x_axe) ;
1512
1516 void fait_d_psi() ;
1517
1521 void extrinsic_curvature() ;
1522
1544 void equilibrium(double ent_c, int mermax, int mermax_potvit,
1545 int mermax_poisson, double relax_poisson,
1546 double relax_potvit, double thres_adapt,
1547 const Tbl& fact, const Tbl* pent_limit,
1548 Tbl& diff) ;
1549
1563 double velocity_potential(int mermax, double precis, double relax) ;
1564
1576 void relaxation(const Star_bin_xcts& star_prev, double relax_ent,
1577 double relax_met, int mer, int fmer_met) ;
1578
1579 friend class Binary_xcts ;
1580
1581};
1582
1583}
1584#endif
Binary systems in eXtended Conformal Thin Sandwich formulation.
Definition binary_xcts.h:59
Binary systems.
Definition binary.h:73
Equation of state base class.
Definition eos.h:190
Basic integer array class.
Definition itbl.h:122
Base class for coordinate mappings.
Definition map.h:670
Flat metric for tensor calculation.
Definition metric.h:261
Metric for tensor calculation.
Definition metric.h:90
Tensor field of valence 0 (or component of a tensorial field).
Definition scalar.h:387
Class for stars in binary system in eXtended Conformal Thin Sandwich formulation.
Definition star.h:1091
Scalar Psi
Total conformal factor .
Definition star.h:1152
const Vector & get_beta_comp() const
Returns the part of the shift vector generated principally by the companion ( ).
Definition star.h:1357
Scalar psi0
Scalar potential of the non-translational part of fluid 4-velocity (in the irrotational case)
Definition star.h:1104
Scalar & set_pot_centri()
Read/write the centrifugal potential.
const Scalar & get_chi() const
Return the function .
Definition star.h:1392
const Vector & get_dcov_chi() const
Returns the covariant derivative of .
Definition star.h:1381
void operator=(const Star_bin_xcts &)
Assignment to another Star_bin_xcts.
const Sym_tensor & get_haij_auto() const
Returns the part of the extrinsic curvature tensor generated by beta_auto.
Definition star.h:1406
Scalar pot_centri
Centrifugal potential.
Definition star.h:1129
Sym_tensor haij_auto
Part of the extrinsic curvature tensor generated by beta_auto.
Definition star.h:1193
Vector & set_beta_auto()
Read/write of .
void fait_d_psi()
Computes the gradient of the total velocity potential .
const Scalar & get_hacar_comp() const
Returns the part of generated by beta_comp.
Definition star.h:1422
void extrinsic_curvature()
Computes haij_auto and hacar_auto from beta_auto, nn and Psi .
const Vector & get_wit_w() const
Returns the spatial projection of the fluid 3-velocity with respect to the co-orbiting observer.
Definition star.h:1331
Scalar psi4
Conformal factor .
Definition star.h:1158
Vector wit_w
Spatial projection of the fluid 3-velocity with respect to the co-orbiting observer.
Definition star.h:1115
virtual double xa_barycenter() const
Absolute coordinate X of the barycenter of the baryon density,.
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Metric_flat flat
Flat metric defined on the mapping (Spherical components with respect to the mapping of the star) .
Definition star.h:1177
Vector & set_beta()
Read/write of .
void equilibrium(double ent_c, int mermax, int mermax_potvit, int mermax_poisson, double relax_poisson, double relax_potvit, double thres_adapt, const Tbl &fact, const Tbl *pent_limit, Tbl &diff)
Computes an equilibrium configuration.
Scalar ssjm1_khi
Effective source at the previous step for the resolution of the Poisson equation for khi.
Definition star.h:1227
Scalar khi
Solution for the scalar part of the vector Poisson equation for .
Definition star.h:1168
Vector ssjm1_wbeta
Effective source at the previous step for wbeta of the vector Poisson equation for wbeta (needed for ...
Definition star.h:1234
Scalar Psi_comp
Scalar field generated principally by the companion star.
Definition star.h:1149
Scalar ssjm1_chi
Effective source at the previous step for the resolution of the Poisson equation for \chi_auto .
Definition star.h:1216
const Scalar & get_Psi() const
Return the conformal factor .
Definition star.h:1389
Vector dcov_Psi
Covariant derivative of the conformal factor .
Definition star.h:1171
Vector beta_comp
Part of the shift vector generated principally by the star (Spherical components with respect to the ...
Definition star.h:1187
const Scalar & get_hacar_auto() const
Returns the part of generated by beta_auto.
Definition star.h:1417
const Scalar & get_chi_auto() const
Returns the scalar field generated principally by the star.
Definition star.h:1372
const Vector & get_dcov_Psi() const
Returns the covariant derivative of the conformal factor .
Definition star.h:1386
Vector beta_auto
Part of the shift vector generated principally by the star (Spherical components with respect to the ...
Definition star.h:1182
Scalar Psi_auto
Scalar field generated principally by the star.
Definition star.h:1144
Scalar hacar_comp
Part of the scalar generated by beta_auto and beta_comp, i.e.
Definition star.h:1211
const Scalar & get_chi_comp() const
Returns the scalar field generated principally by the companion star.
Definition star.h:1377
Scalar chi
Total function .
Definition star.h:1155
void update_metric_der_comp(const Star_bin_xcts &comp)
Computes the derivative of metric functions related to the companion star.
Scalar & set_chi_comp()
Read/write the conformal factor .
const Scalar & get_Psi_comp() const
Returns the scalar field generated principally by the companion star.
Definition star.h:1367
void kinematics(double omega, double x_axe)
Computes the quantities bsn and pot_centri.
const Scalar & get_psi0() const
Returns the non-translational part of the velocity potential.
Definition star.h:1320
const Vector & get_beta_auto() const
Returns the part of the shift vector generated principally by the star ( ).
Definition star.h:1351
const Scalar & get_psi4() const
Return the conformal factor .
Definition star.h:1395
Vector bsn
3-vector shift, divided by N, of the rotating coordinates, .
Definition star.h:1126
virtual void del_hydro_euler()
Sets to ETATNONDEF (undefined state) the hydrodynamical quantities relative to the Eulerian observer.
Scalar loggam
Logarithm of the Lorentz factor between the fluid and the co-orbiting observer.
Definition star.h:1120
Scalar & set_chi_auto()
Read/write the conformal factor .
const Scalar & get_pot_centri() const
Returns the centrifugal potential.
Definition star.h:1345
Sym_tensor haij_comp
Part of the extrinsic curvature tensor generated by beta_comp.
Definition star.h:1199
Scalar & set_Psi_comp()
Read/write the conformal factor .
Scalar & set_Psi_auto()
Read/write the conformal factor .
bool is_irrotational() const
Returns true for an irrotational motion, false for a corotating one.
Definition star.h:1317
const Vector & get_d_psi() const
Returns the covariant derivative of the velocity potential (Spherical components with respect to the ...
Definition star.h:1325
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Vector dcov_chi
Covariant derivative of the function .
Definition star.h:1174
const Metric & get_flat() const
Return the flat metric defined on the mapping (Spherical components with respect to the mapping of th...
Definition star.h:1400
Scalar chi_auto
Scalar field generated principally by the star.
Definition star.h:1134
virtual double mass_g() const
Gravitational mass.
void update_metric(const Star_bin_xcts &comp)
Computes metric coefficients from known potentials, when the companion is another star.
Vector w_beta
Solution for the vector part of the vector Poisson equation for .
Definition star.h:1163
const Vector & get_bsn() const
Returns the shift vector, divided by N, of the rotating coordinates, .
Definition star.h:1342
Scalar chi_comp
Scalar field generated principally by the companion star.
Definition star.h:1139
Scalar ssjm1_psi
Effective source at the previous step for the resolution of the Poisson equation for \Psi_auto .
Definition star.h:1221
const Sym_tensor & get_haij_comp() const
Returns the part of the extrinsic curvature tensor generated by beta_comp.
Definition star.h:1412
virtual void hydro_euler()
Computes the hydrodynamical quantities relative to the Eulerian observer from those in the fluid fram...
void relaxation(const Star_bin_xcts &star_prev, double relax_ent, double relax_met, int mer, int fmer_met)
Performs a relaxation on ent, Psi_auto, chi_auto and beta_auto.
bool irrotational
true for an irrotational star, false for a corotating one
Definition star.h:1099
virtual void del_deriv() const
Deletes all the derived quantities.
Vector d_psi
Gradient of (in the irrotational case) (Spherical components with respect to the mapping of the star...
Definition star.h:1109
virtual void sauve(FILE *) const
Save in a file.
virtual double mass_b() const
Baryon mass.
const Scalar & get_loggam() const
Returns the logarithm of the Lorentz factor between the fluid and the co-orbiting observer.
Definition star.h:1336
virtual ~Star_bin_xcts()
Destructor.
double * p_xa_barycenter
Absolute coordinate X of the barycenter of the baryon density.
Definition star.h:1240
double velocity_potential(int mermax, double precis, double relax)
Computes the non-translational part of the velocity scalar potential by solving the continuity equat...
Scalar hacar_auto
Part of the scalar generated by beta_auto, i.e.
Definition star.h:1205
const Scalar & get_Psi_auto() const
Returns the scalar field generated principally by the star.
Definition star.h:1362
Class for stars in binary system.
Definition star.h:483
Scalar loggam
Logarithm of the Lorentz factor between the fluid and the co-orbiting observer.
Definition star.h:512
Scalar decouple
Function used to construct the part generated by the star from the total .
Definition star.h:676
Scalar lnq_auto
Scalar field generated principally by the star.
Definition star.h:543
Scalar ssjm1_khi
Effective source at the previous step for the resolution of the Poisson equation for khi.
Definition star.h:634
Vector bsn
3-vector shift, divided by N, of the rotating coordinates, .
Definition star.h:518
Vector dcon_logn
Contravariant derivative of the total logarithm of the lapse.
Definition star.h:538
Vector d_psi
Gradient of (in the irrotational case) (Spherical components with respect to the mapping of the star...
Definition star.h:501
Vector dcon_phi
Contravariant derivative of the logarithm of the conformal factor.
Definition star.h:557
double velocity_potential(int mermax, double precis, double relax)
Computes the non-translational part of the velocity scalar potential by solving the continuity equat...
const Scalar & get_logn_auto() const
Returns the part of the lapse logarithm (gravitational potential at the Newtonian limit) generated pr...
Definition star.h:806
Sym_tensor hij_comp
Deviation of the inverse conformal metric from the inverse flat metric generated principally by the ...
Definition star.h:594
double * p_xa_barycenter
Absolute coordinate X of the barycenter of the baryon density.
Definition star.h:687
bool is_irrotational() const
Returns true for an irrotational motion, false for a corotating one.
Definition star.h:773
virtual void del_deriv() const
Deletes all the derived quantities.
Definition star_bin.C:369
Scalar ssjm1_h11
Effective source at the previous step for the resolution of the Poisson equation for h00_auto.
Definition star.h:641
const Scalar & get_lnq_auto() const
Returns the part of the vector field generated principally by the star.
Definition star.h:828
const Vector & get_dcov_logn() const
Returns the covariant derivative of .
Definition star.h:837
Vector beta_comp
Part of the shift vector generated principally by the star (Spherical components with respect to the ...
Definition star.h:575
Scalar ssjm1_logn
Effective source at the previous step for the resolution of the Poisson equation for logn_auto.
Definition star.h:623
Vector dcov_logn
Covariant derivative of the total logarithm of the lapse.
Definition star.h:535
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition star_bin.C:522
void update_metric(const Star_bin &comp, double omega)
Computes metric coefficients from known potentials, when the companion is another star.
Scalar logn_auto
Part of the lapse logarithm (gravitational potential at the Newtonian limit) generated principally by...
Definition star.h:527
Scalar psi0
Scalar potential of the non-translational part of fluid 4-velocity (in the irrotational case)
Definition star.h:496
const Vector & get_bsn() const
Returns the shift vector, divided by N, of the rotating coordinates, .
Definition star.h:798
void kinematics(double omega, double x_axe)
Computes the quantities bsn and pot_centri.
bool irrotational
true for an irrotational star, false for a corotating one
Definition star.h:491
const Vector & get_dcov_phi() const
Returns the covariant derivative of (logarithm of the conformal factor).
Definition star.h:846
virtual double xa_barycenter() const
Absolute coordinate X of the barycenter of the baryon density,.
const Metric & get_flat() const
Return the flat metric defined on the mapping (Spherical components with respect to the mapping of th...
Definition star.h:859
Scalar kcar_comp
Part of the scalar generated by beta_auto and beta_comp, i.e.
Definition star.h:618
Scalar logn_comp
Part of the lapse logarithm (gravitational potential at the Newtonian limit) generated principally by...
Definition star.h:532
Scalar ssjm1_h21
Effective source at the previous step for the resolution of the Poisson equation for h10_auto.
Definition star.h:646
void extrinsic_curvature(double omega)
Computes tkij_auto and akcar_auto from beta_auto, nn and Q.
void fait_d_psi()
Computes the gradient of the total velocity potential .
Definition star_bin.C:645
virtual void hydro_euler()
Computes the hydrodynamical quantities relative to the Eulerian observer from those in the fluid fram...
Sym_tensor tkij_comp
Part of the extrinsic curvature tensor generated by beta_comp.
Definition star.h:606
Scalar ssjm1_h33
Effective source at the previous step for the resolution of the Poisson equation for h22_auto.
Definition star.h:666
const Sym_tensor & get_tkij_comp() const
Returns the part of the extrinsic curvature tensor generated by beta_comp.
Definition star.h:892
void test_K_Hi() const
Test if the gauge conditions we impose are well satisfied.
Definition star_bin.C:727
Scalar lnq_comp
Scalar field generated principally by the companion star.
Definition star.h:548
const Scalar get_decouple() const
Returns the function used to construct beta_auto from beta.
Definition star.h:908
Scalar kcar_auto
Part of the scalar generated by beta_auto, i.e.
Definition star.h:612
void helical(double omega) const
Test of the helical symmetry.
const Sym_tensor & get_hij() const
Return the total deviation of the inverse conformal metric from the inverse flat metric.
Definition star.h:867
const Scalar & get_kcar_auto() const
Returns the part of generated by beta_auto.
Definition star.h:897
Vector beta_auto
Part of the shift vector generated principally by the star (Spherical components with respect to the ...
Definition star.h:570
virtual double mass_b() const
Baryon mass.
const Vector & get_dcon_logn() const
Returns the contravariant derivative of .
Definition star.h:841
Sym_tensor hij
Total deviation of the inverse conformal metric from the inverse flat metric.
Definition star.h:581
void set_lnq_auto(const Scalar &lnq_auto_new)
Assignment of a new lnq_auto.
Definition star.h:755
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition star_bin.C:381
Sym_tensor hij_auto
Deviation of the inverse conformal metric from the inverse flat metric generated principally by the ...
Definition star.h:588
Scalar & set_logn_comp()
Read/write of the logarithm of the lapse generated principally by the companion.
Definition star_bin.C:461
const Scalar & get_psi4() const
Return the conformal factor .
Definition star.h:854
Scalar psi4
Conformal factor .
Definition star.h:552
const Scalar & get_kcar_comp() const
Returns the part of generated by beta_comp.
Definition star.h:902
Scalar ssjm1_h31
Effective source at the previous step for the resolution of the Poisson equation for h20_auto.
Definition star.h:651
void set_conf_flat(bool confflat)
Write if conformally flat.
Definition star.h:765
const Vector & get_dcon_phi() const
Returns the contravariant derivative of (logarithm of the conformal factor).
Definition star.h:851
const Sym_tensor & get_hij_comp() const
Return the deviation of the inverse conformal metric from the inverse flat metric generated principa...
Definition star.h:879
void operator=(const Star_bin &)
Assignment to another Star_bin.
Definition star_bin.C:404
virtual void del_hydro_euler()
Sets to ETATNONDEF (undefined state) the hydrodynamical quantities relative to the Eulerian observer.
Definition star_bin.C:389
const Vector & get_beta_auto() const
Returns the part of the shift vector generated principally by the star.
Definition star.h:817
Vector & set_beta()
Read/write of .
Definition star_bin.C:475
const Scalar & get_lnq_comp() const
Returns the part of the vector field generated principally by the companion star.
Definition star.h:833
void relaxation(const Star_bin &star_prev, double relax_ent, double relax_met, int mer, int fmer_met)
Performs a relaxation on ent, logn_auto, lnq_auto, beta_auto and hij_auto.
Definition star_bin.C:701
const Sym_tensor & get_hij_auto() const
Return the deviation of the inverse conformal metric from the inverse flat metric principally genera...
Definition star.h:873
const Vector & get_d_psi() const
Returns the covariant derivative of the velocity potential (Spherical components with respect to the ...
Definition star.h:781
Sym_tensor tkij_auto
Part of the extrinsic curvature tensor generated by beta_auto.
Definition star.h:600
Vector dcov_phi
Covariant derivative of the logarithm of the conformal factor.
Definition star.h:555
const Sym_tensor & get_tkij_auto() const
Returns the part of the extrinsic curvature tensor generated by beta_auto.
Definition star.h:886
Metric_flat flat
Flat metric defined on the mapping (Spherical components with respect to the mapping of the star) .
Definition star.h:562
Scalar ssjm1_h22
Effective source at the previous step for the resolution of the Poisson equation for h11_auto.
Definition star.h:656
const Vector & get_wit_w() const
Returns the spatial projection of the fluid 3-velocity with respect to the co-orbiting observer.
Definition star.h:787
Vector wit_w
Spatial projection of the fluid 3-velocity with respect to the co-orbiting observer.
Definition star.h:507
void set_logn_auto(const Scalar &logn_auto_new)
Assignment of a new logn_auto.
Definition star.h:751
const Scalar & get_pot_centri() const
Returns the centrifugal potential.
Definition star.h:801
const Scalar & get_psi0() const
Returns the non-translational part of the velocity potential.
Definition star.h:776
Scalar ssjm1_h32
Effective source at the previous step for the resolution of the Poisson equation for h21_auto.
Definition star.h:661
void update_metric_der_comp(const Star_bin &comp, double omega)
Computes the derivative of metric functions related to the companion star.
void equilibrium(double ent_c, int mermax, int mermax_potvit, int mermax_poisson, double relax_poisson, double relax_potvit, double thres_adapt, Tbl &diff, double om)
Computes an equilibrium configuration.
Metric gtilde
Conformal metric .
Definition star.h:565
Scalar pot_centri
Centrifugal potential.
Definition star.h:521
const Vector & get_beta_comp() const
Returns the part of the shift vector generated principally by the star.
Definition star.h:823
Scalar & set_pot_centri()
Read/write the centrifugal potential.
Definition star_bin.C:454
virtual ~Star_bin()
Destructor.
Definition star_bin.C:359
Vector & set_beta_auto()
Read/write of .
Definition star_bin.C:468
const Scalar & get_loggam() const
Returns the logarithm of the Lorentz factor between the fluid and the co-orbiting observer.
Definition star.h:792
const Metric & get_gtilde() const
Return the conformal 3-metric .
Definition star.h:862
virtual void sauve(FILE *) const
Save in a file.
Definition star_bin.C:489
const Scalar & get_logn_comp() const
Returns the part of the lapse logarithm (gravitational potential at the Newtonian limit) generated pr...
Definition star.h:811
Scalar ssjm1_lnq
Effective source at the previous step for the resolution of the Poisson equation for lnq_auto.
Definition star.h:628
virtual double mass_g() const
Gravitational mass.
bool conf_flat
true if the 3-metric is conformally flat, false for a more general metric.
Definition star.h:681
Base class for stars.
Definition star.h:175
virtual ~Star()
Destructor.
Definition star.C:287
const Scalar & get_nn() const
Returns the lapse function N.
Definition star.h:399
Scalar ener
Total energy density in the fluid frame.
Definition star.h:193
virtual double mass_g() const =0
Gravitational mass.
Scalar logn
Logarithm of the lapse N .
Definition star.h:222
Scalar nn
Lapse function N .
Definition star.h:225
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...
Definition star_global.C:63
const Vector & get_u_euler() const
Returns the fluid 3-velocity with respect to the Eulerian observer.
Definition star.h:385
virtual void del_hydro_euler()
Sets to ETATNONDEF (undefined state) the hydrodynamical quantities relative to the Eulerian observer.
Definition star.C:330
const Scalar & get_gam_euler() const
Returns the Lorentz factor between the fluid and Eulerian observers.
Definition star.h:382
const Tensor & get_stress_euler() const
Returns the spatial part of the stress-energy tensor with respect to the Eulerian observer.
Definition star.h:390
double * p_mass_b
Baryon mass.
Definition star.h:268
const Map & get_mp() const
Returns the mapping.
Definition star.h:355
virtual void equilibrium_spher(double ent_c, double precis=1.e-14, const Tbl *pent_limit=0x0)
Computes a spherical static configuration.
double ray_eq_3pis2() const
Coordinate radius at , [r_unit].
const Eos & get_eos() const
Returns the equation of state.
Definition star.h:361
virtual void hydro_euler()
Computes the hydrodynamical quantities relative to the Eulerian observer from those in the fluid fram...
Definition star.C:577
double * p_ray_eq
Coordinate radius at , .
Definition star.h:242
const Eos & eos
Equation of state of the stellar matter.
Definition star.h:185
const Metric & get_gamma() const
Returns the 3-metric .
Definition star.h:409
const Scalar & get_logn() const
Returns the logarithm of the lapse N.
Definition star.h:396
Scalar nbar
Baryon density in the fluid frame.
Definition star.h:192
Scalar ener_euler
Total energy density in the Eulerian frame.
Definition star.h:198
const Tbl & xi_surf() const
Description of the stellar surface: returns a 2-D Tbl containing the values of the radial coordinate ...
Definition star_global.C:89
const Scalar & get_ener_euler() const
Returns the total energy density with respect to the Eulerian observer.
Definition star.h:376
double * p_ray_eq_pi
Coordinate radius at , .
Definition star.h:248
Itbl * p_l_surf
Description of the stellar surface: 2-D Itbl containing the values of the domain index l on the surfa...
Definition star.h:260
void equation_of_state()
Computes the proper baryon and energy density, as well as pressure from the enthalpy.
Definition star.C:462
double * p_ray_eq_pis2
Coordinate radius at , .
Definition star.h:245
double * p_ray_eq_3pis2
Coordinate radius at , .
Definition star.h:251
virtual void del_deriv() const
Deletes all the derived quantities.
Definition star.C:298
const Scalar & get_ener() const
Returns the proper total energy density.
Definition star.h:370
Scalar gam_euler
Lorentz factor between the fluid and Eulerian observers.
Definition star.h:204
Map & set_mp()
Read/write of the mapping.
Definition star.h:322
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition star.C:417
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition star.C:316
double ray_eq() const
Coordinate radius at , [r_unit].
void set_enthalpy(const Scalar &)
Assignment of the enthalpy field.
Definition star.C:379
Scalar s_euler
Trace of the stress scalar in the Eulerian frame.
Definition star.h:201
double * p_mass_g
Gravitational mass.
Definition star.h:269
Metric gamma
3-metric
Definition star.h:235
Tbl * p_xi_surf
Description of the stellar surface: 2-D Tbl containing the values of the radial coordinate on the su...
Definition star.h:266
virtual void sauve(FILE *) const
Save in a file.
Definition star.C:397
double ray_eq_pis2() const
Coordinate radius at , [r_unit].
Sym_tensor stress_euler
Spatial part of the stress-energy tensor with respect to the Eulerian observer.
Definition star.h:212
const Scalar & get_ent() const
Returns the enthalpy field.
Definition star.h:364
const Scalar & get_s_euler() const
Returns the trace of the stress tensor in the Eulerian frame.
Definition star.h:379
Scalar press
Fluid pressure.
Definition star.h:194
double ray_eq_pi() const
Coordinate radius at , [r_unit].
Scalar ent
Log-enthalpy.
Definition star.h:190
Vector u_euler
Fluid 3-velocity with respect to the Eulerian observer.
Definition star.h:207
const Vector & get_beta() const
Returns the shift vector .
Definition star.h:402
Map & mp
Mapping associated with the star.
Definition star.h:180
int get_nzet() const
Returns the number of domains occupied by the star.
Definition star.h:358
const Scalar & get_press() const
Returns the fluid pressure.
Definition star.h:373
void operator=(const Star &)
Assignment to another Star.
Definition star.C:351
const Scalar & get_nbar() const
Returns the proper baryon density.
Definition star.h:367
int nzet
Number of domains of *mp occupied by the star.
Definition star.h:183
double * p_ray_pole
Coordinate radius at .
Definition star.h:254
Vector beta
Shift vector.
Definition star.h:228
friend ostream & operator<<(ostream &, const Star &)
Display.
Definition star.C:412
double ray_pole() const
Coordinate radius at [r_unit].
virtual double mass_b() const =0
Baryon mass.
Class intended to describe valence-2 symmetric tensors.
Definition sym_tensor.h:223
Basic array class.
Definition tbl.h:161
Tensor handling.
Definition tensor.h:288
Tensor field of valence 1.
Definition vector.h:188
Lorene prototypes.
Definition app_hor.h:64