LORENE
star_rot_global.C
1/*
2 * Methods for computing global quantities within the class Star_rot
3 *
4 * (see file star_rot.h for documentation)
5 */
6
7/*
8 * Copyright (c) 2010 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 star_rot_global_C[] = "$Header: /cvsroot/Lorene/C++/Source/Star/star_rot_global.C,v 1.5 2015/05/19 09:30:56 j_novak Exp $" ;
30
31/*
32 * $Id: star_rot_global.C,v 1.5 2015/05/19 09:30:56 j_novak Exp $
33 * $Log: star_rot_global.C,v $
34 * Revision 1.5 2015/05/19 09:30:56 j_novak
35 * New methods for computing the area of the star and its mean radius.
36 *
37 * Revision 1.4 2014/10/13 08:53:39 j_novak
38 * Lorene classes and functions now belong to the namespace Lorene.
39 *
40 * Revision 1.3 2014/10/06 15:13:17 j_novak
41 * Modified #include directives to use c++ syntax.
42 *
43 * Revision 1.2 2010/01/25 22:33:35 e_gourgoulhon
44 * Debugging...
45 *
46 * Revision 1.1 2010/01/25 18:15:52 e_gourgoulhon
47 * First version.
48 *
49 *
50 *
51 * $Header: /cvsroot/Lorene/C++/Source/Star/star_rot_global.C,v 1.5 2015/05/19 09:30:56 j_novak Exp $
52 *
53 */
54
55// Headers C
56#include <cstdlib>
57#include <cmath>
58
59// Headers Lorene
60#include "star_rot.h"
61#include "unites.h"
62
63 //--------------------------//
64 // Stellar surface //
65 //--------------------------//
66
67namespace Lorene {
68 const Itbl& Star_rot::l_surf() const {
69
70 if (p_l_surf == 0x0) { // a new computation is required
71
72 assert(p_xi_surf == 0x0) ; // consistency check
73
74 int np = mp.get_mg()->get_np(0) ;
75 int nt = mp.get_mg()->get_nt(0) ;
76
77 p_l_surf = new Itbl(np, nt) ;
78 p_xi_surf = new Tbl(np, nt) ;
79
80 double ent0 = 0 ; // definition of the surface
81 double precis = 1.e-15 ;
82 int nitermax = 100 ;
83 int niter ;
84
85 (ent.get_spectral_va()).equipot(ent0, nzet, precis, nitermax, niter, *p_l_surf,
86 *p_xi_surf) ;
87
88 }
89
90 return *p_l_surf ;
91
92 }
93
94 //------------------------------//
95 // Baryon mass //
96 //------------------------------//
97
98 double Star_rot::mass_b() const {
99
100 if (p_mass_b == 0x0) { // a new computation is required
101
102 if (relativistic) {
103
105
106 p_mass_b = new double( dens.integrale() ) ;
107 }
108 else{ // Newtonian case
109 assert(nbar.get_etat() == ETATQCQ) ;
110
111 p_mass_b = new double( nbar.integrale() ) ;
112
113 }
114
115 }
116 return *p_mass_b ;
117 }
118
119
120 //----------------------------//
121 // Gravitational mass //
122 //----------------------------//
123
124 double Star_rot::mass_g() const {
125
126 if (p_mass_g == 0x0) { // a new computation is required
127
128 if (relativistic) {
129
131 + 2 * bbb * (ener_euler + press)
132 * tnphi * uuu ;
133 source = a_car * bbb * source ;
134 source.std_spectral_base() ;
135
136 p_mass_g = new double( source.integrale() ) ;
137 }
138 else{ // Newtonian case
139 p_mass_g = new double( mass_b() ) ; // in the Newtonian case
140 // M_g = M_b
141 }
142 }
143 return *p_mass_g ;
144}
145
146 //----------------------------//
147 // Angular momentum //
148 //----------------------------//
149
150 double Star_rot::angu_mom() const {
151
152 if (p_angu_mom == 0x0) { // a new computation is required
153
154 Scalar dens = uuu ;
155
156 dens.mult_r() ; // Multiplication by
157 dens.set_spectral_va() = (dens.get_spectral_va()).mult_st() ; // r sin(theta)
158
159 if (relativistic) {
160 dens = a_car * b_car * (ener_euler + press) * dens ;
161 }
162 else { // Newtonian case
163 dens = nbar * dens ;
164 }
165
166 p_angu_mom = new double( dens.integrale() ) ;
167 }
168 return *p_angu_mom ;
169 }
170
171
172 //----------------------------//
173 // T/W //
174 //----------------------------//
175
176 double Star_rot::tsw() const {
177
178 if (p_tsw == 0x0) { // a new computation is required
179
180 double tcin = 0.5 * omega * angu_mom() ;
181
182 if (relativistic) {
183
185 double mass_p = dens.integrale() ;
186
187 p_tsw = new double( tcin / ( mass_p + tcin - mass_g() ) ) ;
188
189 }
190 else { // Newtonian case
191 Scalar dens = 0.5 * nbar * logn ;
192 double wgrav = dens.integrale() ;
193 p_tsw = new double( tcin / fabs(wgrav) ) ;
194 }
195 }
196 return *p_tsw ;
197}
198
199
200 //----------------------------//
201 // GRV2 //
202 //----------------------------//
203
204double Star_rot::grv2() const {
205
206 using namespace Unites ;
207 if (p_grv2 == 0x0) { // a new computation is required
208
209 Scalar sou_m(mp) ;
210 if (relativistic) {
211 sou_m = 2 * qpig * a_car * (press + (ener_euler+press)
212 * uuu*uuu ) ;
213 }
214 else {
215 sou_m = 2 * qpig * (press + nbar * uuu*uuu ) ;
216 }
217
219 Scalar sou_q = 1.5 * ak_car
220 - dlogn(1)*dlogn(1) - dlogn(2)*dlogn(2) - dlogn(3)*dlogn(3) ;
221
222 p_grv2 = new double( double(1) - lambda_grv2(sou_m, sou_q) ) ;
223
224 }
225
226 return *p_grv2 ;
227
228}
229
230
231 //----------------------------//
232 // GRV3 //
233 //----------------------------//
234
235double Star_rot::grv3(ostream* ost) const {
236
237 using namespace Unites ;
238
239 if (p_grv3 == 0x0) { // a new computation is required
240
241 Scalar source(mp) ;
242
243 // Gravitational term [cf. Eq. (43) of Gourgoulhon & Bonazzola
244 // ------------------ Class. Quantum Grav. 11, 443 (1994)]
245
247
248 if (relativistic) {
249 Scalar alpha = dzeta - logn ;
250 Scalar bet = log( bbb ) ;
251 bet.std_spectral_base() ;
252
254 Vector dbet = bet.derive_cov( mp.flat_met_spher() ) ;
255
256 source = 0.75 * ak_car
257 - dlogn(1)*dlogn(1) - dlogn(2)*dlogn(2) - dlogn(3)*dlogn(3)
258 + 0.5 * ( dalpha(1)*dbet(1) + dalpha(2)*dbet(2) + dalpha(3)*dbet(3) ) ;
259
260 Scalar aa = alpha - 0.5 * bet ;
261 Scalar daadt = aa.srdsdt() ; // 1/r d/dth
262
263 // What follows is valid only for a mapping of class Map_radial :
264 const Map_radial* mpr = dynamic_cast<const Map_radial*>(&mp) ;
265 if (mpr == 0x0) {
266 cout << "Star_rot::grv3: the mapping does not belong"
267 << " to the class Map_radial !" << endl ;
268 abort() ;
269 }
270
271 // Computation of 1/tan(theta) * 1/r daa/dtheta
272 if (daadt.get_etat() == ETATQCQ) {
273 Valeur& vdaadt = daadt.set_spectral_va() ;
274 vdaadt = vdaadt.ssint() ; // division by sin(theta)
275 vdaadt = vdaadt.mult_ct() ; // multiplication by cos(theta)
276 }
277
278 Scalar temp = aa.dsdr() + daadt ;
279 temp = ( bbb - a_car/bbb ) * temp ;
280 temp.std_spectral_base() ;
281
282 // Division by r
283 Valeur& vtemp = temp.set_spectral_va() ;
284 vtemp = vtemp.sx() ; // division by xi in the nucleus
285 // Id in the shells
286 // division by xi-1 in the ZEC
287 vtemp = (mpr->xsr) * vtemp ; // multiplication by xi/r in the nucleus
288 // by 1/r in the shells
289 // by r(xi-1) in the ZEC
290
291 // In the ZEC, a multiplication by r has been performed instead
292 // of the division:
293 temp.set_dzpuis( temp.get_dzpuis() + 2 ) ;
294
295 source = bbb * source + 0.5 * temp ;
296
297 }
298 else{
299 source = - 0.5 * ( dlogn(1)*dlogn(1) + dlogn(2)*dlogn(2) + dlogn(3)*dlogn(3) ) ;
300 }
301
302 source.std_spectral_base() ;
303
304 double int_grav = source.integrale() ;
305
306 // Matter term
307 // -----------
308
309 if (relativistic) {
310 source = qpig * a_car * bbb * s_euler ;
311 }
312 else{
313 source = qpig * ( 3 * press + nbar * uuu * uuu ) ;
314 }
315
316 source.std_spectral_base() ;
317
318 double int_mat = source.integrale() ;
319
320 // Virial error
321 // ------------
322 if (ost != 0x0) {
323 *ost << "Star_rot::grv3 : gravitational term : " << int_grav
324 << endl ;
325 *ost << "Star_rot::grv3 : matter term : " << int_mat
326 << endl ;
327 }
328
329 p_grv3 = new double( (int_grav + int_mat) / int_mat ) ;
330
331 }
332
333 return *p_grv3 ;
334
335}
336
337
338 //----------------------------//
339 // R_circ //
340 //----------------------------//
341
342double Star_rot::r_circ() const {
343
344 if (p_r_circ == 0x0) { // a new computation is required
345
346 // Index of the point at phi=0, theta=pi/2 at the surface of the star:
347 const Mg3d* mg = mp.get_mg() ;
348 assert(mg->get_type_t() == SYM) ;
349 int l_b = nzet - 1 ;
350 int i_b = mg->get_nr(l_b) - 1 ;
351 int j_b = mg->get_nt(l_b) - 1 ;
352 int k_b = 0 ;
353
354 p_r_circ = new double( bbb.val_grid_point(l_b, k_b, j_b, i_b) * ray_eq() ) ;
355
356 }
357
358 return *p_r_circ ;
359
360}
361
362 //----------------------------//
363 // Surface area //
364 //----------------------------//
365
366 double Star_rot::area() const {
367
368 if (p_area == 0x0) { // a new computation is required
369 const Mg3d& mg = *(mp.get_mg()) ;
370 int np = mg.get_np(0) ;
371 int nt = mg.get_nt(0) ;
372 assert(np == 1) ; //Only valid for axisymmetric configurations
373
374 const Map_radial* mp_rad = dynamic_cast<const Map_radial*>(&mp) ;
375 assert(mp_rad != 0x0) ;
376
377 Valeur va_r(mg.get_angu()) ;
378 va_r.annule_hard() ;
379 Itbl lsurf = l_surf() ;
380 Tbl xisurf = xi_surf() ;
381 for (int k=0; k<np; k++) {
382 for (int j=0; j<nt; j++) {
383 int l_star = lsurf(k, j) ;
384 double xi_star = xisurf(k, j) ;
385 va_r.set(0, k, j, 0) = mp_rad->val_r_jk(l_star, xi_star, j, k) ;
386 }
387 }
388 va_r.std_base_scal() ;
389
390 Valeur integ(mg.get_angu()) ;
391 Valeur dr = va_r.dsdt() ;
392 integ = sqrt(va_r*va_r + dr*dr) ;
393 Valeur a2 = get_a_car().get_spectral_va() ; a2.std_base_scal() ;
394 Valeur b = get_bbb().get_spectral_va() ; b.std_base_scal() ;
395 for (int k=0; k<np; k++) {
396 for (int j=0; j<nt; j++) {
397 integ.set(0, k, j, 0) *= sqrt(a2.val_point_jk(lsurf(k, j), xisurf(k, j), j, k))
398 * b.val_point_jk(lsurf(k, j), xisurf(k, j), j, k) * va_r(0, k, j, 0) ;
399 }
400 }
401 integ.std_base_scal() ;
402 Valeur integ2 = integ.mult_st() ;
403 double surftot = 0. ;
404 for (int j=0; j<nt; j++) {
405 surftot += (*integ2.c_cf)(0, 0, j, 0) / double(2*j+1) ;
406 }
407
408 p_area = new double( 4*M_PI*surftot) ;
409
410 }
411
412 return *p_area ;
413
414}
415
416 double Star_rot::mean_radius() const {
417
418 return sqrt(area()/(4*M_PI)) ;
419
420 }
421
422 //----------------------------//
423 // Flattening //
424 //----------------------------//
425
426double Star_rot::aplat() const {
427
428 if (p_aplat == 0x0) { // a new computation is required
429
430 p_aplat = new double( ray_pole() / ray_eq() ) ;
431
432 }
433
434 return *p_aplat ;
435
436}
437
438
439 //----------------------------//
440 // Z_eq_f //
441 //----------------------------//
442
443double Star_rot::z_eqf() const {
444
445 if (p_z_eqf == 0x0) { // a new computation is required
446
447 // Index of the point at phi=0, theta=pi/2 at the surface of the star:
448 const Mg3d* mg = mp.get_mg() ;
449 assert(mg->get_type_t() == SYM) ;
450 int l_b = nzet - 1 ;
451 int i_b = mg->get_nr(l_b) - 1 ;
452 int j_b = mg->get_nt(l_b) - 1 ;
453 int k_b = 0 ;
454
455 double u_eq = uuu.val_grid_point(l_b, k_b, j_b, i_b) ;
456 double n_eq = nn.val_grid_point(l_b, k_b, j_b, i_b) ;
457 double nphi_eq = nphi.val_grid_point(l_b, k_b, j_b, i_b) ;
458
459 p_z_eqf = new double( sqrt((1.-u_eq)/(1.+u_eq))
460 / (n_eq + nphi_eq * r_circ() )
461 - 1. ) ;
462 }
463
464 return *p_z_eqf ;
465
466}
467 //----------------------------//
468 // Z_eq_b //
469 //----------------------------//
470
471double Star_rot::z_eqb() const {
472
473 if (p_z_eqb == 0x0) { // a new computation is required
474
475 // Index of the point at phi=0, theta=pi/2 at the surface of the star:
476 const Mg3d* mg = mp.get_mg() ;
477 assert(mg->get_type_t() == SYM) ;
478 int l_b = nzet - 1 ;
479 int i_b = mg->get_nr(l_b) - 1 ;
480 int j_b = mg->get_nt(l_b) - 1 ;
481 int k_b = 0 ;
482
483 double u_eq = uuu.val_grid_point(l_b, k_b, j_b, i_b) ;
484 double n_eq = nn.val_grid_point(l_b, k_b, j_b, i_b) ;
485 double nphi_eq = nphi.val_grid_point(l_b, k_b, j_b, i_b) ;
486
487 p_z_eqb = new double( sqrt((1.+u_eq)/(1.-u_eq))
488 / (n_eq - nphi_eq * r_circ() )
489 - 1. ) ;
490
491 }
492
493 return *p_z_eqb ;
494
495}
496
497
498 //----------------------------//
499 // Z_pole //
500 //----------------------------//
501
502double Star_rot::z_pole() const {
503
504 if (p_z_pole == 0x0) { // a new computation is required
505
506 double n_pole = nn.val_point(ray_pole(), 0., 0.) ;
507
508 p_z_pole = new double( 1. / n_pole - 1. ) ;
509
510 }
511
512 return *p_z_pole ;
513
514}
515
516
517 //----------------------------//
518 // Quadrupole moment //
519 //----------------------------//
520
521double Star_rot::mom_quad() const {
522
523 using namespace Unites ;
524 if (p_mom_quad == 0x0) { // a new computation is required
525
526 // Source for of the Poisson equation for nu
527 // -----------------------------------------
528
529 Scalar source(mp) ;
530
531 if (relativistic) {
532 Scalar bet = log(bbb) ;
533 bet.std_spectral_base() ;
534
536 Vector dlogn_bet = dlogn + bet.derive_cov( mp.flat_met_spher() ) ;
537
538 source = qpig * a_car *( ener_euler + s_euler ) + ak_car
539 - dlogn(1)*dlogn_bet(1) - dlogn(2)*dlogn_bet(2) - dlogn(3)*dlogn_bet(3) ;
540 }
541 else {
542 source = qpig * nbar ;
543 }
544
545 source.std_spectral_base() ;
546
547 // Multiplication by -r^2 P_2(cos(theta))
548 // [cf Eq.(7) of Salgado et al. Astron. Astrophys. 291, 155 (1994) ]
549 // ------------------------------------------------------------------
550
551 // Multiplication by r^2 :
552 // ----------------------
553 source.mult_r() ;
554 source.mult_r() ;
555 if (source.check_dzpuis(2)) {
556 source.inc_dzpuis(2) ;
557 }
558
559 // Muliplication by cos^2(theta) :
560 // -----------------------------
561 Scalar temp = source ;
562
563 // What follows is valid only for a mapping of class Map_radial :
564 assert( dynamic_cast<const Map_radial*>(&mp) != 0x0 ) ;
565
566 if (temp.get_etat() == ETATQCQ) {
567 Valeur& vtemp = temp.set_spectral_va() ;
568 vtemp = vtemp.mult_ct() ; // multiplication by cos(theta)
569 vtemp = vtemp.mult_ct() ; // multiplication by cos(theta)
570 }
571
572 // Muliplication by -P_2(cos(theta)) :
573 // ----------------------------------
574 source = 0.5 * source - 1.5 * temp ;
575
576 // Final result
577 // ------------
578
579 p_mom_quad = new double( source.integrale() / qpig ) ;
580
581 }
582
583 return *p_mom_quad ;
584
585}
586
587
588
589
590}
Time evolution with partial storage (*** under development ***).
Definition evolution.h:371
Basic integer array class.
Definition itbl.h:122
Base class for pure radial mappings.
Definition map.h:1536
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
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
Tensor field of valence 0 (or component of a tensorial field).
Definition scalar.h:387
const Vector & derive_cov(const Metric &gam) const
Returns the gradient (1-form = covariant vector) of *this
double integrale() const
Computes the integral over all space of *this .
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
const Valeur & get_spectral_va() const
Returns va (read only version)
Definition scalar.h:601
int get_etat() const
Returns the logical state ETATNONDEF (undefined), ETATZERO (null) or ETATQCQ (ordinary).
Definition scalar.h:554
double val_point(double r, double theta, double phi) const
Computes the value of the field at an arbitrary point , by means of the spectral expansion.
Definition scalar.C:890
double * p_angu_mom
Angular momentum.
Definition star_rot.h:232
virtual double mean_radius() const
Mean star radius from the area .
double * p_aplat
Flatening r_pole/r_eq.
Definition star_rot.h:237
double * p_mom_quad
Quadrupole moment.
Definition star_rot.h:242
virtual double mom_quad() const
Quadrupole moment.
double * p_z_eqf
Forward redshift factor at equator.
Definition star_rot.h:239
Scalar tnphi
Component of the shift vector.
Definition star_rot.h:118
virtual double z_pole() const
Redshift factor at North pole.
virtual double z_eqb() const
Backward redshift factor at equator.
Scalar b_car
Square of the metric factor B.
Definition star_rot.h:110
Scalar bbb
Metric factor B.
Definition star_rot.h:107
virtual double mass_b() const
Baryon mass.
double omega
Rotation angular velocity ([f_unit] )
Definition star_rot.h:101
Scalar ak_car
Scalar .
Definition star_rot.h:186
virtual double r_circ() const
Circumferential radius.
Scalar uuu
Norm of u_euler.
Definition star_rot.h:121
virtual double z_eqf() const
Forward redshift factor at equator.
Scalar nphi
Metric coefficient .
Definition star_rot.h:113
virtual double mass_g() const
Gravitational mass.
double * p_area
Integrated surface area.
Definition star_rot.h:238
virtual double aplat() const
Flatening r_pole/r_eq.
double * p_grv3
Error on the virial identity GRV3.
Definition star_rot.h:235
const Scalar & get_a_car() const
Returns the square of the metric factor A.
Definition star_rot.h:319
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_tsw
Ratio T/W.
Definition star_rot.h:233
bool relativistic
Indicator of relativity: true for a relativistic star, false for a Newtonian one.
Definition star_rot.h:94
virtual double tsw() const
Ratio T/W.
virtual double grv3(ostream *ost=0x0) const
Error on the virial identity GRV3.
Scalar dzeta
Metric potential .
Definition star_rot.h:134
virtual double area() const
Integrated surface area in .
Scalar a_car
Square of the metric factor A.
Definition star_rot.h:104
double * p_grv2
Error on the virial identity GRV2.
Definition star_rot.h:234
const Scalar & get_bbb() const
Returns the metric factor B.
Definition star_rot.h:316
double * p_r_circ
Circumferential radius.
Definition star_rot.h:236
virtual double grv2() const
Error on the virial identity GRV2.
virtual double angu_mom() const
Angular momentum.
double * p_z_eqb
Backward redshift factor at equator.
Definition star_rot.h:240
double * p_z_pole
Redshift factor at North pole.
Definition star_rot.h:241
static double lambda_grv2(const Scalar &sou_m, const Scalar &sou_q)
Computes the coefficient which ensures that the GRV2 virial identity is satisfied.
Scalar ener
Total energy density in the fluid frame.
Definition star.h:193
Scalar logn
Logarithm of the lapse N .
Definition star.h:222
Scalar nn
Lapse function N .
Definition star.h:225
double * p_mass_b
Baryon mass.
Definition star.h:268
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
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
Scalar gam_euler
Lorentz factor between the fluid and Eulerian observers.
Definition star.h:204
double ray_eq() const
Coordinate radius at , [r_unit].
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
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
Scalar press
Fluid pressure.
Definition star.h:194
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 ray_pole() const
Coordinate radius at [r_unit].
Basic array class.
Definition tbl.h:161
Values and coefficients of a (real-value) function.
Definition valeur.h:287
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
void std_base_scal()
Sets the bases for spectral expansions (member base ) to the standard ones for a scalar.
Definition valeur.C:824
Tensor field of valence 1.
Definition vector.h:188
Cmp sqrt(const Cmp &)
Square root.
Definition cmp_math.C:220
Cmp log(const Cmp &)
Neperian logarithm.
Definition cmp_math.C:296
Lorene prototypes.
Definition app_hor.h:64
Standard units of space, time and mass.