LORENE
star_global.C
1/*
2 * Methods of class Star to compute global quantities
3 */
4
5/*
6 * Copyright (c) 2004 francois Limousin
7 *
8 * This file is part of LORENE.
9 *
10 * LORENE is free software; you can redistribute it and/or modify
11 * it under the terms of the GNU General Public License as published by
12 * the Free Software Foundation; either version 2 of the License, or
13 * (at your option) any later version.
14 *
15 * LORENE is distributed in the hope that it will be useful,
16 * but WITHOUT ANY WARRANTY; without even the implied warranty of
17 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18 * GNU General Public License for more details.
19 *
20 * You should have received a copy of the GNU General Public License
21 * along with LORENE; if not, write to the Free Software
22 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
23 *
24 */
25
26
27char star_global_C[] = "$Header: /cvsroot/Lorene/C++/Source/Star/star_global.C,v 1.6 2014/10/13 08:53:39 j_novak Exp $" ;
28
29/*
30 * $Id: star_global.C,v 1.6 2014/10/13 08:53:39 j_novak Exp $
31 * $Log: star_global.C,v $
32 * Revision 1.6 2014/10/13 08:53:39 j_novak
33 * Lorene classes and functions now belong to the namespace Lorene.
34 *
35 * Revision 1.5 2013/04/25 15:46:06 j_novak
36 * Added special treatment in the case np = 1, for type_p = NONSYM.
37 *
38 * Revision 1.4 2009/10/26 10:54:33 j_novak
39 * Added the case of a NONSYM base in theta.
40 *
41 * Revision 1.3 2007/06/21 19:55:09 k_taniguchi
42 * Introduction of a method to compute ray_eq_3pis2().
43 *
44 * Revision 1.2 2004/01/20 15:20:48 f_limousin
45 * First version
46 *
47 *
48 * $Header: /cvsroot/Lorene/C++/Source/Star/star_global.C,v 1.6 2014/10/13 08:53:39 j_novak Exp $
49 *
50 */
51
52// Headers C
53#include "math.h"
54
55// Headers Lorene
56#include "star.h"
57
58 //--------------------------//
59 // Stellar surface //
60 //--------------------------//
61
62namespace Lorene {
63const Itbl& Star::l_surf() const {
64
65 if (p_l_surf == 0x0) { // a new computation is required
66
67 assert(p_xi_surf == 0x0) ; // consistency check
68
69 int np = mp.get_mg()->get_np(0) ;
70 int nt = mp.get_mg()->get_nt(0) ;
71
72 p_l_surf = new Itbl(np, nt) ;
73 p_xi_surf = new Tbl(np, nt) ;
74
75 double ent0 = 0 ; // definition of the surface
76 double precis = 1.e-15 ;
77 int nitermax = 100 ;
78 int niter ;
79
80 (ent.get_spectral_va()).equipot(ent0, nzet, precis, nitermax, niter, *p_l_surf,
81 *p_xi_surf) ;
82
83 }
84
85 return *p_l_surf ;
86
87}
88
89const Tbl& Star::xi_surf() const {
90
91 if (p_xi_surf == 0x0) { // a new computation is required
92
93 assert(p_l_surf == 0x0) ; // consistency check
94
95 l_surf() ; // the computation is done by l_surf()
96
97 }
98
99 return *p_xi_surf ;
100
101}
102
103
104 //--------------------------//
105 // Coordinate radii //
106 //--------------------------//
107
108double Star::ray_eq() const {
109
110 if (p_ray_eq == 0x0) { // a new computation is required
111
112 const Mg3d& mg = *(mp.get_mg()) ;
113
114 int type_t = mg.get_type_t() ;
115#ifndef NDEBUG
116 int type_p = mg.get_type_p() ;
117#endif
118 int nt = mg.get_nt(0) ;
119
120 assert( (type_t == SYM) || (type_t == NONSYM) ) ;
121 assert( (type_p == SYM) || (type_p == NONSYM) ) ;
122 int k = 0 ;
123 int j = (type_t == SYM ? nt-1 : nt / 2);
124 int l = l_surf()(k, j) ;
125 double xi = xi_surf()(k, j) ;
126 double theta = M_PI / 2 ;
127 double phi = 0 ;
128
129 p_ray_eq = new double( mp.val_r(l, xi, theta, phi) ) ;
130
131 }
132
133 return *p_ray_eq ;
134
135}
136
137
138double Star::ray_eq_pis2() const {
139
140 if (p_ray_eq_pis2 == 0x0) { // a new computation is required
141
142 const Mg3d& mg = *(mp.get_mg()) ;
143
144 int type_t = mg.get_type_t() ;
145 int type_p = mg.get_type_p() ;
146 int nt = mg.get_nt(0) ;
147 int np = mg.get_np(0) ;
148
149 int j = (type_t == SYM ? nt-1 : nt / 2);
150 double theta = M_PI / 2 ;
151 double phi = M_PI / 2 ;
152
153 switch (type_p) {
154
155 case SYM : {
156 int k = np / 2 ;
157 int l = l_surf()(k, j) ;
158 double xi = xi_surf()(k, j) ;
159 p_ray_eq_pis2 = new double( mp.val_r(l, xi, theta, phi) ) ;
160 break ;
161 }
162
163 case NONSYM : {
164 assert( (np == 1) || (np % 4 == 0) ) ;
165 int k = np / 4 ;
166 int l = l_surf()(k, j) ;
167 double xi = xi_surf()(k, j) ;
168 p_ray_eq_pis2 = new double( mp.val_r(l, xi, theta, phi) ) ;
169 break ;
170 }
171
172 default : {
173 cout << "Star::ray_eq_pis2 : the case type_p = "
174 << type_p << " is not contemplated yet !" << endl ;
175 abort() ;
176 }
177 }
178
179 }
180
181 return *p_ray_eq_pis2 ;
182
183}
184
185
186double Star::ray_eq_pi() const {
187
188 if (p_ray_eq_pi == 0x0) { // a new computation is required
189
190 const Mg3d& mg = *(mp.get_mg()) ;
191
192 int type_t = mg.get_type_t() ;
193 int type_p = mg.get_type_p() ;
194 int nt = mg.get_nt(0) ;
195 int np = mg.get_np(0) ;
196
197 assert ( ( type_t == SYM ) || ( type_t == NONSYM ) ) ;
198
199 switch (type_p) {
200
201 case SYM : {
202 p_ray_eq_pi = new double( ray_eq() ) ;
203 break ;
204 }
205
206 case NONSYM : {
207 int k = np / 2 ;
208 int j = (type_t == SYM ? nt-1 : nt/2 ) ;
209 int l = l_surf()(k, j) ;
210 double xi = xi_surf()(k, j) ;
211 double theta = M_PI / 2 ;
212 double phi = M_PI ;
213
214 p_ray_eq_pi = new double( mp.val_r(l, xi, theta, phi) ) ;
215 break ;
216 }
217
218 default : {
219
220 cout << "Star::ray_eq_pi : the case type_p = " << type_p << endl ;
221 cout << " is not contemplated yet !" << endl ;
222 abort() ;
223 break ;
224 }
225 }
226
227 }
228
229 return *p_ray_eq_pi ;
230
231}
232
233double Star::ray_eq_3pis2() const {
234
235 if (p_ray_eq_3pis2 == 0x0) { // a new computation is required
236
237 const Mg3d& mg = *(mp.get_mg()) ;
238
239 int type_t = mg.get_type_t() ;
240 int type_p = mg.get_type_p() ;
241 int nt = mg.get_nt(0) ;
242 int np = mg.get_np(0) ;
243
244 assert( ( type_t == SYM ) || ( type_t == NONSYM ) ) ;
245
246 int j = (type_t == SYM ? nt-1 : nt/2 );
247 double theta = M_PI / 2 ;
248 double phi = 3. * M_PI / 2 ;
249
250 switch (type_p) {
251
252 case SYM : {
254 break ;
255 }
256
257 case NONSYM : {
258 assert( np % 4 == 0 ) ;
259 int k = 3 * np / 4 ;
260 int l = l_surf()(k, j) ;
261 double xi = xi_surf()(k, j) ;
262 p_ray_eq_3pis2 = new double( mp.val_r(l, xi, theta, phi) ) ;
263 break ;
264 }
265
266 default : {
267 cout << "Star::ray_eq_3pis2 : the case type_p = "
268 << type_p << " is not implemented yet !" << endl ;
269 abort() ;
270 }
271 }
272 }
273
274 return *p_ray_eq_3pis2 ;
275
276}
277
278double Star::ray_pole() const {
279
280 if (p_ray_pole == 0x0) { // a new computation is required
281
282#ifndef NDEBUG
283 const Mg3d& mg = *(mp.get_mg()) ;
284 int type_t = mg.get_type_t() ;
285#endif
286 assert( (type_t == SYM) || (type_t == NONSYM) ) ;
287
288 int k = 0 ;
289 int j = 0 ;
290 int l = l_surf()(k, j) ;
291 double xi = xi_surf()(k, j) ;
292 double theta = 0 ;
293 double phi = 0 ;
294
295 p_ray_pole = new double( mp.val_r(l, xi, theta, phi) ) ;
296
297 }
298
299 return *p_ray_pole ;
300
301}
302
303 //--------------------------//
304 // Baryon mass //
305 //--------------------------//
306
307double Star::mass_b() const {
308
309 if (p_mass_b == 0x0) { // a new computation is required
310
311 cout <<
312 "Star::mass_b : in the relativistic case, the baryon mass"
313 << endl <<
314 "computation cannot be performed by the base class Star !"
315 << endl ;
316 abort() ;
317 }
318
319 return *p_mass_b ;
320
321}
322
323 //----------------------------//
324 // Gravitational mass //
325 //----------------------------//
326
327double Star::mass_g() const {
328
329 if (p_mass_g == 0x0) { // a new computation is required
330
331 cout <<
332 "Star::mass_g : in the relativistic case, the gravitational mass"
333 << endl <<
334 "computation cannot be performed by the base class Star !"
335 << endl ;
336 abort() ;
337 }
338
339 return *p_mass_g ;
340
341}
342
343}
Time evolution with partial storage (*** under development ***).
Definition evolution.h:371
Basic integer array class.
Definition itbl.h:122
virtual double val_r(int l, double xi, double theta, double pphi) const =0
Returns the value of the radial coordinate r for a given in a given domain.
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
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_type_p() const
Returns the type of sampling in the direction: SYM : : symmetry with respect to the transformatio...
Definition grilles.h:495
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
Definition grilles.h:457
const Valeur & get_spectral_va() const
Returns va (read only version)
Definition scalar.h:601
virtual double mass_g() const =0
Gravitational mass.
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
double * p_mass_b
Baryon mass.
Definition star.h:268
double ray_eq_3pis2() const
Coordinate radius at , [r_unit].
double * p_ray_eq
Coordinate radius at , .
Definition star.h:242
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
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
double * p_ray_eq_pis2
Coordinate radius at , .
Definition star.h:245
double * p_ray_eq_3pis2
Coordinate radius at , .
Definition star.h:251
double ray_eq() const
Coordinate radius at , [r_unit].
double * p_mass_g
Gravitational mass.
Definition star.h:269
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
double ray_eq_pis2() const
Coordinate radius at , [r_unit].
double ray_eq_pi() const
Coordinate radius at , [r_unit].
Scalar ent
Log-enthalpy.
Definition star.h:190
Map & mp
Mapping associated with the star.
Definition star.h:180
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
double ray_pole() const
Coordinate radius at [r_unit].
virtual double mass_b() const =0
Baryon mass.
Basic array class.
Definition tbl.h:161
Lorene prototypes.
Definition app_hor.h:64