LORENE
et_bin_vel_pot.C
1/*
2 * Method of class Etoile_bin to compute the velocity scalar potential $\psi$
3 * by solving the continuity equation.
4 *
5 * (see file etoile.h for documentation).
6 *
7 */
8
9/*
10 * Copyright (c) 2000-2001 Eric Gourgoulhon
11 *
12 * This file is part of LORENE.
13 *
14 * LORENE is free software; you can redistribute it and/or modify
15 * it under the terms of the GNU General Public License as published by
16 * the Free Software Foundation; either version 2 of the License, or
17 * (at your option) any later version.
18 *
19 * LORENE is distributed in the hope that it will be useful,
20 * but WITHOUT ANY WARRANTY; without even the implied warranty of
21 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
22 * GNU General Public License for more details.
23 *
24 * You should have received a copy of the GNU General Public License
25 * along with LORENE; if not, write to the Free Software
26 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
27 *
28 */
29
30
31char et_bin_vel_pot_C[] = "$Header: /cvsroot/Lorene/C++/Source/Etoile/et_bin_vel_pot.C,v 1.15 2014/10/13 08:52:56 j_novak Exp $" ;
32
33/*
34 * $Id: et_bin_vel_pot.C,v 1.15 2014/10/13 08:52:56 j_novak Exp $
35 * $Log: et_bin_vel_pot.C,v $
36 * Revision 1.15 2014/10/13 08:52:56 j_novak
37 * Lorene classes and functions now belong to the namespace Lorene.
38 *
39 * Revision 1.14 2007/10/18 14:26:43 e_gourgoulhon
40 * Changed the call to Eos::der_nbar_ent in order to allow for MEos type
41 * of equation of state.
42 *
43 * Revision 1.13 2007/10/16 21:56:26 e_gourgoulhon
44 * Can deal with more than one domain into the star,
45 * thanks to the new function Map_radial::poisson_compact.
46 *
47 * Revision 1.12 2005/10/18 13:12:33 p_grandclement
48 * update of the mixted binary codes
49 *
50 * Revision 1.11 2004/05/25 15:38:38 f_limousin
51 * Minor modifs.
52 *
53 * Revision 1.10 2004/05/10 10:17:27 f_limousin
54 * Add a new member ssjm1_psi of class Etoile for the resolution of the
55 * oisson_interne equation
56 *
57 * Revision 1.9 2004/04/19 11:26:17 f_limousin
58 * Add a new function Etoile_bin::velocity_potential( , , , ) for the
59 * case of strange stars
60 *
61 * Revision 1.8 2004/04/08 17:02:00 f_limousin
62 * Modif to avoid an error in the compilation
63 *
64 * Revision 1.7 2004/04/08 16:52:58 f_limousin
65 * Minor change
66 *
67 * Revision 1.6 2004/04/08 16:36:36 f_limousin
68 * Implement the resolution of the continuity equation for strange
69 * stars.
70 *
71 * Revision 1.5 2003/10/24 11:43:57 e_gourgoulhon
72 * beta is now computed as ln(AN) in the case beta_auto
73 * is undefined (for instance, if the companion is black hole).
74 *
75 * Revision 1.4 2003/01/17 13:38:56 f_limousin
76 * Add comments
77 *
78 * Revision 1.3 2003/01/13 15:31:50 e_gourgoulhon
79 * Suppressed the desaliasing
80 * (did not worked due to missing basis in ylm).
81 *
82 * Revision 1.2 2002/12/10 14:44:21 k_taniguchi
83 * Change the multiplication "*" to "%"
84 * and flat_scalar_prod to flat_scalar_prod_desal.
85 *
86 * Revision 1.1.1.1 2001/11/20 15:19:28 e_gourgoulhon
87 * LORENE
88 *
89 * Revision 2.9 2001/02/23 15:18:59 eric
90 * Modification du calcul de zeta_h pour eviter division par zero
91 * dans les domaines externes a l'etoile.
92 *
93 * Revision 2.8 2001/02/07 09:47:42 eric
94 * zeta_h est desormais donne par Eos::der_nbar_ent.
95 *
96 * Revision 2.7 2000/12/22 13:10:03 eric
97 * Prolongement C^1 de dpsi en dehors de l'etoile.
98 *
99 * Revision 2.6 2000/03/22 12:56:44 eric
100 * Nouveau prototype d'Etoile_bin::velocity_potential : l'erreur est
101 * retournee en double.
102 *
103 * Revision 2.5 2000/02/25 17:35:29 eric
104 * Annulation de la source dans les zones externes avant l'appel a
105 * poisson_compact.
106 *
107 * Revision 2.4 2000/02/22 11:42:55 eric
108 * Test resolution de l'equation.
109 *
110 * Revision 2.3 2000/02/22 10:42:25 eric
111 * Correction erreur dans les termes sources: multiplication par unsurc2 de
112 * termes relativistes.
113 *
114 * Revision 2.2 2000/02/21 15:05:50 eric
115 * Traitement du cas psi0 = 0 .
116 *
117 * Revision 2.1 2000/02/21 13:59:39 eric
118 * Remplacement du membre psi par psi0.
119 * Modif calcul de d_psi a la fin.
120 *
121 * Revision 2.0 2000/02/17 18:50:44 eric
122 * *** empty log message ***
123 *
124 *
125 * $Header: /cvsroot/Lorene/C++/Source/Etoile/et_bin_vel_pot.C,v 1.15 2014/10/13 08:52:56 j_novak Exp $
126 *
127 */
128
129// Headers Lorene
130#include "scalar.h"
131#include "metrique.h"
132#include "etoile.h"
133#include "eos.h"
134#include "param.h"
135#include "et_bin_nsbh.h"
136#include "utilitaires.h"
137
138// Local prototype
139namespace Lorene {
140Cmp raccord_c1(const Cmp& uu, int l1) ;
141
142double Etoile_bin::velocity_potential(int mermax, double precis, double relax) {
143
144 // Which star is that ?
145 const Et_bin_nsbh* pnsbh = dynamic_cast<const Et_bin_nsbh*>(this) ;
146
147 if (eos.identify() == 5 || eos.identify() == 4 ||
148 eos.identify() == 3) {
149
150 // Routine used for binary strange stars.
151
152 int nzm1 = mp.get_mg()->get_nzone() - 1 ;
153
154 //----------------------------------
155 // Specific relativistic enthalpy ---> hhh
156 //----------------------------------
157
158 Tenseur hhh = exp(unsurc2 * ent) ; // = 1 at the Newtonian limit
159 hhh.set_std_base() ;
160
161 //----------------------------------------------
162 // Computation of W^i = - A^2 h Gamma_n B^i/N
163 // See Eq (62) from Gourgoulhon et al. (2001)
164 //----------------------------------------------
165
166 Tenseur www = - a_car * hhh * gam_euler * bsn ;
167
168 www.change_triad( mp.get_bvect_cart() ) ; // components on the mapping
169 // Cartesian basis
170
171 //-------------------------------------------------
172 // Constant value of W^i at the center of the star
173 //-------------------------------------------------
174
175 Tenseur v_orb(mp, 1, CON, mp.get_bvect_cart()) ;
176
177 v_orb.set_etat_qcq() ;
178 for (int i=0; i<3; i++) {
179 v_orb.set(i) = www(i)(0, 0, 0, 0) ;
180 }
181
182 v_orb.annule(nzm1, nzm1) ; // set to zero in the ZEC
183
184
185 v_orb.set_triad( *(www.get_triad()) ) ;
186 v_orb.set_std_base() ;
187
188
189 //-------------------------------------------------
190 // Source and coefficients a,b for poisson_compact (idenpendent from psi0)
191 //-------------------------------------------------
192
193 Cmp dndh_log = eos.der_nbar_ent(ent(), nzet) ;
194
195 // In order to avoid any division by zero in the computation of zeta_h
196 // the value of dndh_log is set to 1 in the external domains:
197 for (int l=nzet; l <= nzm1; l++) {
198 dndh_log.set(l) = 1 ;
199 }
200
201 double erreur ;
202
203 Tenseur zeta_h( ent() / dndh_log ) ;
204 zeta_h.set_std_base() ;
205
206 Scalar zeta_h_scalar (zeta_h()) ;
207 zeta_h_scalar.set_outer_boundary(0, (ent() / dndh_log)(0,0,0,0)) ;
208 for (int l=1; l<=nzm1; l++)
209 zeta_h_scalar.set_domain(l) = 1 ;
210
211 Cmp zeta_h_cmp (zeta_h_scalar) ;
212 zeta_h.set() = zeta_h_cmp ;
213 zeta_h.set_std_base() ;
214
215
216
217 Tenseur beta(mp) ;
218
219 if (pnsbh!=0x0) {
220 beta = log( sqrt(a_car) * nnn ) ;
221 beta.set_std_base() ;
222 }
223 else {
224 beta = beta_auto + beta_comp ;
225 }
226
227 Tenseur tmp_zeta = 1 - unsurc2 * zeta_h ;
228 tmp_zeta.set_std_base() ;
229
230 Tenseur bb = tmp_zeta * ent.gradient_spher()
231 + unsurc2 * zeta_h * beta.gradient_spher() ;
232
233 Tenseur entmb = ent - beta ;
234
235 Tenseur grad_ent (ent.gradient()) ;
236 grad_ent.change_triad(mp.get_bvect_spher()) ;
237
238 // Source for the poisson equation
239 // See Eq (63) from Gourgoulhon et al. (2001)
240 Tenseur source = flat_scalar_prod( www - v_orb, ent.gradient() )
241 + unsurc2 * zeta_h * (
242 flat_scalar_prod( v_orb, entmb.gradient() )
244 / gam_euler ) ;
245
246 for (int l=1; l<=nzm1; l++)
247 source.set().annule(l) ;
248
249 source = (source - flat_scalar_prod(bb, psi0.gradient_spher()))
250 / zeta_h ;
251 source.annule(nzet, nzm1) ;
252
253 Param par ;
254 int niter ;
255 par.add_int(mermax) ;
256 par.add_double(precis, 0) ;
257 par.add_double(relax, 1) ;
258 par.add_int_mod(niter) ;
259
260 par.add_cmp_mod(ssjm1_psi, 0) ;
261
262 if (psi0.get_etat() == ETATZERO) {
263 psi0.set_etat_qcq() ;
264 psi0.set() = 0 ;
265 }
266
267 int nr = mp.get_mg()->get_nr(0);
268 int nt = mp.get_mg()->get_nt(0);
269 int np = mp.get_mg()->get_np(0);
270
271 cout << "nr = " << nr << " nt = " << nt << " np = " << np << endl ;
272
273 cout << "psi0" << endl << norme(psi0()/(nr*nt*np)) << endl ;
274 cout << "d(psi)/dr" << endl << norme(psi0.set().dsdr()/(nr*nt*np)) << endl ;
275
276 Valeur lim(mp.get_mg()->get_angu()) ;
277 lim.annule_hard() ;
278
279 Tenseur normal (mp, 1, CON, mp.get_bvect_cart()) ;
280 Tenseur normal2 (mp, 1, COV, mp.get_bvect_cart()) ;
281 normal.set_etat_qcq() ;
282 normal2.set_etat_qcq() ;
283
284 const Coord& rr0 = mp.r ;
285 Tenseur rr(mp) ;
286 rr.set_etat_qcq() ;
287 rr.set() = rr0 ;
288 rr.set_std_base() ;
289
290 Tenseur_sym plat(mp, 2, COV, mp.get_bvect_cart() ) ;
291 plat.set_etat_qcq() ;
292 for (int i=0; i<3; i++) {
293 for (int j=0; j<i; j++) {
294 plat.set(i,j) = 0 ;
295 }
296 plat.set(i,i) = 1 ;
297 }
298 plat.set_std_base() ;
299
300 Metrique flat(plat, true) ;
301 Tenseur dcov_r = rr.derive_cov(flat) ;
302
303
304 for (int i=0; i<3; i++) {
305 normal.set(i) = dcov_r(i) ;
306 normal2.set(i) = dcov_r(i) ;
307 }
308
310 normal2.change_triad(mp.get_bvect_spher()) ;
311
312
313
314 Tenseur bsn0 (bsn) ;
316 Tenseur aa (mp, 1, CON, mp.get_bvect_cart()) ;
317 aa = - v_orb - a_car * gam_euler * hhh * bsn0 ;
319
320
321 Tenseur dcov_psi = psi0.derive_cov(flat) ;
322 dcov_psi.change_triad(mp.get_bvect_spher()) ;
323
324 Cmp limite (mp) ;
325 limite = ( - dcov_psi(1) * normal(1) - dcov_psi(2) * normal(2)
326 + contract(aa, 0, normal2, 0)())
327 /normal(0) ;
328
329 for (int j=0; j<nt; j++)
330 for (int k=0; k<np; k++)
331 lim.set(0, k, j, 0) = limite(0, k, j, nr-1) ;
332
333// cout << "lim" << endl << lim << endl ;
334
335 lim.std_base_scal() ;
336 Cmp resu (psi0()) ;
337 source().poisson_neumann_interne(lim, par, resu) ;
338 psi0 = resu ;
339
340/*
341 resu.va.ylm() ;
342 Scalar psi00(resu) ;
343 psi00.spectral_display("psi00") ;
344
345 cout << "value of d(psi)/dr at the surface after poisson" << endl ;
346 for (int j=0; j<nt; j++)
347 for (int k=0; k<np; k++)
348 cout << "j = " << j << " ; k = " << k << " : " <<
349 psi0.set().dsdr()(0, k, j, nr-1) << endl ;
350*/
351 for (int l=1; l<=nzm1; l++)
352 psi0.set().annule(l) ;
353
354
355 //---------------------------------------------------
356 // Check of the solution
357 //---------------------------------------------------
358
359 Cmp laplacien_psi0 = psi0().laplacien() ;
360
361 erreur = diffrel(laplacien_psi0, source())(0) ;
362
363 cout << "Check of the resolution of the continuity equation for strange stars: "
364 << endl ;
365 cout << "norme(source) : " << norme(source())(0) << endl
366 << "Error in the solution : " << erreur << endl ;
367
368 //--------------------------------
369 // Computation of grad(psi)
370 //--------------------------------
371
372 // The computation is done component by component because psi0.gradient()
373 // is a covariant vector, whereas v_orb is a contravariant one.
374
376
377 for (int i=0; i<3; i++) {
378 d_psi.set(i) = (psi0.gradient())(i) + v_orb(i) ;
379 }
380
381 d_psi.set_triad( *(v_orb.get_triad()) ) ;
382
383 // C^1 continuation of d_psi outside the star
384 // (to ensure a smooth enthalpy field accross the stellar surface)
385 // ----------------------------------------------------------------
386
387 d_psi.annule(nzet, nzm1) ;
388 for (int i=0; i<3; i++) {
389 d_psi.set(i) = raccord_c1(d_psi(i), nzet) ;
390 }
391
392 assert( d_psi.get_triad() == &(mp.get_bvect_cart()) ) ;
393
395
396 return erreur ;
397
398
399 } // End of strange stars case
400
401//=============================================================================
402
403 else {
404
405 int nzm1 = mp.get_mg()->get_nzone() - 1 ;
406
407 //----------------------------------
408 // Specific relativistic enthalpy ---> hhh
409 //----------------------------------
410
411 Tenseur hhh = exp(unsurc2 * ent) ; // = 1 at the Newtonian limit
412 hhh.set_std_base() ;
413
414 //----------------------------------------------
415 // Computation of W^i = - A^2 h Gamma_n B^i/N
416 // See Eq (62) from Gourgoulhon et al. (2001)
417 //----------------------------------------------
418
419 Tenseur www = - a_car * hhh * gam_euler * bsn ;
420
421 www.change_triad( mp.get_bvect_cart() ) ; // components on the mapping
422 // Cartesian basis
423
424 //-------------------------------------------------
425 // Constant value of W^i at the center of the star
426 //-------------------------------------------------
427
428 Tenseur v_orb(mp, 1, CON, mp.get_bvect_cart()) ;
429
430 v_orb.set_etat_qcq() ;
431 for (int i=0; i<3; i++) {
432 v_orb.set(i) = www(i)(0, 0, 0, 0) ;
433 }
434
435 v_orb.set_triad( *(www.get_triad()) ) ;
436 v_orb.set_std_base() ;
437
438
439 //-------------------------------------------------
440 // Source and coefficients a,b for poisson_compact (idenpendent from psi0)
441 //-------------------------------------------------
442
443 Cmp dndh_log(mp) ;
444 dndh_log = 0 ;
445
446 for (int l=0; l<nzet; l++) {
447
448 Param par ; // Paramater for multi-domain equation of state
449 par.add_int(l) ;
450
451 dndh_log = dndh_log + eos.der_nbar_ent(ent(), 1, l, &par) ;
452
453 }
454
455 // Cmp dndh_log = eos.der_nbar_ent(ent(), nzet) ;
456
457 // In order to avoid any division by zero in the computation of zeta_h
458 // the value of dndh_log is set to 1 in the external domains:
459 for (int l=nzet; l <= nzm1; l++) {
460 dndh_log.set(l) = 1 ;
461 }
462
463 Tenseur zeta_h( ent() / dndh_log ) ;
464 zeta_h.set_std_base() ;
465
466 Tenseur beta(mp) ;
467
468 if (pnsbh!=0x0) {
469 beta = log( sqrt(a_car) * nnn ) ;
470 beta.set_std_base() ;
471 }
472 else {
473 beta = beta_auto + beta_comp ;
474 }
475
476 Tenseur tmp_zeta = 1 - unsurc2 * zeta_h ;
477 tmp_zeta.set_std_base() ;
478
479 Tenseur bb = tmp_zeta * ent.gradient_spher()
480 + unsurc2 * zeta_h * beta.gradient_spher() ;
481
482 Tenseur entmb = ent - beta ;
483
484 // See Eq (63) from Gourgoulhon et al. (2001)
485 Tenseur source = flat_scalar_prod( www - v_orb, ent.gradient() )
486 + unsurc2 * zeta_h * (
487 flat_scalar_prod( v_orb, entmb.gradient() )
489 / gam_euler ) ;
490
491
492 source.annule(nzet, nzm1) ;
493
494 //---------------------------------------------------
495 // Resolution by means of Map_radial::poisson_compact
496 //---------------------------------------------------
497
498 Param par ;
499 int niter ;
500 par.add_int(mermax) ;
501 par.add_double(precis, 0) ;
502 par.add_double(relax, 1) ;
503 par.add_int_mod(niter) ;
504
505
506 if (psi0.get_etat() == ETATZERO) {
507 psi0.set_etat_qcq() ;
508 psi0.set() = 0 ;
509 }
510
511 source.set().va.ylm() ;
512
513 mp.poisson_compact(nzet, source(), zeta_h(), bb, par, psi0.set() ) ;
514
515 //---------------------------------------------------
516 // Check of the solution
517 //---------------------------------------------------
518
519 Tenseur bb_dpsi0 = flat_scalar_prod( bb, psi0.gradient_spher() ) ;
520
521 Cmp oper = zeta_h() * psi0().laplacien() + bb_dpsi0() ;
522
523 source.set().va.ylm_i() ;
524
525 cout << "Check of the resolution of the continuity equation : " << endl ;
526 Tbl terr = diffrel(oper, source()) ;
527 double erreur = 0 ;
528 for (int l=0; l<nzet; l++) {
529 double err = terr(l) ;
530 cout << " domain " << l << " : norme(source) : " << norme(source())(l)
531 << " relative error : " << err << endl ;
532 if (err > erreur) erreur = err ;
533 }
534 // arrete() ;
535
536 //--------------------------------
537 // Computation of grad(psi)
538 //--------------------------------
539
540 // The computation is done component by component because psi0.gradient()
541 // is a covariant vector, whereas v_orb is a contravariant one.
542
544
545 for (int i=0; i<3; i++) {
546 d_psi.set(i) = (psi0.gradient())(i) + v_orb(i) ;
547 }
548
549 d_psi.set_triad( *(v_orb.get_triad()) ) ;
550
551 // C^1 continuation of d_psi outside the star
552 // (to ensure a smooth enthalpy field accross the stellar surface)
553 // ----------------------------------------------------------------
554
555 d_psi.annule(nzet, nzm1) ;
556 for (int i=0; i<3; i++) {
557 d_psi.set(i) = raccord_c1(d_psi(i), nzet) ;
558 }
559
560 assert( d_psi.get_triad() == &(mp.get_bvect_cart()) ) ;
561
563
564 return erreur ;
565
566
567 }
568}
569}
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition cmp.h:446
Valeur va
The numerical value of the Cmp
Definition cmp.h:464
void annule(int l)
Sets the Cmp to zero in a given domain.
Definition cmp.C:348
Tbl & set(int l)
Read/write of the value in a given domain.
Definition cmp.h:724
const Cmp & dsdr() const
Returns of *this .
Definition cmp_deriv.C:84
Active physical coordinates and mapping derivatives.
Definition coord.h:90
virtual int identify() const =0
Returns a number to identify the sub-classe of Eos the object belongs to.
Cmp der_nbar_ent(const Cmp &ent, int nzet, int l_min=0, const Param *par=0x0) const
Computes the logarithmic derivative from the log-enthalpy and extra parameters.
Definition eos.C:407
Class for a star in a NS-BH binary system.
Definition et_bin_nsbh.h:79
Tenseur psi0
Scalar potential of the non-translational part of fluid 4-velocity (in the irrotational case)
Definition etoile.h:833
Tenseur bsn
3-vector shift, divided by N , of the rotating coordinates, .
Definition etoile.h:950
Tenseur d_psi
Gradient of (in the irrotational case) (Cartesian components with respect to ref_triad )
Definition etoile.h:838
double velocity_potential(int mermax, double precis, double relax)
Computes the non-translational part of the velocity scalar potential by solving the continuity equat...
const Base_vect & ref_triad
Reference triad ("absolute frame"), with respect to which the components of all the member Tenseur 's...
Definition etoile.h:828
Tenseur beta_comp
Part of the logarithm of AN generated principaly by the companion star.
Definition etoile.h:874
Cmp ssjm1_psi
Effective source at the previous step for the resolution of the Poisson equation for the scalar by m...
Definition etoile.h:989
int nzet
Number of domains of *mp occupied by the star.
Definition etoile.h:432
Tenseur nnn
Total lapse function.
Definition etoile.h:509
const Eos & eos
Equation of state of the stellar matter.
Definition etoile.h:451
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 ent
Log-enthalpy (relativistic case) or specific enthalpy (Newtonian case)
Definition etoile.h:457
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
double unsurc2
: unsurc2=1 for a relativistic star, 0 for a Newtonian one.
Definition etoile.h:442
const Base_vect_cart & get_bvect_cart() const
Returns the Cartesian basis associated with the coordinates (x,y,z) of the mapping,...
Definition map.h:791
virtual void poisson_compact(const Cmp &source, const Cmp &aa, const Tenseur &bb, const Param &par, Cmp &psi) const =0
Resolution of the elliptic equation in the case where the stellar interior is covered by a single do...
Coord r
r coordinate centered on the grid
Definition map.h:718
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
const Mg3d * get_angu() const
Returns the pointer on the associated angular grid.
Definition mg3d.C:473
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
Parameter storage.
Definition param.h:125
void add_double(const double &x, int position=0)
Adds the the address of a new double to the list.
Definition param.C:315
void add_cmp_mod(Cmp &ti, int position=0)
Adds the address of a new modifiable Cmp to the list.
Definition param.C:1004
void add_int_mod(int &n, int position=0)
Adds the address of a new modifiable int to the list.
Definition param.C:385
void add_int(const int &n, int position=0)
Adds the address of a new int to the list.
Definition param.C:246
Tensor field of valence 0 (or component of a tensorial field).
Definition scalar.h:387
Tbl & set_domain(int l)
Read/write of the value in a given domain.
Definition scalar.h:615
void set_outer_boundary(int l, double x)
Sets the value of the Scalar at the outer boundary of a given domain.
Basic array class.
Definition tbl.h:161
Class intended to describe tensors with a symmetry on the two last indices *** DEPRECATED : use class...
Definition tenseur.h:1253
Tensor handling *** DEPRECATED : use class Tensor instead ***.
Definition tenseur.h:301
const Base_vect * get_triad() const
Returns the vectorial basis (triad) on which the components are defined.
Definition tenseur.h:704
Cmp & set()
Read/write for a scalar (see also operator=(const Cmp&) ).
Definition tenseur.C:824
const Tenseur & derive_cov(const Metrique &met) const
Returns the covariant derivative of *this , with respect to met .
Definition tenseur.C:1554
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition tenseur.C:636
void annule(int l)
Sets the Tenseur to zero in a given domain.
Definition tenseur.C:900
const Tenseur & gradient() const
Returns the gradient of *this (Cartesian coordinates)
Definition tenseur.C:1542
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
void change_triad(const Base_vect &new_triad)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
Definition tenseur.C:668
void set_triad(const Base_vect &new_triad)
Assigns a new vectorial basis (triad) of decomposition.
Definition tenseur.C:674
int get_etat() const
Returns the logical state.
Definition tenseur.h:707
Values and coefficients of a (real-value) function.
Definition valeur.h:287
void ylm()
Computes the coefficients of *this.
Definition valeur_ylm.C:138
Tbl & set(int l)
Read/write of the value in a given domain (configuration space).
Definition valeur.h:363
void ylm_i()
Inverse of ylm()
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 exp(const Cmp &)
Exponential.
Definition cmp_math.C:270
Tbl diffrel(const Cmp &a, const Cmp &b)
Relative difference between two Cmp (norme version).
Definition cmp_math.C:504
Tbl norme(const Cmp &)
Sums of the absolute values of all the values of the Cmp in each domain.
Definition cmp_math.C:481
Cmp log(const Cmp &)
Neperian logarithm.
Definition cmp_math.C:296
Tenseur contract(const Tenseur &, int id1, int id2)
Self contraction of two indices of a Tenseur .
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