LORENE
star_QI_global.C
1 /*
2  * Methods of the class Star_QI to compute global quantities
3  *
4  * (see file compobj.h for documentation).
5  *
6  */
7 
8 /*
9  * Copyright (c) 2012 Claire Some, Eric Gourgoulhon
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 version 2
15  * as published by the Free Software Foundation.
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 
29 
30 /*
31  * $Id: star_QI_global.C,v 1.4 2016/12/05 16:17:49 j_novak Exp $
32  * $Log: star_QI_global.C,v $
33  * Revision 1.4 2016/12/05 16:17:49 j_novak
34  * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
35  *
36  * Revision 1.3 2014/10/13 08:52:49 j_novak
37  * Lorene classes and functions now belong to the namespace Lorene.
38  *
39  * Revision 1.2 2013/06/05 15:10:41 j_novak
40  * Suppression of FINJAC sampling in r. This Jacobi(0,2) base is now
41  * available by setting colloc_r to BASE_JAC02 in the Mg3d constructor.
42  *
43  * Revision 1.1 2012/11/21 14:54:13 c_some
44  * First version
45  *
46  *
47  *
48  * $Header: /cvsroot/Lorene/C++/Source/Compobj/star_QI_global.C,v 1.4 2016/12/05 16:17:49 j_novak Exp $
49  *
50  */
51 
52 
53 // C headers
54 #include <cassert>
55 #include <cstdlib>
56 
57 // Lorene headers
58 #include "compobj.h"
59 #include "unites.h"
60 #include "proto_f77.h"
61 
62  //------------------------//
63  // Gravitational mass //
64  //------------------------//
65 
66 namespace Lorene {
67 double Star_QI::mass_g() const {
68 
69  if (p_mass_g == 0x0) { // a new computation is required
70 
71  Scalar s_euler = stress_euler.trace(gamma) ;
72 
73  //## alternative:
74  // assert(*(stress_euler.get_triad()) == mp.get_bvect_spher()) ;
75  // Scalar s_euler = ( stress_euler(1,1) + stress_euler(2,2) ) / a_car
76  // + stress_euler(3,3) / b_car ;
77 
78  // Cf. Eq. (4.18) of arXiv:1003.5015v2 with (E+p) U = B * mom_euler(3)
79 
80  Scalar source = nn * (ener_euler + s_euler)
81  + 2 * b_car * mom_euler(3)
82  * tnphi ;
83  source = a_car * bbb * source ;
84  source.std_spectral_base() ;
85 
86  p_mass_g = new double( source.integrale() ) ;
87 
88  }
89 
90  return *p_mass_g ;
91 
92 }
93 
94 
95  //------------------------//
96  // Angular momentum //
97  //------------------------//
98 
99 double Star_QI::angu_mom() const {
100 
101  if (p_angu_mom == 0x0) { // a new computation is required
102 
103  // Cf. Eq. (4.39) of arXiv:1003.5015v2 with (E+p) U = B * mom_euler(3)
104 
105  assert(*(mom_euler.get_triad()) == mp.get_bvect_spher()) ;
106 
107  Scalar dens = mom_euler(3) ;
108 
109  dens.mult_r() ; // Multiplication by
110  dens.set_spectral_va() = (dens.get_spectral_va()).mult_st() ; // r sin(theta)
111 
112  dens = a_car * b_car * bbb * dens ;
113 
114  p_angu_mom = new double( dens.integrale() ) ;
115 
116  }
117 
118  return *p_angu_mom ;
119 
120 }
121 
122  //--------------------//
123  // GRV2 //
124  //--------------------//
125 
126 double Star_QI::grv2() const {
127 
128  using namespace Unites ;
129  if (p_grv2 == 0x0) { // a new computation is required
130 
131  assert( *(stress_euler.get_triad()) == mp.get_bvect_spher() ) ;
132  Scalar sou_m = 2 * qpig * a_car * b_car * stress_euler(3,3) ;
133 
134  Vector dlogn = logn.derive_cov( mp.flat_met_spher() ) ;
135  Scalar sou_q = 1.5 * ak_car
136  - dlogn(1)*dlogn(1) - dlogn(2)*dlogn(2) - dlogn(3)*dlogn(3) ;
137 
138  p_grv2 = new double( double(1) - lambda_grv2(sou_m, sou_q) ) ;
139 
140  }
141 
142  return *p_grv2 ;
143 
144 }
145 
146 
147  //--------------------//
148  // GRV3 //
149  //--------------------//
150 
151 double Star_QI::grv3(ostream* ost) const {
152 
153  using namespace Unites ;
154 
155  if (p_grv3 == 0x0) { // a new computation is required
156 
157  Scalar source(mp) ;
158 
159  // Gravitational term [cf. Eq. (43) of Gourgoulhon & Bonazzola
160  // ------------------ Class. Quantum Grav. 11, 443 (1994)]
161 
162  Vector dlogn = logn.derive_cov( mp.flat_met_spher() ) ;
163 
164  Scalar alpha = dzeta - logn ;
165  Scalar bet = log( bbb ) ;
166  bet.std_spectral_base() ;
167 
168  Vector dalpha = alpha.derive_cov( mp.flat_met_spher() ) ;
169  Vector dbet = bet.derive_cov( mp.flat_met_spher() ) ;
170 
171  source = 0.75 * ak_car
172  - dlogn(1)*dlogn(1) - dlogn(2)*dlogn(2) - dlogn(3)*dlogn(3)
173  + 0.5 * ( dalpha(1)*dbet(1) + dalpha(2)*dbet(2) + dalpha(3)*dbet(3) ) ;
174 
175  Scalar aa = alpha - 0.5 * bet ;
176  Scalar daadt = aa.srdsdt() ; // 1/r d/dth
177 
178  // What follows is valid only for a mapping of class Map_radial :
179  const Map_radial* mpr = dynamic_cast<const Map_radial*>(&mp) ;
180  if (mpr == 0x0) {
181  cout << "Star_rot::grv3: the mapping does not belong"
182  << " to the class Map_radial !" << endl ;
183  abort() ;
184  }
185 
186  // Computation of 1/tan(theta) * 1/r daa/dtheta
187  if (daadt.get_etat() == ETATQCQ) {
188  Valeur& vdaadt = daadt.set_spectral_va() ;
189  vdaadt = vdaadt.ssint() ; // division by sin(theta)
190  vdaadt = vdaadt.mult_ct() ; // multiplication by cos(theta)
191  }
192 
193  Scalar temp = aa.dsdr() + daadt ;
194  temp = ( bbb - a_car/bbb ) * temp ;
195  temp.std_spectral_base() ;
196 
197  // Division by r
198  Valeur& vtemp = temp.set_spectral_va() ;
199  vtemp = vtemp.sx() ; // division by xi in the nucleus
200  // Id in the shells
201  // division by xi-1 in the ZEC
202  vtemp = (mpr->xsr) * vtemp ; // multiplication by xi/r in the nucleus
203  // by 1/r in the shells
204  // by r(xi-1) in the ZEC
205 
206  // In the ZEC, a multiplication by r has been performed instead
207  // of the division:
208  temp.set_dzpuis( temp.get_dzpuis() + 2 ) ;
209 
210  source = bbb * source + 0.5 * temp ;
211 
212  source.std_spectral_base() ;
213 
214  double int_grav = source.integrale() ;
215 
216  // Matter term
217  // -----------
218 
219  Scalar s_euler = stress_euler.trace(gamma) ;
220 
221  //## alternative:
222  // assert(*(stress_euler.get_triad()) == mp.get_bvect_spher()) ;
223  // Scalar s_euler = ( stress_euler(1,1) + stress_euler(2,2) ) / a_car
224  // + stress_euler(3,3) / b_car ;
225 
226  source = qpig * a_car * bbb * s_euler ;
227 
228  source.std_spectral_base() ;
229 
230  double int_mat = source.integrale() ;
231 
232  // Virial error
233  // ------------
234  if (ost != 0x0) {
235  *ost << "Star_rot::grv3 : gravitational term : " << int_grav
236  << endl ;
237  *ost << "Star_rot::grv3 : matter term : " << int_mat
238  << endl ;
239  }
240 
241  p_grv3 = new double( (int_grav + int_mat) / int_mat ) ;
242 
243  }
244 
245  return *p_grv3 ;
246 
247 }
248 
249  //----------------------------//
250  // Quadrupole moment //
251  //----------------------------//
252 
253 double Star_QI::mom_quad() const {
254 
255  using namespace Unites ;
256  if (p_mom_quad == 0x0) { // a new computation is required
257 
258  // Source for of the Poisson equation for nu
259  // -----------------------------------------
260 
261  Scalar source(mp) ;
262 
263  Scalar s_euler = stress_euler.trace(gamma) ;
264 
265  //## alternative:
266  // assert(*(stress_euler.get_triad()) == mp.get_bvect_spher()) ;
267  // Scalar s_euler = ( stress_euler(1,1) + stress_euler(2,2) ) / a_car
268  // + stress_euler(3,3) / b_car ;
269 
270  Scalar bet = log(bbb) ;
271  bet.std_spectral_base() ;
272 
273  Vector dlogn = logn.derive_cov( mp.flat_met_spher() ) ;
274  Vector dlogn_bet = dlogn + bet.derive_cov( mp.flat_met_spher() ) ;
275 
276  source = qpig * a_car *( ener_euler + s_euler ) + ak_car
277  - dlogn(1)*dlogn_bet(1) - dlogn(2)*dlogn_bet(2) - dlogn(3)*dlogn_bet(3) ;
278 
279  source.std_spectral_base() ;
280 
281  // Multiplication by -r^2 P_2(cos(theta))
282  // [cf Eq.(7) of Salgado et al. Astron. Astrophys. 291, 155 (1994) ]
283  // ------------------------------------------------------------------
284 
285  // Multiplication by r^2 :
286  // ----------------------
287  source.mult_r() ;
288  source.mult_r() ;
289  if (source.check_dzpuis(2)) {
290  source.inc_dzpuis(2) ;
291  }
292 
293  // Muliplication by cos^2(theta) :
294  // -----------------------------
295  Scalar temp = source ;
296 
297  // What follows is valid only for a mapping of class Map_radial :
298  assert( dynamic_cast<const Map_radial*>(&mp) != 0x0 ) ;
299 
300  if (temp.get_etat() == ETATQCQ) {
301  Valeur& vtemp = temp.set_spectral_va() ;
302  vtemp = vtemp.mult_ct() ; // multiplication by cos(theta)
303  vtemp = vtemp.mult_ct() ; // multiplication by cos(theta)
304  }
305 
306  // Muliplication by -P_2(cos(theta)) :
307  // ----------------------------------
308  source = 0.5 * source - 1.5 * temp ;
309 
310  // Final result
311  // ------------
312 
313  p_mom_quad = new double( source.integrale() / qpig ) ;
314 
315  }
316 
317  return *p_mom_quad ;
318 
319 }
320 
321 
322 // Function Star_QI::lambda_grv2
323 
324 double Star_QI::lambda_grv2(const Scalar& sou_m, const Scalar& sou_q) {
325 
326  const Map_radial* mprad = dynamic_cast<const Map_radial*>( &sou_m.get_mp() ) ;
327 
328  if (mprad == 0x0) {
329  cout << "Star_rot::lambda_grv2: the mapping of sou_m does not"
330  << endl << " belong to the class Map_radial !" << endl ;
331  abort() ;
332  }
333 
334  assert( &sou_q.get_mp() == mprad ) ;
335 
336  sou_q.check_dzpuis(4) ;
337 
338  const Mg3d* mg = mprad->get_mg() ;
339  int nz = mg->get_nzone() ;
340 
341  // Construction of a Map_af which coincides with *mp on the equator
342  // ----------------------------------------------------------------
343 
344  double theta0 = M_PI / 2 ; // Equator
345  double phi0 = 0 ;
346 
347  Map_af mpaff(*mprad) ;
348 
349  for (int l=0 ; l<nz ; l++) {
350  double rmax = mprad->val_r(l, double(1), theta0, phi0) ;
351  switch ( mg->get_type_r(l) ) {
352  case RARE: {
353  double rmin = mprad->val_r(l, double(0), theta0, phi0) ;
354  mpaff.set_alpha(rmax - rmin, l) ;
355  mpaff.set_beta(rmin, l) ;
356  break ;
357  }
358 
359  case FIN: {
360  double rmin = mprad->val_r(l, double(-1), theta0, phi0) ;
361  mpaff.set_alpha( double(.5) * (rmax - rmin), l ) ;
362  mpaff.set_beta( double(.5) * (rmax + rmin), l) ;
363  break ;
364  }
365 
366  case UNSURR: {
367  double rmin = mprad->val_r(l, double(-1), theta0, phi0) ;
368  double umax = double(1) / rmin ;
369  double umin = double(1) / rmax ;
370  mpaff.set_alpha( double(.5) * (umin - umax), l) ;
371  mpaff.set_beta( double(.5) * (umin + umax), l) ;
372  break ;
373  }
374 
375  default: {
376  cout << "Star_rot::lambda_grv2: unknown type_r ! " << endl ;
377  abort () ;
378  break ;
379  }
380 
381  }
382  }
383 
384 
385  // Reduced Jacobian of
386  // the transformation (r,theta,phi) <-> (dzeta,theta',phi')
387  // ------------------------------------------------------------
388 
389  Mtbl jac = 1 / ( (mprad->xsr) * (mprad->dxdr) ) ;
390  // R/x dR/dx in the nucleus
391  // R dR/dx in the shells
392  // - U/(x-1) dU/dx in the ZEC
393  for (int l=0; l<nz; l++) {
394  switch ( mg->get_type_r(l) ) {
395  case RARE: {
396  double a1 = mpaff.get_alpha()[l] ;
397  *(jac.t[l]) = *(jac.t[l]) / (a1*a1) ;
398  break ;
399  }
400 
401  case FIN: {
402  double a1 = mpaff.get_alpha()[l] ;
403  double b1 = mpaff.get_beta()[l] ;
404  assert( jac.t[l]->get_etat() == ETATQCQ ) ;
405  double* tjac = jac.t[l]->t ;
406  double* const xi = mg->get_grille3d(l)->x ;
407  for (int k=0; k<mg->get_np(l); k++) {
408  for (int j=0; j<mg->get_nt(l); j++) {
409  for (int i=0; i<mg->get_nr(l); i++) {
410  *tjac = *tjac /
411  (a1 * (a1 * xi[i] + b1) ) ;
412  tjac++ ;
413  }
414  }
415  }
416 
417  break ;
418  }
419 
420  case UNSURR: {
421  double a1 = mpaff.get_alpha()[l] ;
422  *(jac.t[l]) = - *(jac.t[l]) / (a1*a1) ;
423  break ;
424  }
425 
426  default: {
427  cout << "Star_rot::lambda_grv2: unknown type_r ! " << endl ;
428  abort () ;
429  break ;
430  }
431 
432  }
433 
434  }
435 
436 
437  // Multiplication of the sources by the reduced Jacobian:
438  // -----------------------------------------------------
439 
440  Mtbl s_m(mg) ;
441  if ( sou_m.get_etat() == ETATZERO ) {
442  s_m = 0 ;
443  }
444  else{
445  assert(sou_m.get_spectral_va().get_etat() == ETATQCQ) ;
446  sou_m.get_spectral_va().coef_i() ;
447  s_m = *(sou_m.get_spectral_va().c) ;
448  }
449 
450  Mtbl s_q(mg) ;
451  if ( sou_q.get_etat() == ETATZERO ) {
452  s_q = 0 ;
453  }
454  else{
455  assert(sou_q.get_spectral_va().get_etat() == ETATQCQ) ;
456  sou_q.get_spectral_va().coef_i() ;
457  s_q = *(sou_q.get_spectral_va().c) ;
458  }
459 
460  s_m *= jac ;
461  s_q *= jac ;
462 
463 
464  // Preparations for the call to the Fortran subroutine
465  // ---------------------------------------------------
466 
467  int np1 = 1 ; // Axisymmetry enforced
468  int nt = mg->get_nt(0) ;
469  int nt2 = 2*nt - 1 ; // Number of points for the theta sampling
470  // in [0,Pi], instead of [0,Pi/2]
471 
472  // Array NDL
473  // ---------
474  int* ndl = new int[nz+4] ;
475  ndl[0] = nz ;
476  for (int l=0; l<nz; l++) {
477  ndl[1+l] = mg->get_nr(l) ;
478  }
479  ndl[1+nz] = nt2 ;
480  ndl[2+nz] = np1 ;
481  ndl[3+nz] = nz ;
482 
483  // Parameters NDR, NDT, NDP
484  // ------------------------
485  int nrmax = 0 ;
486  for (int l=0; l<nz ; l++) {
487  nrmax = ( ndl[1+l] > nrmax ) ? ndl[1+l] : nrmax ;
488  }
489  int ndr = nrmax + 5 ;
490  int ndt = nt2 + 2 ;
491  int ndp = np1 + 2 ;
492 
493  // Array ERRE
494  // ----------
495 
496  double* erre = new double [nz*ndr] ;
497 
498  for (int l=0; l<nz; l++) {
499  double a1 = mpaff.get_alpha()[l] ;
500  double b1 = mpaff.get_beta()[l] ;
501  for (int i=0; i<ndl[1+l]; i++) {
502  double xi = mg->get_grille3d(l)->x[i] ;
503  erre[ ndr*l + i ] = a1 * xi + b1 ;
504  }
505  }
506 
507  // Arrays containing the data
508  // --------------------------
509 
510  int ndrt = ndr*ndt ;
511  int ndrtp = ndr*ndt*ndp ;
512  int taille = ndrtp*nz ;
513 
514  double* tsou_m = new double[ taille ] ;
515  double* tsou_q = new double[ taille ] ;
516 
517  // Initialisation to zero :
518  for (int i=0; i<taille; i++) {
519  tsou_m[i] = 0 ;
520  tsou_q[i] = 0 ;
521  }
522 
523  // Copy of s_m into tsou_m
524  // -----------------------
525 
526  for (int l=0; l<nz; l++) {
527  for (int k=0; k<np1; k++) {
528  for (int j=0; j<nt; j++) {
529  for (int i=0; i<mg->get_nr(l); i++) {
530  double xx = s_m(l, k, j, i) ;
531  tsou_m[ndrtp*l + ndrt*k + ndr*j + i] = xx ;
532  // point symetrique par rapport au plan theta = pi/2 :
533  tsou_m[ndrtp*l + ndrt*k + ndr*(nt2-1-j) + i] = xx ;
534  }
535  }
536  }
537  }
538 
539  // Copy of s_q into tsou_q
540  // -----------------------
541 
542  for (int l=0; l<nz; l++) {
543  for (int k=0; k<np1; k++) {
544  for (int j=0; j<nt; j++) {
545  for (int i=0; i<mg->get_nr(l); i++) {
546  double xx = s_q(l, k, j, i) ;
547  tsou_q[ndrtp*l + ndrt*k + ndr*j + i] = xx ;
548  // point symetrique par rapport au plan theta = pi/2 :
549  tsou_q[ndrtp*l + ndrt*k + ndr*(nt2-1-j) + i] = xx ;
550  }
551  }
552  }
553  }
554 
555 
556  // Computation of the integrals
557  // ----------------------------
558 
559  double int_m, int_q ;
560  F77_integrale2d(ndl, &ndr, &ndt, &ndp, erre, tsou_m, &int_m) ;
561  F77_integrale2d(ndl, &ndr, &ndt, &ndp, erre, tsou_q, &int_q) ;
562 
563  // Cleaning
564  // --------
565 
566  delete [] ndl ;
567  delete [] erre ;
568  delete [] tsou_m ;
569  delete [] tsou_q ;
570 
571  // Computation of lambda
572  // ---------------------
573 
574  double lambda ;
575  if ( int_q != double(0) ) {
576  lambda = - int_m / int_q ;
577  }
578  else{
579  lambda = 0 ;
580  }
581 
582  return lambda ;
583 
584 }
585 
586 }
Scalar logn
Logarithm of the lapse N .
Definition: compobj.h:498
const Grille3d * get_grille3d(int l) const
Returns a pointer on the 3D mono-grid for domain no. l.
Definition: grilles.h:517
Cmp log(const Cmp &)
Neperian logarithm.
Definition: cmp_math.C:299
Vector mom_euler
Total 3-momentum density in the Eulerian frame.
Definition: compobj.h:150
const double * get_alpha() const
Returns the pointer on the array alpha.
Definition: map_af.C:607
void mult_r()
Multiplication by r everywhere; dzpuis is not changed.
Scalar tnphi
Component of the shift vector.
Definition: compobj.h:503
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
Definition: grilles.h:479
Multi-domain array.
Definition: mtbl.h:118
Scalar ak_car
Scalar .
Definition: compobj.h:318
const Base_vect_spher & get_bvect_spher() const
Returns the orthonormal vectorial basis associated with the coordinates of the mapping.
Definition: map.h:801
Lorene prototypes.
Definition: app_hor.h:67
Standard units of space, time and mass.
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:783
Tensor field of valence 0 (or component of a tensorial field).
Definition: scalar.h:393
void coef_i() const
Computes the physical value of *this.
virtual double mass_g() const
Gravitational mass.
double integrale() const
Computes the integral over all space of *this .
Definition: scalar_integ.C:64
virtual double mom_quad() const
Quadrupole moment.
const Valeur & sx() const
Returns (r -sampling = RARE ) \ Id (r sampling = FIN ) \ (r -sampling = UNSURR ) ...
Definition: valeur_sx.C:113
virtual void std_spectral_base()
Sets the spectral bases of the Valeur va to the standard ones for a scalar field. ...
Definition: scalar.C:790
Values and coefficients of a (real-value) function.
Definition: valeur.h:297
int get_etat() const
Gives the logical state.
Definition: tbl.h:414
int get_etat() const
Returns the logical state ETATNONDEF (undefined), ETATZERO (null) or ETATQCQ (ordinary).
Definition: scalar.h:560
double * x
Array of values of at the nr collocation points.
Definition: grilles.h:215
Tensor field of valence 1.
Definition: vector.h:188
virtual double grv3(ostream *ost=0x0) const
Error on the virial identity GRV3.
Scalar a_car
Square of the metric factor A.
Definition: compobj.h:290
double * p_grv2
Error on the virial identity GRV2.
Definition: compobj.h:588
void set_dzpuis(int)
Modifies the dzpuis flag.
Definition: scalar.C:814
Scalar b_car
Square of the metric factor B.
Definition: compobj.h:296
void set_beta(double beta0, int l)
Modifies the value of in domain no. l.
Definition: map_af.C:771
int get_etat() const
Returns the logical state.
Definition: valeur.h:760
double * p_angu_mom
Angular momentum.
Definition: compobj.h:324
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 Base_vect * get_triad() const
Returns the vectorial basis (triad) on which the components are defined.
Definition: tensor.h:879
Scalar bbb
Metric factor B.
Definition: compobj.h:293
virtual void inc_dzpuis(int inc=1)
Increases by inc units the value of dzpuis and changes accordingly the values of the Scalar in the co...
static double lambda_grv2(const Scalar &sou_m, const Scalar &sou_q)
Computes the coefficient which ensures that the GRV2 virial identity is satisfied.
const Valeur & ssint() const
Returns of *this.
Definition: valeur_ssint.C:115
Coord dxdr
in the nucleus and in the non-compactified shells; \ in the compactified outer domain.
Definition: map.h:1581
int get_dzpuis() const
Returns dzpuis.
Definition: scalar.h:563
double * t
The array of double.
Definition: tbl.h:176
const double * get_beta() const
Returns the pointer on the array beta.
Definition: map_af.C:611
Mtbl * c
Values of the function at the points of the multi-grid.
Definition: valeur.h:309
Base class for pure radial mappings.
Definition: map.h:1557
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:465
void set_alpha(double alpha0, int l)
Modifies the value of in domain no. l.
Definition: map_af.C:760
Coord xsr
in the nucleus; \ 1/R in the non-compactified shells; \ in the compactified outer domain...
Definition: map.h:1570
Scalar dzeta
Metric potential .
Definition: compobj.h:516
double * p_grv3
Error on the virial identity GRV3.
Definition: compobj.h:589
double * p_mass_g
Gravitational mass (ADM mass as a volume integral)
Definition: compobj.h:591
double * p_mom_quad
Quadrupole moment.
Definition: compobj.h:590
int get_nr(int l) const
Returns the number of points in the radial direction ( ) in domain no. l.
Definition: grilles.h:469
Multi-domain grid.
Definition: grilles.h:279
Affine radial mapping.
Definition: map.h:2048
const Scalar & srdsdt() const
Returns of *this .
Definition: scalar_deriv.C:145
Scalar ener_euler
Total energy density E in the Eulerian frame.
Definition: compobj.h:147
Scalar nn
Lapse function N .
Definition: compobj.h:138
const Scalar & dsdr() const
Returns of *this .
Definition: scalar_deriv.C:113
const Valeur & mult_ct() const
Returns applied to *this.
Metric gamma
3-metric
Definition: compobj.h:144
virtual double grv2() const
Error on the virial identity GRV2.
virtual double angu_mom() const
Angular momentum.
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
Definition: grilles.h:474
Valeur & set_spectral_va()
Returns va (read/write version)
Definition: scalar.h:610
int get_type_r(int l) const
Returns the type of sampling in the radial direction in domain no.
Definition: grilles.h:491
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: scalar.C:879
Tbl ** t
Array (size nzone ) of pointers on the Tbl &#39;s.
Definition: mtbl.h:132
const Vector & derive_cov(const Metric &gam) const
Returns the gradient (1-form = covariant vector) of *this
Definition: scalar_deriv.C:390
const Map & get_mp() const
Returns the mapping.
Definition: tensor.h:874
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:324
const Valeur & get_spectral_va() const
Returns va (read only version)
Definition: scalar.h:607
Tensor trace(int ind1, int ind2) const
Trace on two different type indices.
Map & mp
Mapping describing the coordinate system (r,theta,phi)
Definition: compobj.h:135
Sym_tensor stress_euler
Stress tensor with respect to the Eulerian observer.
Definition: compobj.h:153