LORENE
bin_ns_bh_kij.C
1 /*
2  * Methods for computing the extrinsic curvature tensor for a Bin_ns_bh
3  *
4  */
5 
6 /*
7  * Copyright (c) 2002 Philippe Grandclement, Keisuke Taniguchi,
8  * 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 version 2
14  * as published by the Free Software Foundation.
15  *
16  * LORENE is distributed in the hope that it will be useful,
17  * but WITHOUT ANY WARRANTY; without even the implied warranty of
18  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19  * GNU General Public License for more details.
20  *
21  * You should have received a copy of the GNU General Public License
22  * along with LORENE; if not, write to the Free Software
23  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
24  *
25  */
26 
27 char bin_ns_bh_kij_C[] = "$Header: /cvsroot/Lorene/C++/Source/Bin_ns_bh/bin_ns_bh_kij.C,v 1.11 2014/10/13 08:52:43 j_novak Exp $" ;
28 
29 /*
30  * $Id: bin_ns_bh_kij.C,v 1.11 2014/10/13 08:52:43 j_novak Exp $
31  * $Log: bin_ns_bh_kij.C,v $
32  * Revision 1.11 2014/10/13 08:52:43 j_novak
33  * Lorene classes and functions now belong to the namespace Lorene.
34  *
35  * Revision 1.10 2014/10/06 15:13:01 j_novak
36  * Modified #include directives to use c++ syntax.
37  *
38  * Revision 1.9 2008/09/26 08:44:04 p_grandclement
39  * Mixted binaries with non vanishing spin
40  *
41  * Revision 1.8 2008/08/19 06:41:59 j_novak
42  * Minor modifications to avoid warnings with gcc 4.3. Most of them concern
43  * cast-type operations, and constant strings that must be defined as const char*
44  *
45  * Revision 1.7 2007/04/26 14:14:59 f_limousin
46  * The function fait_tkij now have default values for bound_nn and lim_nn
47  *
48  * Revision 1.6 2007/04/26 13:16:23 f_limousin
49  * Correction of an error in the computation of grad_n_tot and grad_psi_tot
50  *
51  * Revision 1.5 2007/04/24 20:13:53 f_limousin
52  * Implementation of Dirichlet and Neumann BC for the lapse
53  *
54  * Revision 1.4 2005/12/01 12:59:10 p_grandclement
55  * Files for bin_ns_bh project
56  *
57  * Revision 1.3 2005/08/29 15:10:15 p_grandclement
58  * Addition of things needed :
59  * 1) For BBH with different masses
60  * 2) Provisory files for the mixted binaries (Bh and NS) : THIS IS NOT
61  * WORKING YET !!!
62  *
63  * Revision 1.2 2004/05/27 12:41:00 p_grandclement
64  * correction of some shadowed variables
65  *
66  * Revision 1.1 2003/02/13 16:40:25 p_grandclement
67  * Addition of various things for the Bin_ns_bh project, non of them being
68  * completely tested
69  *
70  *
71  * $Header: /cvsroot/Lorene/C++/Source/Bin_ns_bh/bin_ns_bh_kij.C,v 1.11 2014/10/13 08:52:43 j_novak Exp $
72  *
73  */
74 
75 // C++ headers
76 #include "headcpp.h"
77 
78 // C headers
79 #include <cmath>
80 
81 // Lorene headers
82 #include "bin_ns_bh.h"
83 #include "proto.h"
84 #include "utilitaires.h"
85 
86 #include "graphique.h"
87 
88 namespace Lorene {
89 
97 
98  Tenseur copie_shift (shift_auto) ;
99  copie_shift.change_triad(mp.get_bvect_cart()) ;
100 
101  if (shift_auto.get_etat() == ETATZERO)
103  else {
104  Tenseur grad (copie_shift.gradient()) ;
105  Tenseur trace (contract (grad,0, 1)) ;
106 
109  for (int i=0 ; i<3 ; i++)
110  for (int j=i ; j<3 ; j++)
111  taij_auto.set(i, j) = grad(i, j)+grad(j, i) ;
112  for (int i=0 ; i<3 ; i++)
113  taij_auto.set(i, i) -= 2./3.*trace() ;
114  }
115 
117 }
118 
120 
121  int nz_bh = hole.mp.get_mg()->get_nzone() ;
122  int nz_ns = star.mp.get_mg()->get_nzone() ;
123 
124  // On determine R_limite (pour le moment en tout cas...) :
125  double distance = star.mp.get_ori_x()-hole.mp.get_ori_x() ;
126  double lim_ns = distance/2. ;
127  double lim_bh = distance/2. ;
128  double int_ns = lim_ns/3. ;
129  double int_bh = lim_bh/3. ;
130 
131  // Les fonctions totales :
132  Cmp decouple_ns (star.mp) ;
133  decouple_ns.allocate_all() ;
134  Cmp decouple_bh (hole.mp) ;
135  decouple_bh.allocate_all() ;
136 
137  Mtbl xabs_ns (star.mp.xa) ;
138  Mtbl yabs_ns (star.mp.ya) ;
139  Mtbl zabs_ns (star.mp.za) ;
140 
141  Mtbl xabs_bh (hole.mp.xa) ;
142  Mtbl yabs_bh (hole.mp.ya) ;
143  Mtbl zabs_bh (hole.mp.za) ;
144 
145  double xabs, yabs, zabs, air_ns, air_bh, theta, phi ;
146 
147  // On boucle sur les autres zones :
148  for (int l=0 ; l<nz_ns ; l++) {
149  int nr = star.mp.get_mg()->get_nr (l) ;
150 
151  if (l==nz_ns-1)
152  nr -- ;
153 
154  int np = star.mp.get_mg()->get_np (l) ;
155  int nt = star.mp.get_mg()->get_nt (l) ;
156 
157  for (int k=0 ; k<np ; k++)
158  for (int j=0 ; j<nt ; j++)
159  for (int i=0 ; i<nr ; i++) {
160 
161  xabs = xabs_ns (l, k, j, i) ;
162  yabs = yabs_ns (l, k, j, i) ;
163  zabs = zabs_ns (l, k, j, i) ;
164 
165  // les coordonnees du point :
167  (xabs, yabs, zabs, air_ns, theta, phi) ;
169  (xabs, yabs, zabs, air_bh, theta, phi) ;
170 
171  if (air_ns <= lim_ns)
172  if (air_ns < int_ns)
173  decouple_ns.set(l, k, j, i) = 1 ;
174  else
175  decouple_ns.set(l, k, j, i) =
176  0.5*pow(cos((air_ns-int_ns)*M_PI/2./(lim_ns-int_ns)), 2.)+0.5
177  ;
178  else
179  if (air_bh <= lim_bh)
180  if (air_bh < int_bh)
181  decouple_ns.set(l, k, j, i) = 0 ;
182  else
183  decouple_ns.set(l, k, j, i) = 0.5*
184  pow(sin((air_bh-int_bh)*M_PI/2./(lim_bh-int_bh)), 2.)
185  ;
186  else
187  // On est loin des deux trous :
188  decouple_ns.set(l, k, j, i) = 0.5 ;
189  }
190 
191  // Cas infini :
192  if (l==nz_ns-1)
193  for (int k=0 ; k<np ; k++)
194  for (int j=0 ; j<nt ; j++)
195  decouple_ns.set(nz_ns-1, k, j, nr) = 0.5 ;
196  }
197 
198  for (int l=0 ; l<nz_bh ; l++) {
199  int nr = hole.mp.get_mg()->get_nr (l) ;
200 
201  if (l==nz_bh-1)
202  nr -- ;
203 
204  int np = hole.mp.get_mg()->get_np (l) ;
205  int nt = hole.mp.get_mg()->get_nt (l) ;
206 
207  for (int k=0 ; k<np ; k++)
208  for (int j=0 ; j<nt ; j++)
209  for (int i=0 ; i<nr ; i++) {
210 
211  xabs = xabs_bh (l, k, j, i) ;
212  yabs = yabs_bh (l, k, j, i) ;
213  zabs = zabs_bh (l, k, j, i) ;
214 
215  // les coordonnees du point :
217  (xabs, yabs, zabs, air_ns, theta, phi) ;
219  (xabs, yabs, zabs, air_bh, theta, phi) ;
220 
221  if (air_bh <= lim_bh)
222  if (air_bh < int_bh)
223  decouple_bh.set(l, k, j, i) = 1 ;
224  else
225  decouple_bh.set(l, k, j, i) = 0.5*
226  pow(cos((air_bh-int_bh)*M_PI/2./(lim_bh-int_bh)), 2.)+0.5 ;
227  else
228  if (air_ns <= lim_ns)
229  if (air_ns < int_ns)
230  decouple_bh.set(l, k, j, i) = 0 ;
231  else
232  decouple_bh.set(l, k, j, i) = 0.5*
233  pow(sin((air_ns-int_ns)*M_PI/2./(lim_ns-int_ns)), 2.) ;
234 
235  else
236  // On est loin des deux trous :
237  decouple_bh.set(l, k, j, i) = 0.5 ;
238  }
239 
240  // Cas infini :
241  if (l==nz_bh-1)
242  for (int k=0 ; k<np ; k++)
243  for (int j=0 ; j<nt ; j++)
244  decouple_bh.set(nz_bh-1, k, j, nr) = 0.5 ;
245  }
246 
247  decouple_ns.std_base_scal() ;
248  decouple_bh.std_base_scal() ;
249 
250  star.decouple = decouple_ns ;
251  hole.decouple = decouple_bh ;
252 }
253 
254 //********************************************************
255 //calcul de kij total. (la regularisation ayant ete faite)
256 //********************************************************
257 void Bin_ns_bh::fait_tkij (int bound_nn, double lim_nn) {
258 
259  fait_decouple() ;
260 
261  double norme_hole = 0 ;
262  double norme_star = 0 ;
263 
264  for (int i=0 ; i<3 ; i++) {
265  norme_hole += max(norme(hole.get_shift_auto()(i))) ;
266  norme_star += max(norme(star.get_shift_auto()(i))) ;
267  }
268 
269 #ifndef NDEBUG
270  bool zero_shift_hole = (norme_hole <1e-14) ? true : false ;
271 #endif
272  bool zero_shift_star = (norme_star <1e-14) ? true : false ;
273 
274  assert (zero_shift_hole == zero_shift_star) ;
275 
276  if (zero_shift_star == true) {
279 
283 
288  }
289  else {
290 
291  if (bound_nn < 0){
292  int nnt = hole.mp.get_mg()->get_nt(1) ;
293  int nnp = hole.mp.get_mg()->get_np(1) ;
294 
295  for (int k=0; k<nnp; k++)
296  for (int j=0; j<nnt; j++){
297  if (fabs((hole.n_auto+hole.n_comp)()(1, k, j , 0)) < 1e-2){
298  bound_nn = 0 ;
299  lim_nn = 0. ;
300  break ;
301  }
302  }
303  }
304 
305  if (bound_nn != 0 || lim_nn != 0){
306 
307  hole.fait_taij_auto () ;
308  star.fait_taij_auto() ;
309 
310  // On trouve les trucs du compagnon
312  // Pas de membre pour NS
313  Tenseur_sym ns_taij_comp (star.mp, 2, CON, ref_triad) ;
314  ns_taij_comp.set_etat_qcq() ;
315 
316  Tenseur_sym copie_ns (star.taij_auto) ;
317  copie_ns.dec2_dzpuis() ;
318  Tenseur_sym copie_bh (hole.taij_auto) ;
319  copie_bh.dec2_dzpuis() ;
320 
321  // Les importations :
322  if (hole.taij_auto.get_etat() == ETATZERO)
323  ns_taij_comp.set_etat_zero() ;
324  else {
325  ns_taij_comp.set(0, 0).import(copie_bh(0, 0)) ;
326  ns_taij_comp.set(0, 1).import(copie_bh(0, 1)) ;
327  ns_taij_comp.set(0, 2).import(copie_bh(0, 2)) ;
328  ns_taij_comp.set(1, 1).import(copie_bh(1, 1)) ;
329  ns_taij_comp.set(1, 2).import(copie_bh(1, 2)) ;
330  ns_taij_comp.set(2, 2).import(copie_bh(2, 2)) ;
331  ns_taij_comp.set_triad(*copie_bh.get_triad()) ;
332  ns_taij_comp.change_triad(star.ref_triad) ;
333  }
334 
335  if (star.taij_auto.get_etat() == ETATZERO)
337  else {
338  hole.taij_comp.set(0, 0).import(copie_ns(0, 0)) ;
339  hole.taij_comp.set(0, 1).import(copie_ns(0, 1)) ;
340  hole.taij_comp.set(0, 2).import(copie_ns(0, 2)) ;
341  hole.taij_comp.set(1, 1).import(copie_ns(1, 1)) ;
342  hole.taij_comp.set(1, 2).import(copie_ns(1, 2)) ;
343  hole.taij_comp.set(2, 2).import(copie_ns(2, 2)) ;
344  hole.taij_comp.set_triad(*copie_ns.get_triad()) ;
346  }
347 
349  ns_taij_comp.set_std_base() ;
351  ns_taij_comp.inc2_dzpuis() ;
352 
353  // Et enfin les trucs totaux...
355  Tenseur_sym ns_taij_tot (star.taij_auto + ns_taij_comp) ;
356  star.taij_tot = ns_taij_tot ;
357 
358  // Computation of tkij
364 
365  for (int i = 0 ; i<3 ; i++)
366  for (int j = i ; j<3 ; j++) {
367  star.tkij_tot.set(i,j) = 0.5*star.taij_tot(i,j)/star.nnn() ;
368  //star.tkij_auto.set(i,j) = 0.5*star.taij_tot(i,j)/star.nnn() ;
369  //star.tkij_comp.set(i,j) = 0.5*ns_taij_comp(i,j)/star.nnn() ;
370  hole.tkij_tot.set(i,j) = 0.5*hole.taij_tot(i,j)/hole.n_tot() ;
371  //hole.tkij_auto.set(i,j) = 0.5*hole.taij_auto(i,j)/hole.n_tot() ;
372  }
373 
374  for (int lig=0 ; lig<3 ; lig++)
375  for (int col=lig ; col<3 ; col++) {
376  star.tkij_auto.set(lig, col) = star.tkij_tot(lig, col)*
377  star.decouple ;
378  star.tkij_comp.set(lig, col) = star.tkij_tot(lig, col)*
379  (1-star.decouple) ;
380  hole.tkij_auto.set(lig, col) = hole.tkij_tot(lig, col)*
381  hole.decouple ;
382  }
386  }
387  else {
388 
391 
392  // On construit a_ij a partir du shift ...
393  // taij tot doit etre nul sur l'horizon.
394  hole.fait_taij_auto () ;
395  star.fait_taij_auto() ;
396 
397  // On trouve les trucs du compagnon
399  // Pas de membre pour NS
400  Tenseur_sym ns_taij_comp (star.mp, 2, CON, ref_triad) ;
401  ns_taij_comp.set_etat_qcq() ;
402 
403  Tenseur_sym copie_ns (star.taij_auto) ;
404  copie_ns.dec2_dzpuis() ;
405  Tenseur_sym copie_bh (hole.taij_auto) ;
406  copie_bh.dec2_dzpuis() ;
407 
408  // Les importations :
409  if (hole.taij_auto.get_etat() == ETATZERO)
410  ns_taij_comp.set_etat_zero() ;
411  else {
412  ns_taij_comp.set(0, 0).import_asymy(copie_bh(0, 0)) ;
413  ns_taij_comp.set(0, 1).import_symy(copie_bh(0, 1)) ;
414  ns_taij_comp.set(0, 2).import_asymy(copie_bh(0, 2)) ;
415  ns_taij_comp.set(1, 1).import_asymy(copie_bh(1, 1)) ;
416  ns_taij_comp.set(1, 2).import_symy(copie_bh(1, 2)) ;
417  ns_taij_comp.set(2, 2).import_asymy(copie_bh(2, 2)) ;
418  ns_taij_comp.set_triad(*copie_bh.get_triad()) ;
419  ns_taij_comp.change_triad(star.ref_triad) ;
420  }
421 
422  if (star.taij_auto.get_etat() == ETATZERO)
424  else {
425  hole.taij_comp.set(0, 0).import_asymy(copie_ns(0, 0)) ;
426  hole.taij_comp.set(0, 1).import_symy(copie_ns(0, 1)) ;
427  hole.taij_comp.set(0, 2).import_asymy(copie_ns(0, 2)) ;
428  hole.taij_comp.set(1, 1).import_asymy(copie_ns(1, 1)) ;
429  hole.taij_comp.set(1, 2).import_symy(copie_ns(1, 2)) ;
430  hole.taij_comp.set(2, 2).import_asymy(copie_ns(2, 2)) ;
431  hole.taij_comp.set_triad(*copie_ns.get_triad()) ;
433  }
434 
436  ns_taij_comp.set_std_base() ;
438  ns_taij_comp.inc2_dzpuis() ;
439 
440  // Et enfin les trucs totaux...
442  Tenseur_sym ns_taij_tot (star.taij_auto + ns_taij_comp) ;
443  star.taij_tot = ns_taij_tot ;
444 
445  int nz_ns = star.mp.get_mg()->get_nzone() ;
446  Cmp ntot_ns (star.get_nnn()()) ;
447 
448  Cmp ntot_bh (hole.n_tot()) ;
449  //des_coupe_z (ntot_bh, 0, -5, 0, -2.5, 2.5, "N") ;
450  ntot_bh = division_xpun (ntot_bh, 0) ;
451  ntot_bh.raccord(1) ;
452 
453  //des_coupe_z (ntot_bh, 0, -5, 0, -2.5, 2.5, "N/ (x+1)") ;
454 
455  // Boucle sur les composantes :
456  for (int lig = 0 ; lig<3 ; lig++)
457  for (int col = lig ; col<3 ; col++) {
458 
459  // Dans la grille du BH (pas de pb sauf pres horizon :
460  Cmp auxi_bh (hole.taij_tot(lig, col)/2.) ;
461  auxi_bh.dec2_dzpuis() ;
462  auxi_bh = division_xpun (auxi_bh, 0) ;
463 
464  //des_coupe_z (auxi_bh, 0, -10, 20, -7, 7, "Aij/ (x+1)") ;
465  auxi_bh = auxi_bh / ntot_bh ;
466  auxi_bh.raccord(1) ;
467 
468  // Pour la NS on doit faire attention a pas etre pres du trou
469  Cmp auxi_ns (star.mp) ;
470  auxi_ns.allocate_all() ;
471 
472  // copie :
473  Cmp copie_ns_bis (ns_taij_tot(lig, col)) ;
474  copie_ns_bis.dec2_dzpuis() ;
475 
476  // Double le rayon limite :
477  double lim_bh = hole.mp.get_alpha()[1] + hole.mp.get_beta()[1] ;
478 
479  Mtbl xabs_ns (star.mp.xa) ;
480  Mtbl yabs_ns (star.mp.ya) ;
481  Mtbl zabs_ns (star.mp.za) ;
482 
483  double xabs, yabs, zabs, air, theta, phi ;
484 
485  // On boucle sur les zones
486  for (int l=0 ; l<nz_ns ; l++) {
487 
488  int nr = star.mp.get_mg()->get_nr (l) ;
489 
490  if (l==nz_ns-1)
491  nr -- ;
492 
493  int np = star.mp.get_mg()->get_np (l) ;
494  int nt = star.mp.get_mg()->get_nt (l) ;
495 
496  for (int k=0 ; k<np ; k++)
497  for (int j=0 ; j<nt ; j++)
498  for (int i=0 ; i<nr ; i++) {
499 
500  xabs = xabs_ns (l, k, j, i) ;
501  yabs = yabs_ns (l, k, j, i) ;
502  zabs = zabs_ns (l, k, j, i) ;
503 
504  // les coordonnees du point vis a vis BH :
506  (xabs, yabs, zabs, air, theta, phi) ;
507 
508  if (air >= lim_bh)
509  // On est loin du trou 2 : pas de pb :
510  auxi_ns.set(l, k, j, i) =
511  copie_ns_bis(l, k, j, i) / ntot_ns (l, k, j, i)/2. ;
512  else
513  // On est pres du trou (faut pas tomber dedans !) :
514  auxi_ns.set(l, k, j, i) = auxi_bh.val_point (air, theta, phi) ;
515  }
516 
517  // Cas infini :
518  if (l==nz_ns-1)
519  for (int k=0 ; k<np ; k++)
520  for (int j=0 ; j<nt ; j++)
521  auxi_ns.set(nz_ns-1, k, j, nr) = 0 ;
522  }
523 
524 
525  star.tkij_tot.set(lig, col) = auxi_ns ;
526  hole.tkij_tot.set(lig, col) = auxi_bh ;
527  }
528 
532 
534 
535  //Cmp dessin_un (hole.get_tkij_tot()(0,0)) ;
536  //des_coupe_z (dessin_un, 0, -10, 20, -7, 7, "Kij tot BH") ;
537  //Cmp dessin_deux (ns_tkij_tot(0,0)) ;
538  //des_coupe_z (dessin_deux, 0, -10, 20, -7, 7, "Kij tot NS") ;
539 
543 
544  for (int lig=0 ; lig<3 ; lig++)
545  for (int col=lig ; col<3 ; col++) {
546  star.tkij_auto.set(lig, col) = star.tkij_tot(lig, col)*
547  star.decouple ;
548  star.tkij_comp.set(lig, col) = star.tkij_tot(lig, col)*
549  (1-star.decouple) ;
550  hole.tkij_auto.set(lig, col) = hole.tkij_tot(lig, col)*
551  hole.decouple ;
552  }
556 
557  }
558 
559  // On doit mettre a jour les champs akcar de NS :
561  star.akcar_auto.set() = 0 ;
562 
563  for (int i=0; i<3; i++)
564  for (int j=0; j<3; j++)
565  star.akcar_auto.set() += star.tkij_auto(i, j) % star.tkij_auto(i, j) ;
566 
569 
571  star.akcar_comp.set() = 0 ;
572 
573  for (int i=0; i<3; i++)
574  for (int j=0; j<3; j++)
575  star.akcar_comp.set() += star.tkij_auto(i, j) % star.tkij_comp(i, j) ;
576 
579  }
580 }
581 }
Lorene::Etoile_bin::decouple
Cmp decouple
Function used to construct the part of generated by the star from the total .
Definition: etoile.h:1000
Lorene::Etoile_bin::akcar_comp
Tenseur akcar_comp
Part of the scalar generated by shift_auto and shift_comp , i.e.
Definition: etoile.h:944
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::Etoile::a_car
Tenseur a_car
Total conformal factor .
Definition: etoile.h:515
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::Etoile_bin::shift_auto
Tenseur shift_auto
Part of the shift vector generated principaly by the star.
Definition: etoile.h:889
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::Etoile_bin::akcar_auto
Tenseur akcar_auto
Part of the scalar generated by shift_auto , i.e.
Definition: etoile.h:938
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::Tenseur::gradient
const Tenseur & gradient() const
Returns the gradient of *this (Cartesian coordinates)
Definition: tenseur.C:1542
Lorene::Et_bin_nsbh::taij_auto
Tenseur_sym taij_auto
Part of the extrinsic curvature tensor $\tilde A^{ij} = 2 N K^{ij}$ generated by {\tt shift_auto}.
Definition: et_bin_nsbh.h:126
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::Etoile_bin::tkij_comp
Tenseur_sym tkij_comp
Part of the extrinsic curvature tensor generated by shift_comp .
Definition: etoile.h:932
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::Etoile::mp
Map & mp
Mapping associated with the star.
Definition: etoile.h:429
Lorene::norme
Tbl norme(const Cmp &)
Sums of the absolute values of all the values of the Cmp in each domain.
Definition: cmp_math.C:481
Lorene::max
Tbl max(const Cmp &)
Maximum values of a Cmp in each domain.
Definition: cmp_math.C:435
Lorene::Bhole::n_comp
Tenseur n_comp
Part of N generated by the companion hole.
Definition: bhole.h:287
Lorene::Tenseur::get_triad
const Base_vect * get_triad() const
Returns the vectorial basis (triad) on which the components are defined.
Definition: tenseur.h:704
Lorene::Et_bin_nsbh::taij_tot
Tenseur_sym taij_tot
Total extrinsic curvature tensor $\tilde A^{ij} = 2 N K^{ij}$ generated by {\tt shift_auto} and {\tt ...
Definition: et_bin_nsbh.h:140
Lorene::Bhole::n_tot
Tenseur n_tot
Total N .
Definition: bhole.h:288
Lorene::Tenseur::set_triad
void set_triad(const Base_vect &new_triad)
Assigns a new vectorial basis (triad) of decomposition.
Definition: tenseur.C:674
Lorene::pow
Cmp pow(const Cmp &, int)
Power .
Definition: cmp_math.C:348
Lorene::Etoile_bin::ref_triad
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
Lorene::Etoile_bin::get_shift_auto
const Tenseur & get_shift_auto() const
Returns the part of the shift vector generated principaly by the star.
Definition: etoile.h:1152
Lorene::Etoile::nnn
Tenseur nnn
Total lapse function.
Definition: etoile.h:509
Lorene::Tenseur::get_etat
int get_etat() const
Returns the logical state.
Definition: tenseur.h:707
Lorene::Et_bin_nsbh::tkij_tot
Tenseur_sym tkij_tot
Total extrinsic curvature tensor $K^{ij}$ generated by {\tt shift_auto} and {\tt shift_comp}.
Definition: et_bin_nsbh.h:152
Lorene::Bhole::get_shift_auto
const Tenseur & get_shift_auto() const
Returns the part of generated by the hole.
Definition: bhole.h:440
Lorene::Bin_ns_bh::hole
Bhole hole
The black hole.
Definition: bin_ns_bh.h:131
Lorene::Bhole::n_auto
Tenseur n_auto
Part of N generated by the hole.
Definition: bhole.h:286
Lorene::Tenseur_sym
Class intended to describe tensors with a symmetry on the two last indices *** DEPRECATED : use class...
Definition: tenseur.h:1253
Lorene::Cmp
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition: cmp.h:446
Lorene::Et_bin_nsbh::fait_taij_auto
void fait_taij_auto()
Computes (LB)^{ij} auto.
Definition: bin_ns_bh_kij.C:96
Lorene::Bin_ns_bh::fait_tkij
void fait_tkij(int bound_nn=-1, double lim_nn=0)
Computation of the extrinsic curvature tensor for both {\tt star} and {\tt bhole}.
Definition: bin_ns_bh_kij.C:257
Lorene::Tenseur::change_triad
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
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::Map::get_bvect_cart
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
Lorene::Bhole::tkij_auto
Tenseur_sym tkij_auto
Auto .
Definition: bhole.h:307
Lorene::Bin_ns_bh::ref_triad
const Base_vect_cart ref_triad
Cartesian triad of the absolute reference frame.
Definition: bin_ns_bh.h:125
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::Bin_ns_bh::star
Et_bin_nsbh star
The neutron star.
Definition: bin_ns_bh.h:128
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::Bin_ns_bh::fait_decouple
void fait_decouple()
Function used to compute the {\tt decouple} functions for both the NS and the BH.
Definition: bin_ns_bh_kij.C:119
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::Cmp::import
void import(const Cmp &ci)
Assignment to another Cmp defined on a different mapping.
Definition: cmp_import.C:73
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::Etoile::get_nnn
const Tenseur & get_nnn() const
Returns the total lapse function N.
Definition: etoile.h:727
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::Et_bin_nsbh::tkij_auto
Tenseur_sym tkij_auto
Part of the extrinsic curvature tensor $K^{ij}$ generated by {\tt shift_auto}.
Definition: et_bin_nsbh.h:146
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::contract
Tenseur contract(const Tenseur &, int id1, int id2)
Self contraction of two indices of a Tenseur .
Definition: tenseur_operateur.C:279
Lorene::Map::ya
Coord ya
Absolute y coordinate.
Definition: map.h:731