LORENE
bhole_kij.C
1 /*
2  * Copyright (c) 2001 Philippe Grandclement
3  *
4  * This file is part of LORENE.
5  *
6  * LORENE is free software; you can redistribute it and/or modify
7  * it under the terms of the GNU General Public License as published by
8  * the Free Software Foundation; either version 2 of the License, or
9  * (at your option) any later version.
10  *
11  * LORENE is distributed in the hope that it will be useful,
12  * but WITHOUT ANY WARRANTY; without even the implied warranty of
13  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  * GNU General Public License for more details.
15  *
16  * You should have received a copy of the GNU General Public License
17  * along with LORENE; if not, write to the Free Software
18  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
19  *
20  */
21 
22 
23 char bhole_kij_C[] = "$Header: /cvsroot/Lorene/C++/Source/Bhole_binaire/bhole_kij.C,v 1.6 2014/10/13 08:52:40 j_novak Exp $" ;
24 
25 /*
26  * $Id: bhole_kij.C,v 1.6 2014/10/13 08:52:40 j_novak Exp $
27  * $Log: bhole_kij.C,v $
28  * Revision 1.6 2014/10/13 08:52:40 j_novak
29  * Lorene classes and functions now belong to the namespace Lorene.
30  *
31  * Revision 1.5 2014/10/06 15:12:58 j_novak
32  * Modified #include directives to use c++ syntax.
33  *
34  * Revision 1.4 2005/08/29 15:10:14 p_grandclement
35  * Addition of things needed :
36  * 1) For BBH with different masses
37  * 2) Provisory files for the mixted binaries (Bh and NS) : THIS IS NOT
38  * WORKING YET !!!
39  *
40  * Revision 1.3 2004/05/27 07:10:22 p_grandclement
41  * Correction of some shadowed variables
42  *
43  * Revision 1.2 2002/10/16 14:36:33 j_novak
44  * Reorganization of #include instructions of standard C++, in order to
45  * use experimental version 3 of gcc.
46  *
47  * Revision 1.1.1.1 2001/11/20 15:19:28 e_gourgoulhon
48  * LORENE
49  *
50  * Revision 1.6 2001/05/07 09:12:02 phil
51  * *** empty log message ***
52  *
53  * Revision 1.5 2001/04/30 09:22:52 phil
54  * cas omega = 0
55  *
56  * Revision 1.4 2001/04/27 11:44:01 phil
57  * correction devant assurer la symetrie entre les deux trous
58  *
59  * Revision 1.3 2001/04/26 12:04:44 phil
60  * *** empty log message ***
61  *
62  * Revision 1.2 2001/04/25 15:54:30 phil
63  * corrections diverses rien de bien mechant
64  *
65  * Revision 1.1 2001/04/25 15:10:23 phil
66  * Initial revision
67  *
68  *
69  * $Header: /cvsroot/Lorene/C++/Source/Bhole_binaire/bhole_kij.C,v 1.6 2014/10/13 08:52:40 j_novak Exp $
70  *
71  */
72 
73 
74 //standard
75 #include <cstdlib>
76 #include <cmath>
77 
78 // Lorene
79 #include "nbr_spx.h"
80 #include "tenseur.h"
81 #include "bhole.h"
82 #include "proto.h"
83 #include "utilitaires.h"
84 #include "graphique.h"
85 
86 namespace Lorene {
87 //calcul de kij total. (la regularisation ayant ete faite)
89 
92 
93  // On construit a_ij a partir du shift ...
94  // taij tot doit etre nul sur les deux horizons.
97 
98  // On trouve les trucs du compagnon
101 
102  Tenseur sans_dz_un (hole1.taij_auto) ;
103  sans_dz_un.dec2_dzpuis() ;
104  Tenseur sans_dz_deux (hole2.taij_auto) ;
105  sans_dz_deux.dec2_dzpuis() ;
106 
107  // ON DOIT VERIFIER SI LES DEUX Aij sont alignes ou non !
108  // Les bases des deux vecteurs doivent etre alignees ou non alignees :
109  double orientation_un = hole1.taij_auto.get_mp()->get_rot_phi() ;
110  assert ((orientation_un==0) || (orientation_un==M_PI)) ;
111  double orientation_deux = hole2.taij_auto.get_mp()->get_rot_phi() ;
112  assert ((orientation_deux==0) || (orientation_deux==M_PI)) ;
113  int same_orient = (orientation_un == orientation_deux) ? 1 : -1 ;
114 
115  // Les importations :
116  if (hole2.taij_auto.get_etat() == ETATZERO)
118  else {
119  hole1.taij_comp.set(0, 0).import_asymy(sans_dz_deux(0, 0)) ;
120  hole1.taij_comp.set(0, 1).import_symy(sans_dz_deux(0, 1)) ;
121  hole1.taij_comp.set(0, 2).import_asymy(same_orient*sans_dz_deux(0, 2)) ;
122  hole1.taij_comp.set(1, 1).import_asymy(sans_dz_deux(1, 1)) ;
123  hole1.taij_comp.set(1, 2).import_symy(same_orient*sans_dz_deux(1, 2)) ;
124  hole1.taij_comp.set(2, 2).import_asymy(sans_dz_deux(2, 2)) ;
125  }
126 
127  if (hole1.taij_auto.get_etat() == ETATZERO)
129  else {
130  hole2.taij_comp.set(0, 0).import_asymy(sans_dz_un(0, 0)) ;
131  hole2.taij_comp.set(0, 1).import_symy(sans_dz_un(0, 1)) ;
132  hole2.taij_comp.set(0, 2).import_asymy(same_orient*sans_dz_un(0, 2)) ;
133  hole2.taij_comp.set(1, 1).import_asymy(sans_dz_un(1, 1)) ;
134  hole2.taij_comp.set(1, 2).import_symy(same_orient*sans_dz_un(1, 2)) ;
135  hole2.taij_comp.set(2, 2).import_asymy(sans_dz_un(2, 2)) ;
136  }
137 
142 
143  // Et enfin les trucs totaux...
146 
147  if ((hole1.taij_tot.get_etat() == ETATZERO) &&
148  (hole2.taij_tot.get_etat() == ETATZERO)) {
149 
154  }
155  else {
156  int nz_un = hole1.mp.get_mg()->get_nzone() ;
157  int nz_deux = hole2.mp.get_mg()->get_nzone() ;
158 
159  Cmp ntot_un (hole1.n_tot()) ;
160  ntot_un = division_xpun (ntot_un, 0) ;
161  ntot_un.raccord(1) ;
162 
163  Cmp ntot_deux (hole2.n_tot()) ;
164  ntot_deux = division_xpun (ntot_deux, 0) ;
165  ntot_deux.raccord(1) ;
166 
167  // Boucle sur les composantes :
168  for (int lig = 0 ; lig<3 ; lig++)
169  for (int col = lig ; col<3 ; col++) {
170 
171  // Le sens d orientation
172  int ind = 1 ;
173  if (lig !=2)
174  ind *= -1 ;
175  if (col != 2)
176  ind *= -1 ;
177  if (same_orient == 1)
178  ind = 1 ;
179 
180  // Pres de H1 :
181  Cmp auxi_un (hole1.taij_tot(lig, col)/2.) ;
182  auxi_un.dec2_dzpuis() ;
183  auxi_un = division_xpun (auxi_un, 0) ;
184  auxi_un = auxi_un / ntot_un ;
185  auxi_un.raccord(1) ;
186 
187  // Pres de H2 :
188  Cmp auxi_deux (hole2.taij_tot(lig, col)/2.) ;
189  auxi_deux.dec2_dzpuis() ;
190  auxi_deux = division_xpun (auxi_deux, 0) ;
191  auxi_deux = auxi_deux / ntot_deux ;
192  auxi_deux.raccord(1) ;
193 
194  // copie :
195  Cmp copie_un (hole1.taij_tot(lig, col)) ;
196  copie_un.dec2_dzpuis() ;
197 
198  Cmp copie_deux (hole2.taij_tot(lig, col)) ;
199  copie_deux.dec2_dzpuis() ;
200 
201  // Double les rayons limites :
202  double lim_un = hole1.mp.get_alpha()[1] + hole1.mp.get_beta()[1] ;
203  double lim_deux = hole2.mp.get_alpha()[1] + hole2.mp.get_beta()[1] ;
204 
205  Mtbl xabs_un (hole1.mp.xa) ;
206  Mtbl yabs_un (hole1.mp.ya) ;
207  Mtbl zabs_un (hole1.mp.za) ;
208 
209  Mtbl xabs_deux (hole2.mp.xa) ;
210  Mtbl yabs_deux (hole2.mp.ya) ;
211  Mtbl zabs_deux (hole2.mp.za) ;
212 
213  double xabs, yabs, zabs, air, theta, phi ;
214 
215  // On boucle sur les autres zones :
216  for (int l=2 ; l<nz_un ; l++) {
217 
218  int nr = hole1.mp.get_mg()->get_nr (l) ;
219 
220  if (l==nz_un-1)
221  nr -- ;
222 
223  int np = hole1.mp.get_mg()->get_np (l) ;
224  int nt = hole1.mp.get_mg()->get_nt (l) ;
225 
226  for (int k=0 ; k<np ; k++)
227  for (int j=0 ; j<nt ; j++)
228  for (int i=0 ; i<nr ; i++) {
229 
230  xabs = xabs_un (l, k, j, i) ;
231  yabs = yabs_un (l, k, j, i) ;
232  zabs = zabs_un (l, k, j, i) ;
233 
234  // les coordonnees du point en 2 :
236  (xabs, yabs, zabs, air, theta, phi) ;
237 
238  if (air >= lim_deux)
239  // On est loin du trou 2 : pas de pb :
240  auxi_un.set(l, k, j, i) =
241  copie_un(l, k, j, i) / ntot_un (l, k, j, i)/2. ;
242  else
243  // On est pres du trou deux :
244  auxi_un.set(l, k, j, i) =
245  ind * auxi_deux.val_point (air, theta, phi) ;
246  }
247 
248  // Cas infini :
249  if (l==nz_un-1)
250  for (int k=0 ; k<np ; k++)
251  for (int j=0 ; j<nt ; j++)
252  auxi_un.set(nz_un-1, k, j, nr-1) = 0 ;
253  }
254 
255  // Le second trou :
256  for (int l=2 ; l<nz_deux ; l++) {
257 
258  int nr = hole2.mp.get_mg()->get_nr (l) ;
259 
260  if (l==nz_deux-1)
261  nr -- ;
262 
263  int np = hole2.mp.get_mg()->get_np (l) ;
264  int nt = hole2.mp.get_mg()->get_nt (l) ;
265 
266  for (int k=0 ; k<np ; k++)
267  for (int j=0 ; j<nt ; j++)
268  for (int i=0 ; i<nr ; i++) {
269 
270  xabs = xabs_deux (l, k, j, i) ;
271  yabs = yabs_deux (l, k, j, i) ;
272  zabs = zabs_deux (l, k, j, i) ;
273 
274  // les coordonnees du point en 1 :
276  (xabs, yabs, zabs, air, theta, phi) ;
277 
278  if (air >= lim_un)
279  // On est loin du trou 1 : pas de pb :
280  auxi_deux.set(l, k, j, i) =
281  copie_deux(l, k, j, i) / ntot_deux (l, k, j, i) /2.;
282  else
283  // On est pres du trou deux :
284  auxi_deux.set(l, k, j, i) =
285  ind * auxi_un.val_point (air, theta, phi) ;
286  }
287  // Cas infini :
288  if (l==nz_deux-1)
289  for (int k=0 ; k<np ; k++)
290  for (int j=0 ; j<nt ; j++)
291  auxi_deux.set(nz_deux-1, k, j, nr-1) = 0 ;
292  }
293 
294  auxi_un.inc2_dzpuis() ;
295  auxi_deux.inc2_dzpuis() ;
296 
297  hole1.tkij_tot.set(lig, col) = auxi_un ;
298  hole2.tkij_tot.set(lig, col) = auxi_deux ;
299  }
300 
303 
304  for (int lig=0 ; lig<3 ; lig++)
305  for (int col=lig ; col<3 ; col++) {
306  hole1.tkij_auto.set(lig, col) = hole1.tkij_tot(lig, col)*
307  hole1.decouple ;
308  hole2.tkij_auto.set(lig, col) = hole2.tkij_tot(lig, col)*
309  hole2.decouple ;
310  }
311  }
312 }
313 
315 
316  int nz_un = hole1.mp.get_mg()->get_nzone() ;
317  int nz_deux = hole2.mp.get_mg()->get_nzone() ;
318 
319  // On determine R_limite (pour le moment en tout cas...) :
320  double distance = hole1.mp.get_ori_x() - hole2.mp.get_ori_x() ;
321  double lim_un = distance/2. ;
322  double lim_deux = distance/2. ;
323  double int_un = distance/6. ;
324  double int_deux = distance/6. ;
325 
326  // Les fonctions de base
327  Cmp fonction_f_un (hole1.mp) ;
328  fonction_f_un = 0.5*pow(
329  cos((hole1.mp.r-int_un)*M_PI/2./(lim_un-int_un)), 2.)+0.5 ;
330  fonction_f_un.std_base_scal();
331 
332  Cmp fonction_g_un (hole1.mp) ;
333  fonction_g_un = 0.5*pow
334  (sin((hole1.mp.r-int_un)*M_PI/2./(lim_un-int_un)), 2.) ;
335  fonction_g_un.std_base_scal();
336 
337  Cmp fonction_f_deux (hole2.mp) ;
338  fonction_f_deux = 0.5*pow(
339  cos((hole2.mp.r-int_deux)*M_PI/2./(lim_deux-int_deux)), 2.)+0.5 ;
340  fonction_f_deux.std_base_scal();
341 
342  Cmp fonction_g_deux (hole2.mp) ;
343  fonction_g_deux = 0.5*pow(
344  sin((hole2.mp.r-int_deux)*M_PI/2./(lim_un-int_deux)), 2.) ;
345  fonction_g_deux.std_base_scal();
346 
347  // Les fonctions totales :
348  Cmp decouple_un (hole1.mp) ;
349  decouple_un.allocate_all() ;
350  Cmp decouple_deux (hole2.mp) ;
351  decouple_deux.allocate_all() ;
352 
353  Mtbl xabs_un (hole1.mp.xa) ;
354  Mtbl yabs_un (hole1.mp.ya) ;
355  Mtbl zabs_un (hole1.mp.za) ;
356 
357  Mtbl xabs_deux (hole2.mp.xa) ;
358  Mtbl yabs_deux (hole2.mp.ya) ;
359  Mtbl zabs_deux (hole2.mp.za) ;
360 
361  double xabs, yabs, zabs, air_un, air_deux, theta, phi ;
362 
363  // On boucle sur les autres zones :
364  for (int l=0 ; l<nz_un ; l++) {
365  int nr = hole1.mp.get_mg()->get_nr (l) ;
366 
367  if (l==nz_un-1)
368  nr -- ;
369 
370  int np = hole1.mp.get_mg()->get_np (l) ;
371  int nt = hole1.mp.get_mg()->get_nt (l) ;
372 
373  for (int k=0 ; k<np ; k++)
374  for (int j=0 ; j<nt ; j++)
375  for (int i=0 ; i<nr ; i++) {
376 
377  xabs = xabs_un (l, k, j, i) ;
378  yabs = yabs_un (l, k, j, i) ;
379  zabs = zabs_un (l, k, j, i) ;
380 
381  // les coordonnees du point :
383  (xabs, yabs, zabs, air_un, theta, phi) ;
385  (xabs, yabs, zabs, air_deux, theta, phi) ;
386 
387  if (air_un <= lim_un)
388  if (air_un < int_un)
389  decouple_un.set(l, k, j, i) = 1 ;
390  else
391  // pres du trou un :
392  decouple_un.set(l, k, j, i) =
393  fonction_f_un (l, k, j, i) ;
394  else
395  if (air_deux <= lim_deux)
396  if (air_deux < int_deux)
397  decouple_un.set(l, k, j, i) = 0 ;
398  else
399  // On est pres du trou deux :
400  decouple_un.set(l, k, j, i) =
401  fonction_g_deux.val_point (air_deux, theta, phi) ;
402 
403  else
404  // On est loin des deux trous :
405  decouple_un.set(l, k, j, i) = 0.5 ;
406  }
407 
408  // Cas infini :
409  if (l==nz_un-1)
410  for (int k=0 ; k<np ; k++)
411  for (int j=0 ; j<nt ; j++)
412  decouple_un.set(nz_un-1, k, j, nr) = 0.5 ;
413  }
414 
415  for (int l=0 ; l<nz_deux ; l++) {
416 
417  int nr = hole2.mp.get_mg()->get_nr (l) ;
418 
419  if (l==nz_deux-1)
420  nr -- ;
421 
422  int np = hole2.mp.get_mg()->get_np (l) ;
423  int nt = hole2.mp.get_mg()->get_nt (l) ;
424 
425  for (int k=0 ; k<np ; k++)
426  for (int j=0 ; j<nt ; j++)
427  for (int i=0 ; i<nr ; i++) {
428 
429  xabs = xabs_deux (l, k, j, i) ;
430  yabs = yabs_deux (l, k, j, i) ;
431  zabs = zabs_deux (l, k, j, i) ;
432 
433  // les coordonnees du point :
435  (xabs, yabs, zabs, air_un, theta, phi) ;
437  (xabs, yabs, zabs, air_deux, theta, phi) ;
438 
439  if (air_deux <= lim_deux)
440  if (air_deux < int_deux)
441  decouple_deux.set(l, k, j, i) = 1 ;
442  else
443  // pres du trou deux :
444  decouple_deux.set(l, k, j, i) =
445  fonction_f_deux (l, k, j, i) ;
446  else
447  if (air_un <= lim_un)
448  if (air_un < int_un)
449  decouple_deux.set(l, k, j, i) = 0 ;
450  else
451  // On est pres du trou un :
452  decouple_deux.set(l, k, j, i) =
453  fonction_g_un.val_point (air_un, theta, phi) ;
454 
455  else
456  // On est loin des deux trous :
457  decouple_deux.set(l, k, j, i) = 0.5 ;
458  }
459 
460  // Cas infini :
461  if (l==nz_deux-1)
462  for (int k=0 ; k<np ; k++)
463  for (int j=0 ; j<nt ; j++)
464  decouple_deux.set(nz_deux-1, k, j, nr) = 0.5 ;
465  }
466 
467  hole1.decouple = decouple_un ;
468  hole2.decouple = decouple_deux ;
469 }
470 }
Lorene::Cmp::import_asymy
void import_asymy(const Cmp &ci)
Assignment to another Cmp defined on a different mapping.
Definition: cmp_import_asymy.C:67
Lorene::Bhole::taij_auto
Tenseur_sym taij_auto
Part of generated by the hole.
Definition: bhole.h:299
Lorene::Bhole::tkij_tot
Tenseur_sym tkij_tot
Total .
Definition: bhole.h:308
Lorene::Bhole::mp
Map_af & mp
Affine mapping.
Definition: bhole.h:273
Lorene::Tenseur::set_etat_zero
void set_etat_zero()
Sets the logical state to ETATZERO (zero state).
Definition: tenseur.C:645
Lorene::Mg3d::get_np
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
Definition: grilles.h:462
Lorene::Map::get_ori_x
double get_ori_x() const
Returns the x coordinate of the origin.
Definition: map.h:768
Lorene::Tenseur
Tensor handling *** DEPRECATED : use class Tensor instead ***.
Definition: tenseur.h:301
Lorene::Bhole_binaire::fait_decouple
void fait_decouple()
Calculates {tt decouple} which is used to obtain tkij_auto by the formula : tkij_auto = decouple * tk...
Definition: bhole_kij.C:314
Lorene::Bhole::fait_taij_auto
void fait_taij_auto()
Calculates the part of generated by shift_auto .
Definition: bhole.C:379
Lorene::Tenseur::set
Cmp & set()
Read/write for a scalar (see also operator=(const Cmp&) ).
Definition: tenseur.C:824
Lorene::Map::get_rot_phi
double get_rot_phi() const
Returns the angle between the x –axis and X –axis.
Definition: map.h:775
Lorene::Tenseur::inc2_dzpuis
void inc2_dzpuis()
dzpuis += 2 ;
Definition: tenseur.C:1143
Lorene::Mg3d::get_nt
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
Definition: grilles.h:457
Lorene
Lorene prototypes.
Definition: app_hor.h:64
Lorene::Bhole_binaire::hole2
Bhole hole2
Black hole two.
Definition: bhole.h:763
Lorene::Mtbl
Multi-domain array.
Definition: mtbl.h:118
Lorene::Map::get_mg
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:765
Lorene::Bhole_binaire::fait_tkij
void fait_tkij()
Calculation af the extrinsic curvature tensor.
Definition: bhole_kij.C:88
Lorene::Cmp::val_point
double val_point(double r, double theta, double phi) const
Computes the value of the field represented by *this at an arbitrary point , by means of the spectral...
Definition: cmp.C:732
Lorene::Map::r
Coord r
r coordinate centered on the grid
Definition: map.h:718
Lorene::Bhole::n_tot
Tenseur n_tot
Total N .
Definition: bhole.h:288
Lorene::pow
Cmp pow(const Cmp &, int)
Power .
Definition: cmp_math.C:348
Lorene::Bhole_binaire::hole1
Bhole hole1
Black hole one.
Definition: bhole.h:762
Lorene::Tenseur::get_etat
int get_etat() const
Returns the logical state.
Definition: tenseur.h:707
Lorene::Cmp::inc2_dzpuis
void inc2_dzpuis()
Increases by 2 the value of dzpuis and changes accordingly the values of the Cmp in the external comp...
Definition: cmp_r_manip.C:192
Lorene::Cmp
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition: cmp.h:446
Lorene::Mg3d::get_nzone
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:448
Lorene::Mg3d::get_nr
int get_nr(int l) const
Returns the number of points in the radial direction ( ) in domain no. l.
Definition: grilles.h:452
Lorene::Cmp::set
Tbl & set(int l)
Read/write of the value in a given domain.
Definition: cmp.h:724
Lorene::Bhole::tkij_auto
Tenseur_sym tkij_auto
Auto .
Definition: bhole.h:307
Lorene::Tenseur::get_mp
const Map * get_mp() const
Returns pointer on the mapping.
Definition: tenseur.h:699
Lorene::Bhole::decouple
Cmp decouple
Function used to construct the part of generated by the hole from the total .
Definition: bhole.h:318
Lorene::Cmp::dec2_dzpuis
void dec2_dzpuis()
Decreases by 2 the value of dzpuis and changes accordingly the values of the Cmp in the external comp...
Definition: cmp_r_manip.C:180
Lorene::Map::xa
Coord xa
Absolute x coordinate.
Definition: map.h:730
Lorene::Tenseur::set_etat_qcq
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: tenseur.C:636
Lorene::Cmp::import_symy
void import_symy(const Cmp &ci)
Assignment to another Cmp defined on a different mapping.
Definition: cmp_import_symy.C:67
Lorene::cos
Cmp cos(const Cmp &)
Cosine.
Definition: cmp_math.C:94
Lorene::Map_af::get_beta
const double * get_beta() const
Returns the pointer on the array beta.
Definition: map_af.C:481
Lorene::Cmp::allocate_all
void allocate_all()
Sets the logical state to ETATQCQ (ordinary state) and performs the memory allocation of all the elem...
Definition: cmp.C:323
Lorene::Map::convert_absolute
void convert_absolute(double xx, double yy, double zz, double &rr, double &theta, double &pphi) const
Determines the coordinates corresponding to given absolute Cartesian coordinates (X,...
Definition: map.C:302
Lorene::Bhole::taij_comp
Tenseur_sym taij_comp
Part of generated by the companion hole.
Definition: bhole.h:300
Lorene::Cmp::raccord
void raccord(int n)
Performs the matching of the nucleus with respect to the first shell.
Definition: cmp_raccord.C:170
Lorene::Map::za
Coord za
Absolute z coordinate.
Definition: map.h:732
Lorene::Map_af::get_alpha
const double * get_alpha() const
Returns the pointer on the array alpha.
Definition: map_af.C:477
Lorene::Tenseur::dec2_dzpuis
void dec2_dzpuis()
dzpuis -= 2 ;
Definition: tenseur.C:1130
Lorene::Bhole::taij_tot
Tenseur_sym taij_tot
Total , which must be zero on the horizon of the regularisation on the shift has been done.
Definition: bhole.h:305
Lorene::Cmp::std_base_scal
void std_base_scal()
Sets the spectral bases of the Valeur va to the standard ones for a scalar.
Definition: cmp.C:644
Lorene::Tenseur::set_std_base
void set_std_base()
Set the standard spectal basis of decomposition for each component.
Definition: tenseur.C:1170
Lorene::sin
Cmp sin(const Cmp &)
Sine.
Definition: cmp_math.C:69
Lorene::Map::ya
Coord ya
Absolute y coordinate.
Definition: map.h:731