LORENE
et_rot_global.C
1/*
2 * Methods for computing global quantities within the class Etoile_rot
3 *
4 * (see file etoile.h for documentation)
5 */
6
7/*
8 * Copyright (c) 2000-2001 Eric Gourgoulhon
9 *
10 * This file is part of LORENE.
11 *
12 * LORENE is free software; you can redistribute it and/or modify
13 * it under the terms of the GNU General Public License as published by
14 * the Free Software Foundation; either version 2 of the License, or
15 * (at your option) any later version.
16 *
17 * LORENE is distributed in the hope that it will be useful,
18 * but WITHOUT ANY WARRANTY; without even the implied warranty of
19 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
20 * GNU General Public License for more details.
21 *
22 * You should have received a copy of the GNU General Public License
23 * along with LORENE; if not, write to the Free Software
24 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
25 *
26 */
27
28
29char et_rot_global_C[] = "$Header: /cvsroot/Lorene/C++/Source/Etoile/et_rot_global.C,v 1.9 2015/06/12 12:38:25 j_novak Exp $" ;
30
31/*
32 * $Id: et_rot_global.C,v 1.9 2015/06/12 12:38:25 j_novak Exp $
33 * $Log: et_rot_global.C,v $
34 * Revision 1.9 2015/06/12 12:38:25 j_novak
35 * Implementation of the corrected formula for the quadrupole momentum.
36 *
37 * Revision 1.8 2015/06/10 14:37:44 a_sourie
38 * Corrected the formula for the quadrupole.
39 *
40 * Revision 1.7 2014/10/13 08:52:57 j_novak
41 * Lorene classes and functions now belong to the namespace Lorene.
42 *
43 * Revision 1.6 2014/10/06 15:13:09 j_novak
44 * Modified #include directives to use c++ syntax.
45 *
46 * Revision 1.5 2012/08/12 17:48:35 p_cerda
47 * Magnetstar: New classes for magnetstar. Allowing for non-equatorial symmetry in Etoile et al. Adding B_phi in Et_rot_mag.
48 *
49 * Revision 1.4 2004/03/25 10:29:06 j_novak
50 * All LORENE's units are now defined in the namespace Unites (in file unites.h).
51 *
52 * Revision 1.3 2003/11/03 16:47:13 e_gourgoulhon
53 * Corrected error in grv2() in the Newtonian case.
54 *
55 * Revision 1.2 2002/04/05 09:09:37 j_novak
56 * The inversion of the EOS for 2-fluids polytrope has been modified.
57 * Some errors in the determination of the surface were corrected.
58 *
59 * Revision 1.1.1.1 2001/11/20 15:19:28 e_gourgoulhon
60 * LORENE
61 *
62 * Revision 1.5 2000/11/19 18:52:09 eric
63 * grv2() operationnelle.
64 *
65 * Revision 1.4 2000/10/12 15:34:55 eric
66 * Calcul de la masse grav, de GRV3 et du moment quadrupolaire.
67 *
68 * Revision 1.3 2000/08/31 11:25:58 eric
69 * *** empty log message ***
70 *
71 * Revision 1.2 2000/08/25 12:28:16 eric
72 * *** empty log message ***
73 *
74 * Revision 1.1 2000/07/20 15:32:56 eric
75 * Initial revision
76 *
77 *
78 * $Header: /cvsroot/Lorene/C++/Source/Etoile/et_rot_global.C,v 1.9 2015/06/12 12:38:25 j_novak Exp $
79 *
80 */
81
82// Headers C
83#include <cstdlib>
84#include <cmath>
85
86// Headers Lorene
87#include "etoile.h"
88#include "unites.h"
89
90 //--------------------------//
91 // Stellar surface //
92 //--------------------------//
93
94namespace Lorene {
95const Itbl& Etoile_rot::l_surf() const {
96
97 if (p_l_surf == 0x0) { // a new computation is required
98
99 assert(p_xi_surf == 0x0) ; // consistency check
100
101 int np = mp.get_mg()->get_np(0) ;
102 int nt = mp.get_mg()->get_nt(0) ;
103
104 p_l_surf = new Itbl(np, nt) ;
105 p_xi_surf = new Tbl(np, nt) ;
106
107 double ent0 = 0 ; // definition of the surface
108 double precis = 1.e-15 ;
109 int nitermax = 100 ;
110 int niter ;
111
112 (ent().va).equipot(ent0, nzet, precis, nitermax, niter, *p_l_surf,
113 *p_xi_surf) ;
114
115 }
116
117 return *p_l_surf ;
118
119}
120
121 //--------------------------//
122 // Baryon mass //
123 //--------------------------//
124
125double Etoile_rot::mass_b() const {
126
127 if (p_mass_b == 0x0) { // a new computation is required
128
129 if (relativistic) {
130
131 Cmp dens = a_car() * bbb() * gam_euler() * nbar() ;
132
133 dens.std_base_scal() ;
134
135 p_mass_b = new double( dens.integrale() ) ;
136
137
138 }
139 else{ // Newtonian case
140 assert(nbar.get_etat() == ETATQCQ) ;
141
142 p_mass_b = new double( nbar().integrale() ) ;
143
144 }
145
146 }
147
148 return *p_mass_b ;
149
150}
151
152
153 //----------------------------//
154 // Gravitational mass //
155 //----------------------------//
156
157double Etoile_rot::mass_g() const {
158
159 if (p_mass_g == 0x0) { // a new computation is required
160
161 if (relativistic) {
162
163 Tenseur source = nnn * (ener_euler + s_euler)
164 + 2 * bbb * (ener_euler + press)
165 * tnphi * uuu ;
166 source = a_car * bbb * source ;
167
168 source.set_std_base() ;
169
170 p_mass_g = new double( source().integrale() ) ;
171
172
173 }
174 else{ // Newtonian case
175 p_mass_g = new double( mass_b() ) ; // in the Newtonian case
176 // M_g = M_b
177 }
178 }
179
180 return *p_mass_g ;
181
182}
183
184 //----------------------------//
185 // Angular momentum //
186 //----------------------------//
187
188double Etoile_rot::angu_mom() const {
189
190 if (p_angu_mom == 0x0) { // a new computation is required
191
192 Cmp dens = uuu() ;
193
194 dens.mult_r() ; // Multiplication by
195 dens.va = (dens.va).mult_st() ; // r sin(theta)
196
197 if (relativistic) {
198 dens = a_car() * b_car() * (ener_euler() + press())
199 * dens ;
200 }
201 else { // Newtonian case
202 dens = nbar() * dens ;
203 }
204
205 dens.std_base_scal() ;
206
207 p_angu_mom = new double( dens.integrale() ) ;
208
209 }
210
211 return *p_angu_mom ;
212
213}
214
215
216 //----------------------------//
217 // T/W //
218 //----------------------------//
219
220double Etoile_rot::tsw() const {
221
222 if (p_tsw == 0x0) { // a new computation is required
223
224 double tcin = 0.5 * omega * angu_mom() ;
225
226 if (relativistic) {
227
228 Cmp dens = a_car() * bbb() * gam_euler() * ener() ;
229 dens.std_base_scal() ;
230 double mass_p = dens.integrale() ;
231
232 p_tsw = new double( tcin / ( mass_p + tcin - mass_g() ) ) ;
233
234 }
235 else { // Newtonian case
236 Cmp dens = 0.5 * nbar() * logn() ;
237 dens.std_base_scal() ;
238 double wgrav = dens.integrale() ;
239 p_tsw = new double( tcin / fabs(wgrav) ) ;
240 }
241
242
243 }
244
245 return *p_tsw ;
246
247}
248
249
250 //----------------------------//
251 // GRV2 //
252 //----------------------------//
253
254double Etoile_rot::grv2() const {
255
256 using namespace Unites ;
257 if (p_grv2 == 0x0) { // a new computation is required
258
259 Tenseur sou_m(mp) ;
260 if (relativistic) {
261 sou_m = 2 * qpig * a_car * (press + (ener_euler+press)
262 * uuu*uuu ) ;
263 }
264 else {
265 sou_m = 2 * qpig * (press + nbar * uuu*uuu ) ;
266 }
267
268 Tenseur sou_q = 1.5 * ak_car
270 logn.gradient_spher() ) ;
271
272 p_grv2 = new double( double(1) - lambda_grv2(sou_m(), sou_q()) ) ;
273
274 }
275
276 return *p_grv2 ;
277
278}
279
280
281 //----------------------------//
282 // GRV3 //
283 //----------------------------//
284
285double Etoile_rot::grv3(ostream* ost) const {
286
287 using namespace Unites ;
288
289 if (p_grv3 == 0x0) { // a new computation is required
290
291 Tenseur source(mp) ;
292
293 // Gravitational term [cf. Eq. (43) of Gourgoulhon & Bonazzola
294 // ------------------ Class. Quantum Grav. 11, 443 (1994)]
295
296 if (relativistic) {
297 Tenseur alpha = dzeta - logn ;
298 Tenseur beta = log( bbb ) ;
299 beta.set_std_base() ;
300
301 source = 0.75 * ak_car
304 + 0.5 * flat_scalar_prod(alpha.gradient_spher(),
305 beta.gradient_spher() ) ;
306
307 Cmp aa = alpha() - 0.5 * beta() ;
308 Cmp daadt = aa.srdsdt() ; // 1/r d/dth
309
310 // What follows is valid only for a mapping of class Map_radial :
311 const Map_radial* mpr = dynamic_cast<const Map_radial*>(&mp) ;
312 if (mpr == 0x0) {
313 cout << "Etoile_rot::grv3: the mapping does not belong"
314 << " to the class Map_radial !" << endl ;
315 abort() ;
316 }
317
318 // Computation of 1/tan(theta) * 1/r daa/dtheta
319 if (daadt.get_etat() == ETATQCQ) {
320 Valeur& vdaadt = daadt.va ;
321 vdaadt = vdaadt.ssint() ; // division by sin(theta)
322 vdaadt = vdaadt.mult_ct() ; // multiplication by cos(theta)
323 }
324
325 Cmp temp = aa.dsdr() + daadt ;
326 temp = ( bbb() - a_car()/bbb() ) * temp ;
327 temp.std_base_scal() ;
328
329 // Division by r
330 Valeur& vtemp = temp.va ;
331 vtemp = vtemp.sx() ; // division by xi in the nucleus
332 // Id in the shells
333 // division by xi-1 in the ZEC
334 vtemp = (mpr->xsr) * vtemp ; // multiplication by xi/r in the nucleus
335 // by 1/r in the shells
336 // by r(xi-1) in the ZEC
337
338 // In the ZEC, a multiplication by r has been performed instead
339 // of the division:
340 temp.set_dzpuis( temp.get_dzpuis() + 2 ) ;
341
342 source = bbb() * source() + 0.5 * temp ;
343
344 }
345 else{
346 source = - 0.5 * flat_scalar_prod(logn.gradient_spher(),
347 logn.gradient_spher() ) ;
348 }
349
350 source.set_std_base() ;
351
352 double int_grav = source().integrale() ;
353
354 // Matter term
355 // -----------
356
357 if (relativistic) {
358 source = qpig * a_car * bbb * s_euler ;
359 }
360 else{
361 source = qpig * ( 3 * press + nbar * uuu * uuu ) ;
362 }
363
364 source.set_std_base() ;
365
366 double int_mat = source().integrale() ;
367
368 // Virial error
369 // ------------
370 if (ost != 0x0) {
371 *ost << "Etoile_rot::grv3 : gravitational term : " << int_grav
372 << endl ;
373 *ost << "Etoile_rot::grv3 : matter term : " << int_mat
374 << endl ;
375 }
376
377 p_grv3 = new double( (int_grav + int_mat) / int_mat ) ;
378
379 }
380
381 return *p_grv3 ;
382
383}
384
385
386 //----------------------------//
387 // R_circ //
388 //----------------------------//
389
390double Etoile_rot::r_circ() const {
391
392 if (p_r_circ == 0x0) { // a new computation is required
393
394 // Index of the point at phi=0, theta=pi/2 at the surface of the star:
395 const Mg3d* mg = mp.get_mg() ;
396
397
398 int l_b = nzet - 1 ;
399 int i_b = mg->get_nr(l_b) - 1 ;
400 int j_b;
401 if (mg->get_type_t() == SYM) {
402 j_b = mg->get_nt(l_b) - 1 ;
403 }else{
404 j_b = (mg->get_nt(l_b) - 1)/2 ;
405 }
406 int k_b = 0 ;
407
408 p_r_circ = new double( bbb()(l_b, k_b, j_b, i_b) * ray_eq() ) ;
409
410 }
411
412 return *p_r_circ ;
413
414}
415
416 //----------------------------//
417 // Surface area //
418 //----------------------------//
419
420 double Etoile_rot::area() const {
421
422 if (p_area == 0x0) { // a new computation is required
423 const Mg3d& mg = *(mp.get_mg()) ;
424 int np = mg.get_np(0) ;
425 int nt = mg.get_nt(0) ;
426 assert(np == 1) ; //Only valid for axisymmetric configurations
427
428 const Map_radial* mp_rad = dynamic_cast<const Map_radial*>(&mp) ;
429 assert(mp_rad != 0x0) ;
430
431 Valeur va_r(mg.get_angu()) ;
432 va_r.annule_hard() ;
433 Itbl lsurf = l_surf() ;
434 Tbl xisurf = xi_surf() ;
435 for (int k=0; k<np; k++) {
436 for (int j=0; j<nt; j++) {
437 int l_star = lsurf(k, j) ;
438 double xi_star = xisurf(k, j) ;
439 va_r.set(0, k, j, 0) = mp_rad->val_r_jk(l_star, xi_star, j, k) ;
440 }
441 }
442 va_r.std_base_scal() ;
443
444 Valeur integ(mg.get_angu()) ;
445 Valeur dr = va_r.dsdt() ;
446 integ = sqrt(va_r*va_r + dr*dr) ;
447 Cmp aaaa = get_a_car()() ;
448 Valeur a2 = aaaa.va ; a2.std_base_scal() ;
449 Cmp bbbb = get_bbb()() ;
450 Valeur b = bbbb.va ; b.std_base_scal() ;
451 for (int k=0; k<np; k++) {
452 for (int j=0; j<nt; j++) {
453 integ.set(0, k, j, 0) *= sqrt(a2.val_point_jk(lsurf(k, j), xisurf(k, j), j, k))
454 * b.val_point_jk(lsurf(k, j), xisurf(k, j), j, k) * va_r(0, k, j, 0) ;
455 }
456 }
457 integ.std_base_scal() ;
458 Valeur integ2 = integ.mult_st() ;
459 double surftot = 0. ;
460 for (int j=0; j<nt; j++) {
461 surftot += (*integ2.c_cf)(0, 0, j, 0) / double(2*j+1) ;
462 }
463
464 p_area = new double( 4*M_PI*surftot) ;
465
466 }
467
468 return *p_area ;
469
470}
471
472 double Etoile_rot::mean_radius() const {
473
474 return sqrt(area()/(4*M_PI)) ;
475
476 }
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492 //----------------------------//
493 // Flattening //
494 //----------------------------//
495
496double Etoile_rot::aplat() const {
497
498 if (p_aplat == 0x0) { // a new computation is required
499
500 p_aplat = new double( ray_pole() / ray_eq() ) ;
501
502 }
503
504 return *p_aplat ;
505
506}
507
508
509 //----------------------------//
510 // Z_eq_f //
511 //----------------------------//
512
513double Etoile_rot::z_eqf() const {
514
515 if (p_z_eqf == 0x0) { // a new computation is required
516
517 // Index of the point at phi=0, theta=pi/2 at the surface of the star:
518 const Mg3d* mg = mp.get_mg() ;
519 int l_b = nzet - 1 ;
520 int i_b = mg->get_nr(l_b) - 1 ;
521 int j_b;
522 if (mg->get_type_t() == SYM) {
523 j_b = mg->get_nt(l_b) - 1 ;
524 }else{
525 j_b = (mg->get_nt(l_b) - 1)/2 ;
526 }
527 int k_b = 0 ;
528
529 double u_eq = uuu()(l_b, k_b, j_b, i_b) ;
530 double n_eq = nnn()(l_b, k_b, j_b, i_b) ;
531 double nphi_eq = nphi()(l_b, k_b, j_b, i_b) ;
532
533 p_z_eqf = new double( sqrt((1.-u_eq)/(1.+u_eq))
534 / (n_eq + nphi_eq * r_circ() )
535 - 1. ) ;
536 }
537
538 return *p_z_eqf ;
539
540}
541 //----------------------------//
542 // Z_eq_b //
543 //----------------------------//
544
545double Etoile_rot::z_eqb() const {
546
547 if (p_z_eqb == 0x0) { // a new computation is required
548
549 // Index of the point at phi=0, theta=pi/2 at the surface of the star:
550 const Mg3d* mg = mp.get_mg() ;
551 int l_b = nzet - 1 ;
552 int i_b = mg->get_nr(l_b) - 1 ;
553 int j_b;
554 if (mg->get_type_t() == SYM) {
555 j_b = mg->get_nt(l_b) - 1 ;
556 }else{
557 j_b = (mg->get_nt(l_b) - 1) / 2 ;
558 }
559 int k_b = 0 ;
560
561 double u_eq = uuu()(l_b, k_b, j_b, i_b) ;
562 double n_eq = nnn()(l_b, k_b, j_b, i_b) ;
563 double nphi_eq = nphi()(l_b, k_b, j_b, i_b) ;
564
565 p_z_eqb = new double( sqrt((1.+u_eq)/(1.-u_eq))
566 / (n_eq - nphi_eq * r_circ() )
567 - 1. ) ;
568
569 }
570
571 return *p_z_eqb ;
572
573}
574
575
576 //----------------------------//
577 // Z_pole //
578 //----------------------------//
579
580double Etoile_rot::z_pole() const {
581
582 if (p_z_pole == 0x0) { // a new computation is required
583
584 double n_pole = nnn().val_point(ray_pole(), 0., 0.) ;
585
586 p_z_pole = new double( 1. / n_pole - 1. ) ;
587
588 }
589
590 return *p_z_pole ;
591
592}
593
594
595 //----------------------------//
596 // Quadrupole moment //
597 //----------------------------//
598
599
600double Etoile_rot::mom_quad() const {
601
602 using namespace Unites ;
603
604 if (p_mom_quad == 0x0) { // a new computation is required
605
606 p_mom_quad = new double( mom_quad_old() ) ;
607 if (relativistic) {
608 double b = mom_quad_Bo() / ( mass_g() * mass_g() ) ;
609 *p_mom_quad -= 4./3. * ( 1./4. + b ) * pow(mass_g(), 3) * ggrav * ggrav ;
610 }
611 }
612
613 return *p_mom_quad ;
614
615}
616
617
619
620 using namespace Unites ;
621
622 if (p_mom_quad_Bo == 0x0) { // a new computation is required
623
624 Cmp dens(mp) ;
625
626 dens = press() ;
627 dens = a_car() * bbb() * nnn() * dens ;
628 dens.mult_rsint() ;
629 dens.std_base_scal() ;
630
631 p_mom_quad_Bo = new double( - 32. * dens.integrale() / qpig ) ;
632
633 }
634
635 return *p_mom_quad_Bo ;
636
637}
638
639
640
642
643 using namespace Unites ;
644
645 if (p_mom_quad_old == 0x0) { // a new computation is required
646
647 // Source for of the Poisson equation for nu
648 // -----------------------------------------
649
650 Tenseur source(mp) ;
651
652 if (relativistic) {
653 Tenseur beta = log(bbb) ;
654 beta.set_std_base() ;
655 source = qpig * a_car *( ener_euler + s_euler )
657 logn.gradient_spher() + beta.gradient_spher()) ;
658 }
659 else {
660 source = qpig * nbar ;
661 }
662 source.set_std_base() ;
663
664 // Multiplication by -r^2 P_2(cos(theta))
665 // [cf Eq.(7) of Salgado et al. Astron. Astrophys. 291, 155 (1994) ]
666 // ------------------------------------------------------------------
667
668 // Multiplication by r^2 :
669 // ----------------------
670 Cmp& csource = source.set() ;
671 csource.mult_r() ;
672 csource.mult_r() ;
673 if (csource.check_dzpuis(2)) {
674 csource.inc2_dzpuis() ;
675 }
676
677 // Muliplication by cos^2(theta) :
678 // -----------------------------
679 Cmp temp = csource ;
680
681 // What follows is valid only for a mapping of class Map_radial :
682 assert( dynamic_cast<const Map_radial*>(&mp) != 0x0 ) ;
683
684 if (temp.get_etat() == ETATQCQ) {
685 Valeur& vtemp = temp.va ;
686 vtemp = vtemp.mult_ct() ; // multiplication by cos(theta)
687 vtemp = vtemp.mult_ct() ; // multiplication by cos(theta)
688 }
689
690 // Muliplication by -P_2(cos(theta)) :
691 // ----------------------------------
692 source = 0.5 * source() - 1.5 * temp ;
693
694 // Final result
695 // ------------
696
697 p_mom_quad_old = new double(- source().integrale() / qpig ) ;
698
699 }
700
701 return *p_mom_quad_old ;
702
703}
704
705
706
707
708}
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition cmp.h:446
void mult_rsint()
Multiplication by .
int get_etat() const
Returns the logical state.
Definition cmp.h:899
Valeur va
The numerical value of the Cmp
Definition cmp.h:464
void std_base_scal()
Sets the spectral bases of the Valeur va to the standard ones for a scalar.
Definition cmp.C:644
int get_dzpuis() const
Returns dzpuis.
Definition cmp.h:903
void mult_r()
Multiplication by r everywhere.
Definition cmp_r_manip.C:91
void inc2_dzpuis()
Increases by 2 the value of dzpuis and changes accordingly the values of the Cmp in the external comp...
void set_dzpuis(int)
Set a value to dzpuis.
Definition cmp.C:654
double integrale() const
Computes the integral over all space of *this .
Definition cmp_integ.C:55
const Cmp & srdsdt() const
Returns of *this .
Definition cmp_deriv.C:105
bool check_dzpuis(int dzi) const
Returns false if the last domain is compactified and *this is not zero in this domain and dzpuis is n...
Definition cmp.C:715
const Cmp & dsdr() const
Returns of *this .
Definition cmp_deriv.C:84
Tenseur uuu
Norm of u_euler.
Definition etoile.h:1518
double omega
Rotation angular velocity ([f_unit] )
Definition etoile.h:1501
double * p_z_pole
Redshift factor at North pole.
Definition etoile.h:1640
virtual double mom_quad() const
Quadrupole moment.
double * p_z_eqf
Forward redshift factor at equator.
Definition etoile.h:1638
virtual double r_circ() const
Circumferential radius.
Tenseur & logn
Metric potential = logn_auto.
Definition etoile.h:1521
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...
double * p_mom_quad_old
Part of the quadrupole moment.
Definition etoile.h:1642
double * p_z_eqb
Backward redshift factor at equator.
Definition etoile.h:1639
virtual double mass_g() const
Gravitational mass.
double * p_r_circ
Circumferential radius.
Definition etoile.h:1635
virtual double mom_quad_Bo() const
Part of the quadrupole moment.
virtual double aplat() const
Flatening r_pole/r_eq.
virtual double grv3(ostream *ost=0x0) const
Error on the virial identity GRV3.
const Tenseur & get_bbb() const
Returns the metric factor B.
Definition etoile.h:1712
virtual double tsw() const
Ratio T/W.
Tenseur nphi
Metric coefficient .
Definition etoile.h:1510
virtual double mass_b() const
Baryon mass.
static double lambda_grv2(const Cmp &sou_m, const Cmp &sou_q)
Computes the coefficient which ensures that the GRV2 virial identity is satisfied.
Tenseur bbb
Metric factor B.
Definition etoile.h:1504
double * p_area
Surface area.
Definition etoile.h:1636
Tenseur ak_car
Scalar .
Definition etoile.h:1586
Tenseur & dzeta
Metric potential = beta_auto.
Definition etoile.h:1534
virtual double mean_radius() const
Mean radius.
double * p_grv3
Error on the virial identity GRV3.
Definition etoile.h:1634
double * p_grv2
Error on the virial identity GRV2.
Definition etoile.h:1633
double * p_aplat
Flatening r_pole/r_eq.
Definition etoile.h:1637
double * p_mom_quad_Bo
Part of the quadrupole moment.
Definition etoile.h:1643
virtual double grv2() const
Error on the virial identity GRV2.
double * p_mom_quad
Quadrupole moment.
Definition etoile.h:1641
virtual double area() const
Surface area.
double * p_angu_mom
Angular momentum.
Definition etoile.h:1631
virtual double angu_mom() const
Angular momentum.
virtual double z_eqf() const
Forward redshift factor at equator.
virtual double mom_quad_old() const
Part of the quadrupole moment.
double * p_tsw
Ratio T/W.
Definition etoile.h:1632
Tenseur b_car
Square of the metric factor B.
Definition etoile.h:1507
virtual double z_eqb() const
Backward redshift factor at equator.
Tenseur tnphi
Component of the shift vector.
Definition etoile.h:1515
virtual double z_pole() const
Redshift factor at North pole.
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
double ray_eq() const
Coordinate radius at , [r_unit].
Tenseur nnn
Total lapse function.
Definition etoile.h:509
Tenseur nbar
Baryon density in the fluid frame.
Definition etoile.h:459
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
Tenseur ener
Total energy density in the fluid frame.
Definition etoile.h:460
const Tenseur & get_a_car() const
Returns the total conformal factor .
Definition etoile.h:733
Tenseur press
Fluid pressure.
Definition etoile.h:461
bool relativistic
Indicator of relativity: true for a relativistic star, false for a Newtonian one.
Definition etoile.h:437
Tenseur ener_euler
Total energy density in the Eulerian frame.
Definition etoile.h:465
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
const Tbl & xi_surf() const
Description of the stellar surface: returns a 2-D Tbl containing the values of the radial coordinat...
Tenseur a_car
Total conformal factor .
Definition etoile.h:515
double ray_pole() const
Coordinate radius at [r_unit].
Basic integer array class.
Definition itbl.h:122
Base class for pure radial mappings.
Definition map.h:1536
Coord xsr
in the nucleus; \ 1/R in the non-compactified shells; \ in the compactified outer domain.
Definition map.h:1549
virtual double val_r_jk(int l, double xi, int j, int k) const =0
Returns the value of the radial coordinate r for a given and a given collocation point in in a give...
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 Mg3d * get_angu() const
Returns the pointer on the associated angular grid.
Definition mg3d.C:473
int get_type_t() const
Returns the type of sampling in the direction: SYM : : symmetry with respect to the equatorial pl...
Definition grilles.h:485
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_nr(int l) const
Returns the number of points in the radial direction ( ) in domain no. l.
Definition grilles.h:452
Basic array class.
Definition tbl.h:161
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 set_std_base()
Set the standard spectal basis of decomposition for each component.
Definition tenseur.C:1170
const Tenseur & gradient_spher() const
Returns the gradient of *this (Spherical coordinates) (scalar field only).
Definition tenseur.C:1548
int get_etat() const
Returns the logical state.
Definition tenseur.h:707
Values and coefficients of a (real-value) function.
Definition valeur.h:287
const Valeur & mult_ct() const
Returns applied to *this.
double val_point_jk(int l, double x, int j, int k) const
Computes the value of the field represented by *this at an arbitrary point in , but collocation point...
Definition valeur.C:900
const Valeur & sx() const
Returns (r -sampling = RARE ) \ Id (r sampling = FIN ) \ (r -sampling = UNSURR )
Definition valeur_sx.C:110
const Valeur & dsdt() const
Returns of *this.
Tbl & set(int l)
Read/write of the value in a given domain (configuration space).
Definition valeur.h:363
const Valeur & ssint() const
Returns of *this.
Mtbl_cf * c_cf
Coefficients of the spectral expansion of the function.
Definition valeur.h:302
const Valeur & mult_st() const
Returns applied to *this.
void std_base_scal()
Sets the bases for spectral expansions (member base ) to the standard ones for a scalar.
Definition valeur.C:824
void annule_hard()
Sets the Valeur to zero in a hard way.
Definition valeur.C:723
Cmp sqrt(const Cmp &)
Square root.
Definition cmp_math.C:220
Cmp pow(const Cmp &, int)
Power .
Definition cmp_math.C:348
Cmp log(const Cmp &)
Neperian logarithm.
Definition cmp_math.C:296
Tenseur flat_scalar_prod(const Tenseur &t1, const Tenseur &t2)
Scalar product of two Tenseur when the metric is : performs the contraction of the last index of t1 w...
Lorene prototypes.
Definition app_hor.h:64
Standard units of space, time and mass.