LORENE
etoile.C
1/*
2 * Methods of class Etoile
3 *
4 * (see file etoile.h for documentation)
5 */
6
7/*
8 * Copyright (c) 2000-2001 Eric Gourgoulhon
9 * Copyright (c) 2000-2001 Keisuke Taniguchi
10 *
11 * This file is part of LORENE.
12 *
13 * LORENE is free software; you can redistribute it and/or modify
14 * it under the terms of the GNU General Public License as published by
15 * the Free Software Foundation; either version 2 of the License, or
16 * (at your option) any later version.
17 *
18 * LORENE is distributed in the hope that it will be useful,
19 * but WITHOUT ANY WARRANTY; without even the implied warranty of
20 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
21 * GNU General Public License for more details.
22 *
23 * You should have received a copy of the GNU General Public License
24 * along with LORENE; if not, write to the Free Software
25 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
26 *
27 */
28
29
30char etoile_C[] = "$Header: /cvsroot/Lorene/C++/Source/Etoile/etoile.C,v 1.10 2014/10/13 08:52:58 j_novak Exp $" ;
31
32/*
33 * $Id: etoile.C,v 1.10 2014/10/13 08:52:58 j_novak Exp $
34 * $Log: etoile.C,v $
35 * Revision 1.10 2014/10/13 08:52:58 j_novak
36 * Lorene classes and functions now belong to the namespace Lorene.
37 *
38 * Revision 1.9 2012/08/12 17:48:35 p_cerda
39 * Magnetstar: New classes for magnetstar. Allowing for non-equatorial symmetry in Etoile et al. Adding B_phi in Et_rot_mag.
40 *
41 * Revision 1.8 2005/01/18 22:36:50 k_taniguchi
42 * Delete a pointer for ray_eq(int kk).
43 *
44 * Revision 1.7 2005/01/18 20:35:05 k_taniguchi
45 * Addition of ray_eq(int kk).
46 *
47 * Revision 1.6 2005/01/17 20:40:25 k_taniguchi
48 * Addition of ray_eq_3pis2().
49 *
50 * Revision 1.5 2004/03/25 10:29:06 j_novak
51 * All LORENE's units are now defined in the namespace Unites (in file unites.h).
52 *
53 * Revision 1.4 2003/10/13 15:23:56 f_limousin
54 * *** empty log message ***
55 *
56 * Revision 1.3 2002/04/09 14:32:15 e_gourgoulhon
57 * 1/ Added extra parameters in EOS computational functions (argument par)
58 * 2/ New class MEos for multi-domain EOS
59 *
60 * Revision 1.2 2001/12/04 21:27:53 e_gourgoulhon
61 *
62 * All writing/reading to a binary file are now performed according to
63 * the big endian convention, whatever the system is big endian or
64 * small endian, thanks to the functions fwrite_be and fread_be
65 *
66 * Revision 1.1.1.1 2001/11/20 15:19:28 e_gourgoulhon
67 * LORENE
68 *
69 * Revision 2.14 2000/11/24 13:27:44 eric
70 * Dans eqution_of_state(): changement leger de ent dans le cas ou l'on a
71 * deux domaine avant d'appeler l'EOS.
72 *
73 * Revision 2.13 2000/09/25 12:22:02 keisuke
74 * *** empty log message ***
75 *
76 * Revision 2.12 2000/09/22 15:50:58 keisuke
77 * Ajout du membre d_logn_auto_div.
78 *
79 * Revision 2.11 2000/09/07 14:34:09 keisuke
80 * Ajout du membre logn_auto_regu.
81 *
82 * Revision 2.10 2000/08/31 15:36:54 eric
83 * Bases spectrales standards pour nnn, a_car et gam_euler dans le
84 * constructeur (initialisation a la metrique plate).
85 *
86 * Revision 2.9 2000/08/29 11:37:49 eric
87 * Ajout des membres k_div et logn_auto_div.
88 *
89 * Revision 2.8 2000/07/21 12:01:11 eric
90 * Modif dans Etoile::del_deriv() :
91 * appel de Etoile::set_der_0x0() et non de la fonction virtuelle set_der_0x0().
92 *
93 * Revision 2.7 2000/03/21 12:39:34 eric
94 * Le constructeur standard teste la compatibilite de l'EOS avec le
95 * caractere relativiste de l'etoile.
96 *
97 * Revision 2.6 2000/02/21 14:32:40 eric
98 * gam_euler est initialise a 1 dans le constructeur standard.
99 * Suppression de l'appel a del_hydro_euler dans equation_of_state().
100 *
101 * Revision 2.5 2000/02/09 19:30:47 eric
102 * La triade de decomposition doit desormais figurer en argument des
103 * constructeurs de Tenseur.
104 *
105 * Revision 2.4 2000/02/02 09:23:34 eric
106 * Affichage de la masse.
107 *
108 * Revision 2.3 2000/01/28 17:18:10 eric
109 * Modifs noms des quantites globales.
110 * Affichage.
111 *
112 * Revision 2.2 2000/01/24 17:13:36 eric
113 * Le mapping mp n'est plus constant.
114 *
115 * Revision 2.1 2000/01/24 13:37:22 eric
116 * *** empty log message ***
117 *
118 * Revision 2.0 2000/01/20 17:04:45 eric
119 * *** empty log message ***
120 *
121 *
122 * $Header: /cvsroot/Lorene/C++/Source/Etoile/etoile.C,v 1.10 2014/10/13 08:52:58 j_novak Exp $
123 *
124 */
125
126// Headers C
127#include "math.h"
128
129// Headers Lorene
130#include "etoile.h"
131#include "eos.h"
132#include "utilitaires.h"
133#include "param.h"
134#include "unites.h"
135
136 //--------------//
137 // Constructors //
138 //--------------//
139
140// Standard constructor
141// --------------------
142namespace Lorene {
143Etoile::Etoile(Map& mpi, int nzet_i, bool relat, const Eos& eos_i)
144 : mp(mpi),
145 nzet(nzet_i),
146 relativistic(relat),
147 k_div(0),
148 eos(eos_i),
149 ent(mpi),
150 nbar(mpi),
151 ener(mpi),
152 press(mpi),
153 ener_euler(mpi),
154 s_euler(mpi),
155 gam_euler(mpi),
156 u_euler(mpi, 1, CON, mp.get_bvect_cart()),
157 logn_auto(mpi),
158 logn_auto_regu(mpi),
159 logn_auto_div(mpi),
160 d_logn_auto_div(mpi, 1, COV, mp.get_bvect_spher()),
161 beta_auto(mpi),
162 nnn(mpi),
163 shift(mpi, 1, CON, mp.get_bvect_cart()),
164 a_car(mpi) {
165
166 // Check of the EOS
167 const Eos_poly* p_eos_poly = dynamic_cast<const Eos_poly*>( &eos ) ;
168 const Eos_poly_newt* p_eos_poly_newt =
169 dynamic_cast<const Eos_poly_newt*>( &eos ) ;
170 const Eos_incomp* p_eos_incomp = dynamic_cast<const Eos_incomp*>( &eos ) ;
171 const Eos_incomp_newt* p_eos_incomp_newt =
172 dynamic_cast<const Eos_incomp_newt*>( &eos ) ;
173
174 if (relativistic) {
175
176 if (p_eos_poly_newt != 0x0) {
177 cout <<
178 "Etoile::Etoile : the EOS Eos_poly_newt must not be employed"
179 << " for a relativistic star ! " << endl ;
180 cout << "(Use Eos_poly instead)" << endl ;
181 abort() ;
182 }
183 if (p_eos_incomp_newt != 0x0) {
184 cout <<
185 "Etoile::Etoile : the EOS Eos_incomp_newt must not be employed"
186 << " for a relativistic star ! " << endl ;
187 cout << "(Use Eos_incomp instead)" << endl ;
188 abort() ;
189 }
190
191 }
192 else{
193
194 if ( (p_eos_poly != 0x0) && (p_eos_poly_newt == 0x0) ) {
195 cout <<
196 "Etoile::Etoile : the EOS Eos_poly must not be employed"
197 << " for a Newtonian star ! " << endl ;
198 cout << "(Use Eos_poly_newt instead)" << endl ;
199 abort() ;
200 }
201 if ( (p_eos_incomp != 0x0) && (p_eos_incomp_newt == 0x0) ) {
202 cout <<
203 "Etoile::Etoile : the EOS Eos_incomp must not be employed"
204 << " for a relativistic star ! " << endl ;
205 cout << "(Use Eos_incomp_newt instead)" << endl ;
206 abort() ;
207 }
208
209 }
210
211
212 // Parameter 1/c^2
213 unsurc2 = relativistic ? double(1) : double(0) ;
214
215 // Pointers of derived quantities initialized to zero :
216 set_der_0x0() ;
217
218 // All the matter quantities are initialized to zero :
219 nbar = 0 ;
220 ener = 0 ;
221 press = 0 ;
222 ent = 0 ;
223 ener_euler = 0 ;
224 s_euler = 0 ;
225 gam_euler = 1 ;
227 u_euler = 0 ;
228
229 // The metric is initialized to the flat one :
230 logn_auto = 0 ;
231 logn_auto_regu = 0 ;
232 logn_auto_div = 0 ;
233 d_logn_auto_div = 0 ;
234 beta_auto = 0 ;
235 nnn = 1 ;
236 nnn.set_std_base() ;
237 shift = 0 ;
238 a_car = 1 ;
240
241}
242
243// Copy constructor
244// ----------------
246 : mp(et.mp),
247 nzet(et.nzet),
248 relativistic(et.relativistic),
249 unsurc2(et.unsurc2),
250 k_div(et.k_div),
251 eos(et.eos),
252 ent(et.ent),
253 nbar(et.nbar),
254 ener(et.ener),
255 press(et.press),
256 ener_euler(et.ener_euler),
257 s_euler(et.s_euler),
258 gam_euler(et.gam_euler),
259 u_euler(et.u_euler),
260 logn_auto(et.logn_auto),
261 logn_auto_regu(et.logn_auto_regu),
262 logn_auto_div(et.logn_auto_div),
263 d_logn_auto_div(et.d_logn_auto_div),
264 beta_auto(et.beta_auto),
265 nnn(et.nnn),
266 shift(et.shift),
267 a_car(et.a_car) {
268
269 set_der_0x0() ;
270
271}
272
273// Constructor from a file
274// -----------------------
275Etoile::Etoile(Map& mpi, const Eos& eos_i, FILE* fich)
276 : mp(mpi),
277 eos(eos_i),
278 ent(mpi),
279 nbar(mpi),
280 ener(mpi),
281 press(mpi),
282 ener_euler(mpi),
283 s_euler(mpi),
284 gam_euler(mpi),
285 u_euler(mpi, 1, CON, mp.get_bvect_cart()),
286 logn_auto(mpi),
287 logn_auto_regu(mpi),
288 logn_auto_div(mpi),
289 d_logn_auto_div(mpi, 1, COV, mp.get_bvect_spher()),
290 beta_auto(mpi),
291 nnn(mpi),
292 shift(mpi, 1, CON, mp.get_bvect_cart()),
293 a_car(mpi) {
294
295 // Etoile parameters
296 // -----------------
297
298 // nzet and relativistic are read in the file:
299 int xx ;
300 fread_be(&xx, sizeof(int), 1, fich) ;
301 k_div = xx / 1000 ; // integer part
302 nzet = xx - k_div * 1000 ;
303
304 fread(&relativistic, sizeof(bool), 1, fich) ;
305
306 // Parameter 1/c^2 is deduced from relativistic:
307 unsurc2 = relativistic ? double(1) : double(0) ;
308
309
310 // Equation of state
311 // -----------------
312
313 // Read of the saved EOS
314 Eos* p_eos_file = Eos::eos_from_file(fich) ;
315
316 // Comparison with the assigned EOS:
317 if (eos != *p_eos_file) {
318 cout <<
319 "Etoile::Etoile(const Map&, const Eos&, FILE*) : the EOS given in "
320 << endl <<
321 " argument and that read in the file are different !" << endl ;
322 abort() ;
323 }
324
325 // p_eos_file is no longer required (it was used only for checking the
326 // EOS compatibility)
327 delete p_eos_file ;
328
329 // Read of the saved fields:
330 // ------------------------
331 Tenseur ent_file(mp, fich) ;
332 ent = ent_file ;
333
334 Tenseur logn_auto_file(mp, fich) ;
335 logn_auto = logn_auto_file ;
336
337 Tenseur beta_auto_file(mp, fich) ;
338 beta_auto = beta_auto_file ;
339
340 if (k_div == 0) {
341 logn_auto_div = 0 ;
342 d_logn_auto_div = 0 ;
343 }
344 else {
345
346 Tenseur logn_auto_div_file(mp, fich) ;
347 logn_auto_div = logn_auto_div_file ;
348
349 Tenseur d_logn_auto_div_file(mp, mp.get_bvect_spher(), fich) ;
350 d_logn_auto_div = d_logn_auto_div_file ;
351 }
352
354
355 shift = 0 ;
356
357 // Pointers of derived quantities initialized to zero
358 // --------------------------------------------------
359 set_der_0x0() ;
360
361}
362
363 //------------//
364 // Destructor //
365 //------------//
366
368
369 del_deriv() ;
370
371}
372
373
374 //----------------------------------//
375 // Management of derived quantities //
376 //----------------------------------//
377
378void Etoile::del_deriv() const {
379
380 if (p_mass_b != 0x0) delete p_mass_b ;
381 if (p_mass_g != 0x0) delete p_mass_g ;
382 if (p_ray_eq != 0x0) delete p_ray_eq ;
383 if (p_ray_eq_pis2 != 0x0) delete p_ray_eq_pis2 ;
384 if (p_ray_eq_pi != 0x0) delete p_ray_eq_pi ;
385 if (p_ray_eq_3pis2 != 0x0) delete p_ray_eq_3pis2 ;
386 if (p_ray_pole != 0x0) delete p_ray_pole ;
387 if (p_l_surf != 0x0) delete p_l_surf ;
388 if (p_xi_surf != 0x0) delete p_xi_surf ;
389
391}
392
393
394
395
397
398 p_mass_b = 0x0 ;
399 p_mass_g = 0x0 ;
400 p_ray_eq = 0x0 ;
401 p_ray_eq_pis2 = 0x0 ;
402 p_ray_eq_pi = 0x0 ;
403 p_ray_eq_3pis2 = 0x0 ;
404 p_ray_pole = 0x0 ;
405 p_l_surf = 0x0 ;
406 p_xi_surf = 0x0 ;
407
408}
409
420
421
422
423
424 //--------------//
425 // Assignment //
426 //--------------//
427
428// Assignment to another Etoile
429// ----------------------------
430void Etoile::operator=(const Etoile& et) {
431
432 assert( &(et.mp) == &mp ) ; // Same mapping
433 assert( &(et.eos) == &eos ) ; // Same EOS
434
435 nzet = et.nzet ;
437 k_div = et.k_div ;
438 unsurc2 = et.unsurc2 ;
439
440 ent = et.ent ;
441 nbar = et.nbar ;
442 ener = et.ener ;
443 press = et.press ;
445 s_euler = et.s_euler ;
446 gam_euler = et.gam_euler ;
447 u_euler = et.u_euler ;
448 logn_auto = et.logn_auto ;
452 beta_auto = et.beta_auto ;
453 nnn = et.nnn ;
454 shift = et.shift ;
455 a_car = et.a_car ;
456
457
458 del_deriv() ; // Deletes all derived quantities
459
460}
461
462// Assignment of the enthalpy field
463// --------------------------------
464
465void Etoile::set_enthalpy(const Cmp& ent_i) {
466
467 ent = ent_i ;
468
469 // Update of (nbar, ener, press) :
471
472 // The derived quantities are obsolete:
473 del_deriv() ;
474
475}
476
477 //--------------//
478 // Outputs //
479 //--------------//
480
481// Save in a file
482// --------------
483void Etoile::sauve(FILE* fich) const {
484
485 int xx = nzet + k_div * 1000 ;
486 fwrite_be(&xx, sizeof(int), 1, fich) ;
487
488 fwrite(&relativistic, sizeof(bool), 1, fich) ;
489
490 eos.sauve(fich) ;
491
492 ent.sauve(fich) ;
493 logn_auto.sauve(fich) ;
494 beta_auto.sauve(fich) ;
495
496 if (k_div != 0) {
497 logn_auto_div.sauve(fich) ;
498 d_logn_auto_div.sauve(fich) ;
499 }
500
501}
502
503// Printing
504// --------
505
506ostream& operator<<(ostream& ost, const Etoile& et) {
507 et >> ost ;
508 return ost ;
509}
510
511ostream& Etoile::operator>>(ostream& ost) const {
512
513 using namespace Unites ;
514
515 ost << endl ;
516 if (relativistic) {
517 ost << "Relativistic star" << endl ;
518 ost << "-----------------" << endl ;
519 }
520 else {
521 ost << "Newtonian star" << endl ;
522 ost << "--------------" << endl ;
523 }
524
525 ost << "Number of domains occupied by the star : " << nzet << endl ;
526
527 ost << "Equation of state : " << endl ;
528 ost << eos << endl ;
529
530 ost << endl << "Central enthalpy : " << ent()(0,0,0,0) << " c^2" << endl ;
531 ost << "Central proper baryon density : " << nbar()(0,0,0,0)
532 << " x 0.1 fm^-3" << endl ;
533 ost << "Central proper energy density : " << ener()(0,0,0,0)
534 << " rho_nuc c^2" << endl ;
535 ost << "Central pressure : " << press()(0,0,0,0)
536 << " rho_nuc c^2" << endl ;
537
538 ost << endl
539 << "Regularization index of the gravitational potential : k_div = "
540 << k_div << endl ;
541 ost << "Central lapse N : " << nnn()(0,0,0,0) << endl ;
542 ost << "Central value of A^2 : " << a_car()(0,0,0,0) << endl ;
543
544 ost << endl
545 << "Coordinate equatorial radius (phi=0) a1 = "
546 << ray_eq()/km << " km" << endl ;
547 ost << "Coordinate equatorial radius (phi=pi/2) a2 = "
548 << ray_eq_pis2()/km << " km" << endl ;
549 ost << "Coordinate equatorial radius (phi=pi): "
550 << ray_eq_pi()/km << " km" << endl ;
551 ost << "Coordinate polar radius a3 = "
552 << ray_pole()/km << " km" << endl ;
553 ost << "Axis ratio a2/a1 = " << ray_eq_pis2() / ray_eq()
554 << " a3/a1 = " << ray_pole() / ray_eq() << endl ;
555
556 ost << endl << "Baryon mass : " << mass_b() / msol << " M_sol" << endl ;
557 ost << "Gravitational mass : " << mass_g() / msol << " M_sol" << endl ;
558
559 return ost ;
560}
561
562 //-----------------------------------------//
563 // Computation of hydro quantities //
564 //-----------------------------------------//
565
567
568 Cmp ent_eos = ent() ;
569
570
571 // Slight rescale of the enthalpy field in case of 2 domains inside the
572 // star
573
574
575 double epsilon = 1.e-12 ;
576
577 const Mg3d* mg = mp.get_mg() ;
578 int nz = mg->get_nzone() ;
579
580 Mtbl xi(mg) ;
581 xi.set_etat_qcq() ;
582 for (int l=0; l<nz; l++) {
583 xi.t[l]->set_etat_qcq() ;
584 for (int k=0; k<mg->get_np(l); k++) {
585 for (int j=0; j<mg->get_nt(l); j++) {
586 for (int i=0; i<mg->get_nr(l); i++) {
587 xi.set(l,k,j,i) =
588 mg->get_grille3d(l)->x[i] ;
589 }
590 }
591 }
592
593 }
594
595 Cmp fact_ent(mp) ;
596 fact_ent.allocate_all() ;
597
598 fact_ent.set(0) = 1 + epsilon * xi(0) * xi(0) ;
599 fact_ent.set(1) = 1 - 0.25 * epsilon * (xi(1) - 1) * (xi(1) - 1) ;
600
601 for (int l=nzet; l<nz; l++) {
602 fact_ent.set(l) = 1 ;
603 }
604
605 if (nzet > 1) {
606
607 if(nzet == 3) {
608 fact_ent.set(1) = 1 - 0.5 * epsilon * (xi(1) - 0.5) * (xi(1) - 0.5) ;
609 fact_ent.set(2) = 1 - 0.25 * epsilon * (xi(2) - 1) * (xi(2) - 1) ;
610 }
611
612 if (nzet > 3) {
613
614 cout << "Etoile::equation_of_state: not ready yet for nzet > 3 !"
615 << endl ;
616 }
617
618 ent_eos = fact_ent * ent_eos ;
619 ent_eos.std_base_scal() ;
620 }
621
622
623
624
625
626 // Call to the EOS (the EOS is called domain by domain in order to
627 // allow for the use of MEos)
628
629 Cmp tempo(mp) ;
630
632 nbar.set() = 0 ;
633 for (int l=0; l<nzet; l++) {
634
635 Param par ; // Paramater for multi-domain equation of state
636 par.add_int(l) ;
637
638 tempo = eos.nbar_ent(ent_eos, 1, l, &par) ;
639
640 nbar = nbar() + tempo ;
641
642 }
643
645 ener.set() = 0 ;
646 for (int l=0; l<nzet; l++) {
647
648 Param par ; // Paramater for multi-domain equation of state
649 par.add_int(l) ;
650
651 tempo = eos.ener_ent(ent_eos, 1, l, &par) ;
652
653 ener = ener() + tempo ;
654
655 }
656
658 press.set() = 0 ;
659 for (int l=0; l<nzet; l++) {
660
661 Param par ; // Paramater for multi-domain equation of state
662 par.add_int(l) ;
663
664 tempo = eos.press_ent(ent_eos, 1, l, &par) ;
665
666 press = press() + tempo ;
667
668 }
669
670
671 // Set the bases for spectral expansion
672 nbar.set_std_base() ;
673 ener.set_std_base() ;
675
676 // The Eulerian quantities are obsolete
677 //## del_hydro_euler() ;
678
679 // The derived quantities are obsolete
680 del_deriv() ;
681
682}
683
685
686 cout <<
687 "Etoile::hydro_euler : hydro_euler must be called via a derived class"
688 << endl << " of Etoile !" << endl ;
689
690 abort() ;
691
692}
693}
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition cmp.h:446
void allocate_all()
Sets the logical state to ETATQCQ (ordinary state) and performs the memory allocation of all the elem...
Definition cmp.C:323
void std_base_scal()
Sets the spectral bases of the Valeur va to the standard ones for a scalar.
Definition cmp.C:644
Tbl & set(int l)
Read/write of the value in a given domain.
Definition cmp.h:724
Equation of state of incompressible matter (Newtonian case).
Definition eos.h:1376
Equation of state of incompressible matter (relativistic case).
Definition eos.h:1206
Polytropic equation of state (Newtonian case).
Definition eos.h:1044
Polytropic equation of state (relativistic case).
Definition eos.h:757
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
Base class for stars *** DEPRECATED : use class Star instead ***.
Definition etoile.h:424
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
double * p_mass_b
Baryon mass.
Definition etoile.h:547
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].
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
Tenseur logn_auto
Total of the logarithm of the part of the lapse N generated principaly by the star.
Definition etoile.h:484
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
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
double * p_ray_eq_pi
Coordinate radius at , .
Definition etoile.h:527
virtual double mass_b() const
Baryon mass.
Tenseur ener
Total energy density in the fluid frame.
Definition etoile.h:460
Etoile(Map &mp_i, int nzet_i, bool relat, const Eos &eos_i)
Standard constructor.
Definition etoile.C:143
double ray_eq_pis2() const
Coordinate radius at , [r_unit].
virtual ~Etoile()
Destructor.
Definition etoile.C:367
double * p_ray_pole
Coordinate radius at .
Definition etoile.h:533
Tenseur d_logn_auto_div
Gradient of logn_auto_div (if k_div!=0 )
Definition etoile.h:501
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition etoile.C:511
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
double * p_ray_eq_pis2
Coordinate radius at , .
Definition etoile.h:524
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
Tenseur a_car
Total conformal factor .
Definition etoile.h:515
void set_enthalpy(const Cmp &)
Assignment of the enthalpy field.
Definition etoile.C:465
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
double ray_pole() const
Coordinate radius at [r_unit].
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 Base_vect_spher & get_bvect_spher() const
Returns the orthonormal vectorial basis associated with the coordinates of the mapping.
Definition map.h:783
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition map.h:765
Multi-domain grid.
Definition grilles.h:273
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_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
Definition grilles.h:462
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 & set(int l)
Read/write of the Tbl in a given domain.
Definition mtbl.h:221
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
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
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition tbl.C:361
Tensor handling *** DEPRECATED : use class Tensor instead ***.
Definition tenseur.h:301
Cmp & set()
Read/write for a scalar (see also operator=(const Cmp&) ).
Definition tenseur.C:824
void sauve(FILE *) const
Save in a file.
Definition tenseur.C:1325
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition tenseur.C:636
void set_std_base()
Set the standard spectal basis of decomposition for each component.
Definition tenseur.C:1170
void set_etat_nondef()
Sets the logical state to ETATNONDEF (undefined state).
Definition tenseur.C:650
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
Lorene prototypes.
Definition app_hor.h:64
Standard units of space, time and mass.