LORENE
star.C
1 /*
2  * Methods of class Star
3  *
4  * (see file star.h for documentation)
5  */
6 
7 /*
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 char star_C[] = "$Header: /cvsroot/Lorene/C++/Source/Star/star.C,v 1.21 2014/10/13 08:53:37 j_novak Exp $" ;
33 
34 /*
35  * $Id: star.C,v 1.21 2014/10/13 08:53:37 j_novak Exp $
36  * $Log: star.C,v $
37  * Revision 1.21 2014/10/13 08:53:37 j_novak
38  * Lorene classes and functions now belong to the namespace Lorene.
39  *
40  * Revision 1.20 2013/04/20 20:56:15 m_bejger
41  * Fix for three domains in star in Star::equation_of_state from Etoile/etoile.C
42  *
43  * Revision 1.19 2010/02/02 12:45:16 e_gourgoulhon
44  * Improved the display (operator>>)
45  *
46  * Revision 1.18 2010/01/26 16:49:03 e_gourgoulhon
47  * Commented the test on the relativistic character of the EOS: the
48  * relativity parameter is not defined (yet !) in the base class Star.
49  *
50  * Revision 1.17 2007/11/06 16:22:03 j_novak
51  * The data member stress_euler is now a Sym_tensor instead of a Tensor.
52  *
53  * Revision 1.16 2007/06/21 19:53:47 k_taniguchi
54  * Addition of p_ray_eq_3pis2
55  *
56  * Revision 1.15 2005/09/14 12:30:52 f_limousin
57  * Saving of fields lnq and logn in class Star.
58  *
59  * Revision 1.14 2005/09/13 19:38:31 f_limousin
60  * Reintroduction of the resolution of the equations in cartesian coordinates.
61  *
62  * Revision 1.13 2005/02/17 17:29:04 f_limousin
63  * Change the name of some quantities to be consistent with other classes
64  * (for instance nnn is changed to nn, shift to beta, beta to lnq...)
65  *
66  * Revision 1.12 2005/01/05 17:43:03 f_limousin
67  * u_euler is now constructed in the spherical triad to be consistent
68  * with all the others vectors ans tensors.
69  *
70  * Revision 1.11 2004/11/11 16:29:49 j_novak
71  * set_der_0x0 is no longer virtual (to be coherent with Tensor/Scalar classes).
72  *
73  * Revision 1.10 2004/06/22 12:48:08 f_limousin
74  * Change qq, qq_auto and qq_comp to beta, beta_auto and beta_comp.
75  *
76  * Revision 1.9 2004/06/07 16:21:35 f_limousin
77  * Add outputs
78  *
79  * Revision 1.8 2004/04/08 16:32:10 f_limousin
80  * The new variable is ln(Q) instead of Q=psi^2*N. It improves the
81  * convergence of the code.
82  *
83  * Revision 1.7 2004/03/25 10:29:26 j_novak
84  * All LORENE's units are now defined in the namespace Unites (in file unites.h).
85  *
86  * Revision 1.6 2004/03/08 11:48:00 f_limousin
87  * Error in del_deriv() and set_der_0x0() : p_mass_b and p_mass_g were
88  * missing. And so they were never recomputed.
89  *
90  * Revision 1.5 2004/02/27 09:36:46 f_limousin
91  * u_euler is now constructed on a cartesian basis instead
92  * of a spherical one.
93  *
94  * Revision 1.4 2004/02/21 17:05:13 e_gourgoulhon
95  * Method Scalar::point renamed Scalar::val_grid_point.
96  * Method Scalar::set_point renamed Scalar::set_grid_point.
97  *
98  * Revision 1.3 2004/01/20 15:16:58 f_limousin
99  * First version
100  *
101  *
102  * $Header: /cvsroot/Lorene/C++/Source/Star/star.C,v 1.21 2014/10/13 08:53:37 j_novak Exp $
103  *
104  */
105 
106 // Headers C
107 #include "math.h"
108 
109 // Headers Lorene
110 #include "star.h"
111 #include "eos.h"
112 #include "utilitaires.h"
113 #include "param.h"
114 #include "unites.h"
115 
116 
117  //--------------//
118  // Constructors //
119  //--------------//
120 
121 // Standard constructor
122 // --------------------
123 namespace Lorene {
124 Star::Star(Map& mpi, int nzet_i, const Eos& eos_i)
125  : mp(mpi),
126  nzet(nzet_i),
127  eos(eos_i),
128  ent(mpi),
129  nbar(mpi),
130  ener(mpi),
131  press(mpi),
132  ener_euler(mpi),
133  s_euler(mpi),
134  gam_euler(mpi),
135  u_euler(mpi, CON, mp.get_bvect_spher()),
136  stress_euler(mpi, CON, mp.get_bvect_spher()),
137  logn(mpi),
138  nn(mpi),
139  beta(mpi, CON, mp.get_bvect_spher()),
140  lnq(mpi),
141  gamma(mp.flat_met_spher()){
142 
143 
144  // Check of the EOS
145 // const Eos_poly_newt* p_eos_poly_newt =
146 // dynamic_cast<const Eos_poly_newt*>( &eos ) ;
147 //
148 // const Eos_incomp_newt* p_eos_incomp_newt =
149 // dynamic_cast<const Eos_incomp_newt*>( &eos ) ;
150 //
151 //
152 //
153 // if (p_eos_poly_newt != 0x0) {
154 // cout <<
155 // "Star::Star : the EOS Eos_poly_newt must not be employed"
156 // << " for a relativistic star ! " << endl ;
157 // cout << "(Use Eos_poly instead)" << endl ;
158 // abort() ;
159 // }
160 // if (p_eos_incomp_newt != 0x0) {
161 // cout <<
162 // "Star::Star : the EOS Eos_incomp_newt must not be employed"
163 // << " for a relativistic star ! " << endl ;
164 // cout << "(Use Eos_incomp instead)" << endl ;
165 // abort() ;
166 // }
167 //
168  // Pointers of derived quantities initialized to zero :
169  set_der_0x0() ;
170 
171  // All the matter quantities are initialized to zero :
172  nbar = 0 ;
173  ener = 0 ;
174  press = 0 ;
175  ent = 0 ;
176  ener_euler = 0 ;
177  s_euler = 0 ;
178  gam_euler = 1. ;
182 
183  // The metric is initialized to the flat one :
184  Metric flat(mp.flat_met_spher()) ;
185  flat.cov() ;
186  gamma = flat ;
187 
188  logn = 0 ;
189  nn = 1. ;
190  nn.std_spectral_base() ;
191  beta.set_etat_zero() ;
192  lnq = 0 ;
193 
194 }
195 
196 // Copy constructor
197 // ----------------
198 Star::Star(const Star& et)
199  : mp(et.mp),
200  nzet(et.nzet),
201  eos(et.eos),
202  ent(et.ent),
203  nbar(et.nbar),
204  ener(et.ener),
205  press(et.press),
206  ener_euler(et.ener_euler),
207  s_euler(et.s_euler),
208  gam_euler(et.gam_euler),
209  u_euler(et.u_euler),
210  stress_euler(et.stress_euler),
211  logn(et.logn),
212  nn(et.nn),
213  beta(et.beta),
214  lnq(et.lnq),
215  gamma(et.gamma){
216 
217  set_der_0x0() ;
218 
219 }
220 
221 // Constructor from a file
222 // -----------------------
223 Star::Star(Map& mpi, const Eos& eos_i, FILE* fich)
224  : mp(mpi),
225  eos(eos_i),
226  ent(mpi),
227  nbar(mpi),
228  ener(mpi),
229  press(mpi),
230  ener_euler(mpi),
231  s_euler(mpi),
232  gam_euler(mpi),
233  u_euler(mpi, CON, mp.get_bvect_spher()),
234  stress_euler(mpi, CON, mp.get_bvect_spher()),
235  logn(mpi, *(mpi.get_mg()), fich),
236  nn(mpi),
237  beta(mpi, CON, mp.get_bvect_spher()),
238  lnq(mpi, *(mpi.get_mg()), fich),
239  gamma(mpi.flat_met_spher()){
240 
241  // Star parameters
242  // -----------------
243 
244  // nzet is read in the file:
245  int xx ;
246  fread_be(&xx, sizeof(int), 1, fich) ;
247  nzet = xx ;
248 
249  // Equation of state
250  // -----------------
251 
252  // Read of the saved EOS
253  Eos* p_eos_file = Eos::eos_from_file(fich) ;
254 
255  // Comparison with the assigned EOS:
256  if (eos != *p_eos_file) {
257  cout <<
258  "Star::Star(const Map&, const Eos&, FILE*) : the EOS given in "
259  << endl <<
260  " argument and that read in the file are different !" << endl ;
261  abort() ;
262  }
263 
264  // p_eos_file is no longer required (it was used only for checking the
265  // EOS compatibility)
266  delete p_eos_file ;
267 
268  // Read of the saved fields:
269  // ------------------------
270  Scalar ent_file(mp, *(mp.get_mg()), fich) ;
271  ent = ent_file ;
274  nn = 1 ;
275  beta.set_etat_zero() ;
276 
277  // Pointers of derived quantities initialized to zero
278  // --------------------------------------------------
279  set_der_0x0() ;
280 
281 }
282 
283  //------------//
284  // Destructor //
285  //------------//
286 
288 
289  del_deriv() ;
290 
291 }
292 
293 
294  //----------------------------------//
295  // Management of derived quantities //
296  //----------------------------------//
297 
298 void Star::del_deriv() const {
299 
300  if (p_mass_b != 0x0) delete p_mass_b ;
301  if (p_mass_g != 0x0) delete p_mass_g ;
302  if (p_ray_eq != 0x0) delete p_ray_eq ;
303  if (p_ray_eq_pis2 != 0x0) delete p_ray_eq_pis2 ;
304  if (p_ray_eq_pi != 0x0) delete p_ray_eq_pi ;
305  if (p_ray_eq_3pis2 != 0x0) delete p_ray_eq_3pis2 ;
306  if (p_ray_pole != 0x0) delete p_ray_pole ;
307  if (p_l_surf != 0x0) delete p_l_surf ;
308  if (p_xi_surf != 0x0) delete p_xi_surf ;
309 
310  Star::set_der_0x0() ;
311 }
312 
313 
314 
315 
316 void Star::set_der_0x0() const {
317 
318  p_mass_b = 0x0 ;
319  p_mass_g = 0x0 ;
320  p_ray_eq = 0x0 ;
321  p_ray_eq_pis2 = 0x0 ;
322  p_ray_eq_pi = 0x0 ;
323  p_ray_eq_3pis2 = 0x0 ;
324  p_ray_pole = 0x0 ;
325  p_l_surf = 0x0 ;
326  p_xi_surf = 0x0 ;
327 
328 }
329 
331 
337 
338  del_deriv() ;
339 
340 }
341 
342 
343 
344 
345  //--------------//
346  // Assignment //
347  //--------------//
348 
349 // Assignment to another Star
350 // ----------------------------
351 void Star::operator=(const Star& et) {
352 
353  assert( &(et.mp) == &mp ) ; // Same mapping
354  assert( &(et.eos) == &eos ) ; // Same EOS
355 
356  nzet = et.nzet ;
357  ent = et.ent ;
358  nbar = et.nbar ;
359  ener = et.ener ;
360  press = et.press ;
361  ener_euler = et.ener_euler ;
362  s_euler = et.s_euler ;
363  gam_euler = et.gam_euler ;
364  u_euler = et.u_euler ;
366  logn = et.logn ;
367  nn = et.nn ;
368  beta = et.beta ;
369  lnq = et.lnq ;
370  gamma = et.gamma ;
371 
372  del_deriv() ; // Deletes all derived quantities
373 
374 }
375 
376 // Assignment of the enthalpy field
377 // --------------------------------
378 
379 void Star::set_enthalpy(const Scalar& ent_i) {
380 
381  ent = ent_i ;
382 
383  // Update of (nbar, ener, press) :
384  equation_of_state() ;
385 
386  // The derived quantities are obsolete:
387  del_deriv() ;
388 
389 }
390 
391  //--------------//
392  // Outputs //
393  //--------------//
394 
395 // Save in a file
396 // --------------
397 void Star::sauve(FILE* fich) const {
398 
399  logn.sauve(fich) ;
400  lnq.sauve(fich) ;
401 
402  int xx = nzet ;
403  fwrite_be(&xx, sizeof(int), 1, fich) ;
404 
405  eos.sauve(fich) ;
406  ent.sauve(fich) ;
407 }
408 
409 // Printing
410 // --------
411 
412 ostream& operator<<(ostream& ost, const Star& et) {
413  et >> ost ;
414  return ost ;
415 }
416 
417 ostream& Star::operator>>(ostream& ost) const {
418 
419  using namespace Unites ;
420 
421  ost << endl ;
422 
423  ost << "Number of domains occupied by the star : " << nzet << endl ;
424 
425  ost << "Equation of state : " << endl ;
426  ost << eos << endl ;
427 
428  ost << endl << "Central enthalpy : " << ent.val_grid_point(0,0,0,0) << " c^2" << endl ;
429  ost << "Central proper baryon density : " << nbar.val_grid_point(0,0,0,0)
430  << " x 0.1 fm^-3" << endl ;
431  ost << "Central proper energy density : " << ener.val_grid_point(0,0,0,0)
432  << " rho_nuc c^2" << endl ;
433  ost << "Central pressure : " << press.val_grid_point(0,0,0,0)
434  << " rho_nuc c^2" << endl ;
435 
436  ost << endl ;
437  ost << "Central lapse N : " << nn.val_grid_point(0,0,0,0) << endl ;
438 // ost << "Central value of lnq : " << lnq.val_grid_point(0,0,0,0) << endl ;
439 
440  ost << endl
441  << "Coordinate equatorial radius (phi=0) a1 = "
442  << ray_eq()/km << " km" << endl ;
443  ost << "Coordinate equatorial radius (phi=pi/2) a2 = "
444  << ray_eq_pis2()/km << " km" << endl ;
445  ost << "Coordinate equatorial radius (phi=pi): "
446  << ray_eq_pi()/km << " km" << endl ;
447  ost << "Coordinate polar radius a3 = "
448  << ray_pole()/km << " km" << endl ;
449  ost << "Axis ratio a2/a1 = " << ray_eq_pis2() / ray_eq()
450  << " a3/a1 = " << ray_pole() / ray_eq() << endl ;
451  ost << endl << "Baryon mass : " << mass_b() / msol << " M_sol" << endl ;
452  ost << "Gravitational mass : " << mass_g() / msol << " M_sol" << endl ;
453 
454 
455  return ost ;
456 }
457 
458  //-----------------------------------------//
459  // Computation of hydro quantities //
460  //-----------------------------------------//
461 
463 
464  Scalar ent_eos = ent ;
465 
466 
467  // Slight rescale of the enthalpy field in case of 2 domains inside the
468  // star
469 
470 
471  double epsilon = 1.e-12 ;
472 
473  const Mg3d* mg = mp.get_mg() ;
474  int nz = mg->get_nzone() ;
475 
476  Mtbl xi(mg) ;
477  xi.set_etat_qcq() ;
478  for (int l=0; l<nz; l++) {
479  xi.t[l]->set_etat_qcq() ;
480  for (int k=0; k<mg->get_np(l); k++) {
481  for (int j=0; j<mg->get_nt(l); j++) {
482  for (int i=0; i<mg->get_nr(l); i++) {
483  xi.set(l,k,j,i) =
484  mg->get_grille3d(l)->x[i] ;
485  }
486  }
487  }
488 
489  }
490 
491  Scalar fact_ent(mp) ;
492  fact_ent.allocate_all() ;
493 
494  fact_ent.set_domain(0) = 1 + epsilon * xi(0) * xi(0) ;
495  fact_ent.set_domain(1) = 1 - 0.25 * epsilon * (xi(1) - 1) * (xi(1) - 1) ;
496 
497  for (int l=nzet; l<nz; l++) {
498  fact_ent.set_domain(l) = 1 ;
499  }
500 
501  if (nzet > 1) {
502 
503  if(nzet == 3) {
504  fact_ent.set_domain(1) = 1 - 0.5 * epsilon * (xi(1) - 0.5) * (xi(1) - 0.5) ;
505  fact_ent.set_domain(2) = 1 - 0.25 * epsilon * (xi(2) - 1) * (xi(2) - 1) ;
506  }
507 
508  if (nzet > 3) {
509 
510  cout << "Star::equation_of_state: not ready yet for nzet > 3 !"
511  << endl ;
512  }
513 
514  ent_eos = fact_ent * ent_eos ;
515  ent_eos.std_spectral_base() ;
516  }
517 
518 
519 
520 
521 
522  // Call to the EOS (the EOS is called domain by domain in order to
523  // allow for the use of MEos)
524 
525  Scalar tempo(mp) ;
526 
527  nbar.set_etat_qcq() ;
528  nbar = 0 ;
529  for (int l=0; l<nzet; l++) {
530 
531  Param par ; // Paramater for multi-domain equation of state
532  par.add_int(l) ;
533 
534  tempo = eos.nbar_ent(ent_eos, 1, l, &par) ;
535 
536  nbar = nbar + tempo ;
537 
538  }
539 
540  ener.set_etat_qcq() ;
541  ener = 0 ;
542  for (int l=0; l<nzet; l++) {
543 
544  Param par ; // Paramater for multi-domain equation of state
545  par.add_int(l) ;
546 
547  tempo = eos.ener_ent(ent_eos, 1, l, &par) ;
548 
549  ener = ener + tempo ;
550 
551  }
552 
553  press.set_etat_qcq() ;
554  press = 0 ;
555  for (int l=0; l<nzet; l++) {
556 
557  Param par ; // Paramater for multi-domain equation of state
558  par.add_int(l) ;
559 
560  tempo = eos.press_ent(ent_eos, 1, l, &par) ;
561 
562  press = press + tempo ;
563 
564  }
565 
566 
567  // Set the bases for spectral expansion
571 
572  // The derived quantities are obsolete
573  del_deriv() ;
574 
575 }
576 
578 
579  cout <<
580  "Star::hydro_euler : hydro_euler must be called via a derived class"
581  << endl << " of Star !" << endl ;
582 
583  abort() ;
584 
585 }
586 }
Equation of state base class.
Definition: eos.h:190
Cmp nbar_ent(const Cmp &ent, int nzet, int l_min=0, const Param *par=0x0) const
Computes the baryon density field from the log-enthalpy field and extra parameters.
Definition: eos.C:338
static Eos * eos_from_file(FILE *)
Construction of an EOS from a binary file.
Cmp ener_ent(const Cmp &ent, int nzet, int l_min=0, const Param *par=0x0) const
Computes the total energy density from the log-enthalpy and extra parameters.
Definition: eos.C:363
virtual void sauve(FILE *) const
Save in a file.
Definition: eos.C:179
Cmp press_ent(const Cmp &ent, int nzet, int l_min=0, const Param *par=0x0) const
Computes the pressure from the log-enthalpy and extra parameters.
Definition: eos.C:385
double * x
Array of values of at the nr collocation points.
Definition: grilles.h:209
Base class for coordinate mappings.
Definition: map.h:670
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:765
const Metric_flat & flat_met_spher() const
Returns the flat metric associated with the spherical coordinates and with components expressed in th...
Definition: map.C:321
Metric for tensor calculation.
Definition: metric.h:90
virtual const Sym_tensor & cov() const
Read-only access to the covariant representation.
Definition: metric.C:280
Multi-domain grid.
Definition: grilles.h:273
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
Definition: grilles.h:462
const Grille3d * get_grille3d(int l) const
Returns a pointer on the 3D mono-grid for domain no. l.
Definition: grilles.h:500
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
Definition: grilles.h:457
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:448
int get_nr(int l) const
Returns the number of points in the radial direction ( ) in domain no. l.
Definition: grilles.h:452
Multi-domain array.
Definition: mtbl.h:118
Tbl ** t
Array (size nzone ) of pointers on the Tbl 's.
Definition: mtbl.h:132
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: mtbl.C:299
Tbl & set(int l)
Read/write of the Tbl in a given domain.
Definition: mtbl.h:221
Parameter storage.
Definition: param.h:125
void add_int(const int &n, int position=0)
Adds the address of a new int to the list.
Definition: param.C:246
Tensor field of valence 0 (or component of a tensorial field).
Definition: scalar.h:387
virtual void sauve(FILE *) const
Save in a file.
Definition: scalar.C:686
virtual void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: scalar.C:353
Tbl & set_domain(int l)
Read/write of the value in a given domain.
Definition: scalar.h:615
virtual void std_spectral_base()
Sets the spectral bases of the Valeur va to the standard ones for a scalar field.
Definition: scalar.C:784
double val_grid_point(int l, int k, int j, int i) const
Returns the value of the field at a specified grid point.
Definition: scalar.h:637
virtual void allocate_all()
Sets the logical state to ETATQCQ (ordinary state) and performs the memory allocation of all the elem...
Definition: scalar.C:365
virtual void set_etat_nondef()
Sets the logical state to ETATNONDEF (undefined).
Definition: scalar.C:344
Base class for stars.
Definition: star.h:175
virtual ~Star()
Destructor.
Definition: star.C:287
Star(Map &mp_i, int nzet_i, const Eos &eos_i)
Standard constructor.
Definition: star.C:124
Scalar ener
Total energy density in the fluid frame.
Definition: star.h:193
virtual double mass_g() const =0
Gravitational mass.
Definition: star_global.C:327
Scalar logn
Logarithm of the lapse N .
Definition: star.h:222
Scalar nn
Lapse function N .
Definition: star.h:225
virtual void del_hydro_euler()
Sets to ETATNONDEF (undefined state) the hydrodynamical quantities relative to the Eulerian observer.
Definition: star.C:330
double * p_mass_b
Baryon mass.
Definition: star.h:268
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
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
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
Scalar gam_euler
Lorentz factor between the fluid and Eulerian observers.
Definition: star.h:204
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].
Definition: star_global.C:108
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].
Definition: star_global.C:138
Sym_tensor stress_euler
Spatial part of the stress-energy tensor with respect to the Eulerian observer.
Definition: star.h:212
Scalar press
Fluid pressure.
Definition: star.h:194
double ray_eq_pi() const
Coordinate radius at , [r_unit].
Definition: star_global.C:186
Scalar ent
Log-enthalpy.
Definition: star.h:190
Vector u_euler
Fluid 3-velocity with respect to the Eulerian observer.
Definition: star.h:207
Map & mp
Mapping associated with the star.
Definition: star.h:180
void operator=(const Star &)
Assignment to another Star.
Definition: star.C:351
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
double ray_pole() const
Coordinate radius at [r_unit].
Definition: star_global.C:278
virtual double mass_b() const =0
Baryon mass.
Definition: star_global.C:307
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: tbl.C:361
int fread_be(int *aa, int size, int nb, FILE *fich)
Reads integer(s) from a binary file according to the big endian convention.
Definition: fread_be.C:69
int fwrite_be(const int *aa, int size, int nb, FILE *fich)
Writes integer(s) into a binary file according to the big endian convention.
Definition: fwrite_be.C:70
virtual void set_etat_nondef()
Sets the logical state of all components to ETATNONDEF (undefined state).
Definition: tensor.C:489
virtual void set_etat_zero()
Sets the logical state of all components to ETATZERO (zero state).
Definition: tensor.C:497
Lorene prototypes.
Definition: app_hor.h:64
Standard units of space, time and mass.