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 
28 
29 /*
30  * $Id: bin_ns_bh_kij.C,v 1.12 2016/12/05 16:17:46 j_novak Exp $
31  * $Log: bin_ns_bh_kij.C,v $
32  * Revision 1.12 2016/12/05 16:17:46 j_novak
33  * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
34  *
35  * Revision 1.11 2014/10/13 08:52:43 j_novak
36  * Lorene classes and functions now belong to the namespace Lorene.
37  *
38  * Revision 1.10 2014/10/06 15:13:01 j_novak
39  * Modified #include directives to use c++ syntax.
40  *
41  * Revision 1.9 2008/09/26 08:44:04 p_grandclement
42  * Mixted binaries with non vanishing spin
43  *
44  * Revision 1.8 2008/08/19 06:41:59 j_novak
45  * Minor modifications to avoid warnings with gcc 4.3. Most of them concern
46  * cast-type operations, and constant strings that must be defined as const char*
47  *
48  * Revision 1.7 2007/04/26 14:14:59 f_limousin
49  * The function fait_tkij now have default values for bound_nn and lim_nn
50  *
51  * Revision 1.6 2007/04/26 13:16:23 f_limousin
52  * Correction of an error in the computation of grad_n_tot and grad_psi_tot
53  *
54  * Revision 1.5 2007/04/24 20:13:53 f_limousin
55  * Implementation of Dirichlet and Neumann BC for the lapse
56  *
57  * Revision 1.4 2005/12/01 12:59:10 p_grandclement
58  * Files for bin_ns_bh project
59  *
60  * Revision 1.3 2005/08/29 15:10:15 p_grandclement
61  * Addition of things needed :
62  * 1) For BBH with different masses
63  * 2) Provisory files for the mixted binaries (Bh and NS) : THIS IS NOT
64  * WORKING YET !!!
65  *
66  * Revision 1.2 2004/05/27 12:41:00 p_grandclement
67  * correction of some shadowed variables
68  *
69  * Revision 1.1 2003/02/13 16:40:25 p_grandclement
70  * Addition of various things for the Bin_ns_bh project, non of them being
71  * completely tested
72  *
73  *
74  * $Header: /cvsroot/Lorene/C++/Source/Bin_ns_bh/bin_ns_bh_kij.C,v 1.12 2016/12/05 16:17:46 j_novak Exp $
75  *
76  */
77 
78 // C++ headers
79 #include "headcpp.h"
80 
81 // C headers
82 #include <cmath>
83 
84 // Lorene headers
85 #include "bin_ns_bh.h"
86 #include "proto.h"
87 #include "utilitaires.h"
88 
89 #include "graphique.h"
90 
91 namespace Lorene {
92 
100 
101  Tenseur copie_shift (shift_auto) ;
102  copie_shift.change_triad(mp.get_bvect_cart()) ;
103 
104  if (shift_auto.get_etat() == ETATZERO)
106  else {
107  Tenseur grad (copie_shift.gradient()) ;
108  Tenseur trace (contract (grad,0, 1)) ;
109 
112  for (int i=0 ; i<3 ; i++)
113  for (int j=i ; j<3 ; j++)
114  taij_auto.set(i, j) = grad(i, j)+grad(j, i) ;
115  for (int i=0 ; i<3 ; i++)
116  taij_auto.set(i, i) -= 2./3.*trace() ;
117  }
118 
120 }
121 
123 
124  int nz_bh = hole.mp.get_mg()->get_nzone() ;
125  int nz_ns = star.mp.get_mg()->get_nzone() ;
126 
127  // On determine R_limite (pour le moment en tout cas...) :
128  double distance = star.mp.get_ori_x()-hole.mp.get_ori_x() ;
129  double lim_ns = distance/2. ;
130  double lim_bh = distance/2. ;
131  double int_ns = lim_ns/3. ;
132  double int_bh = lim_bh/3. ;
133 
134  // Les fonctions totales :
135  Cmp decouple_ns (star.mp) ;
136  decouple_ns.allocate_all() ;
137  Cmp decouple_bh (hole.mp) ;
138  decouple_bh.allocate_all() ;
139 
140  Mtbl xabs_ns (star.mp.xa) ;
141  Mtbl yabs_ns (star.mp.ya) ;
142  Mtbl zabs_ns (star.mp.za) ;
143 
144  Mtbl xabs_bh (hole.mp.xa) ;
145  Mtbl yabs_bh (hole.mp.ya) ;
146  Mtbl zabs_bh (hole.mp.za) ;
147 
148  double xabs, yabs, zabs, air_ns, air_bh, theta, phi ;
149 
150  // On boucle sur les autres zones :
151  for (int l=0 ; l<nz_ns ; l++) {
152  int nr = star.mp.get_mg()->get_nr (l) ;
153 
154  if (l==nz_ns-1)
155  nr -- ;
156 
157  int np = star.mp.get_mg()->get_np (l) ;
158  int nt = star.mp.get_mg()->get_nt (l) ;
159 
160  for (int k=0 ; k<np ; k++)
161  for (int j=0 ; j<nt ; j++)
162  for (int i=0 ; i<nr ; i++) {
163 
164  xabs = xabs_ns (l, k, j, i) ;
165  yabs = yabs_ns (l, k, j, i) ;
166  zabs = zabs_ns (l, k, j, i) ;
167 
168  // les coordonnees du point :
170  (xabs, yabs, zabs, air_ns, theta, phi) ;
172  (xabs, yabs, zabs, air_bh, theta, phi) ;
173 
174  if (air_ns <= lim_ns)
175  if (air_ns < int_ns)
176  decouple_ns.set(l, k, j, i) = 1 ;
177  else
178  decouple_ns.set(l, k, j, i) =
179  0.5*pow(cos((air_ns-int_ns)*M_PI/2./(lim_ns-int_ns)), 2.)+0.5
180  ;
181  else
182  if (air_bh <= lim_bh)
183  if (air_bh < int_bh)
184  decouple_ns.set(l, k, j, i) = 0 ;
185  else
186  decouple_ns.set(l, k, j, i) = 0.5*
187  pow(sin((air_bh-int_bh)*M_PI/2./(lim_bh-int_bh)), 2.)
188  ;
189  else
190  // On est loin des deux trous :
191  decouple_ns.set(l, k, j, i) = 0.5 ;
192  }
193 
194  // Cas infini :
195  if (l==nz_ns-1)
196  for (int k=0 ; k<np ; k++)
197  for (int j=0 ; j<nt ; j++)
198  decouple_ns.set(nz_ns-1, k, j, nr) = 0.5 ;
199  }
200 
201  for (int l=0 ; l<nz_bh ; l++) {
202  int nr = hole.mp.get_mg()->get_nr (l) ;
203 
204  if (l==nz_bh-1)
205  nr -- ;
206 
207  int np = hole.mp.get_mg()->get_np (l) ;
208  int nt = hole.mp.get_mg()->get_nt (l) ;
209 
210  for (int k=0 ; k<np ; k++)
211  for (int j=0 ; j<nt ; j++)
212  for (int i=0 ; i<nr ; i++) {
213 
214  xabs = xabs_bh (l, k, j, i) ;
215  yabs = yabs_bh (l, k, j, i) ;
216  zabs = zabs_bh (l, k, j, i) ;
217 
218  // les coordonnees du point :
220  (xabs, yabs, zabs, air_ns, theta, phi) ;
222  (xabs, yabs, zabs, air_bh, theta, phi) ;
223 
224  if (air_bh <= lim_bh)
225  if (air_bh < int_bh)
226  decouple_bh.set(l, k, j, i) = 1 ;
227  else
228  decouple_bh.set(l, k, j, i) = 0.5*
229  pow(cos((air_bh-int_bh)*M_PI/2./(lim_bh-int_bh)), 2.)+0.5 ;
230  else
231  if (air_ns <= lim_ns)
232  if (air_ns < int_ns)
233  decouple_bh.set(l, k, j, i) = 0 ;
234  else
235  decouple_bh.set(l, k, j, i) = 0.5*
236  pow(sin((air_ns-int_ns)*M_PI/2./(lim_ns-int_ns)), 2.) ;
237 
238  else
239  // On est loin des deux trous :
240  decouple_bh.set(l, k, j, i) = 0.5 ;
241  }
242 
243  // Cas infini :
244  if (l==nz_bh-1)
245  for (int k=0 ; k<np ; k++)
246  for (int j=0 ; j<nt ; j++)
247  decouple_bh.set(nz_bh-1, k, j, nr) = 0.5 ;
248  }
249 
250  decouple_ns.std_base_scal() ;
251  decouple_bh.std_base_scal() ;
252 
253  star.decouple = decouple_ns ;
254  hole.decouple = decouple_bh ;
255 }
256 
257 //********************************************************
258 //calcul de kij total. (la regularisation ayant ete faite)
259 //********************************************************
260 void Bin_ns_bh::fait_tkij (int bound_nn, double lim_nn) {
261 
262  fait_decouple() ;
263 
264  double norme_hole = 0 ;
265  double norme_star = 0 ;
266 
267  for (int i=0 ; i<3 ; i++) {
268  norme_hole += max(norme(hole.get_shift_auto()(i))) ;
269  norme_star += max(norme(star.get_shift_auto()(i))) ;
270  }
271 
272 #ifndef NDEBUG
273  bool zero_shift_hole = (norme_hole <1e-14) ? true : false ;
274 #endif
275  bool zero_shift_star = (norme_star <1e-14) ? true : false ;
276 
277  assert (zero_shift_hole == zero_shift_star) ;
278 
279  if (zero_shift_star == true) {
282 
286 
291  }
292  else {
293 
294  if (bound_nn < 0){
295  int nnt = hole.mp.get_mg()->get_nt(1) ;
296  int nnp = hole.mp.get_mg()->get_np(1) ;
297 
298  for (int k=0; k<nnp; k++)
299  for (int j=0; j<nnt; j++){
300  if (fabs((hole.n_auto+hole.n_comp)()(1, k, j , 0)) < 1e-2){
301  bound_nn = 0 ;
302  lim_nn = 0. ;
303  break ;
304  }
305  }
306  }
307 
308  if (bound_nn != 0 || lim_nn != 0){
309 
310  hole.fait_taij_auto () ;
311  star.fait_taij_auto() ;
312 
313  // On trouve les trucs du compagnon
315  // Pas de membre pour NS
316  Tenseur_sym ns_taij_comp (star.mp, 2, CON, ref_triad) ;
317  ns_taij_comp.set_etat_qcq() ;
318 
319  Tenseur_sym copie_ns (star.taij_auto) ;
320  copie_ns.dec2_dzpuis() ;
321  Tenseur_sym copie_bh (hole.taij_auto) ;
322  copie_bh.dec2_dzpuis() ;
323 
324  // Les importations :
325  if (hole.taij_auto.get_etat() == ETATZERO)
326  ns_taij_comp.set_etat_zero() ;
327  else {
328  ns_taij_comp.set(0, 0).import(copie_bh(0, 0)) ;
329  ns_taij_comp.set(0, 1).import(copie_bh(0, 1)) ;
330  ns_taij_comp.set(0, 2).import(copie_bh(0, 2)) ;
331  ns_taij_comp.set(1, 1).import(copie_bh(1, 1)) ;
332  ns_taij_comp.set(1, 2).import(copie_bh(1, 2)) ;
333  ns_taij_comp.set(2, 2).import(copie_bh(2, 2)) ;
334  ns_taij_comp.set_triad(*copie_bh.get_triad()) ;
335  ns_taij_comp.change_triad(star.ref_triad) ;
336  }
337 
338  if (star.taij_auto.get_etat() == ETATZERO)
340  else {
341  hole.taij_comp.set(0, 0).import(copie_ns(0, 0)) ;
342  hole.taij_comp.set(0, 1).import(copie_ns(0, 1)) ;
343  hole.taij_comp.set(0, 2).import(copie_ns(0, 2)) ;
344  hole.taij_comp.set(1, 1).import(copie_ns(1, 1)) ;
345  hole.taij_comp.set(1, 2).import(copie_ns(1, 2)) ;
346  hole.taij_comp.set(2, 2).import(copie_ns(2, 2)) ;
347  hole.taij_comp.set_triad(*copie_ns.get_triad()) ;
349  }
350 
352  ns_taij_comp.set_std_base() ;
354  ns_taij_comp.inc2_dzpuis() ;
355 
356  // Et enfin les trucs totaux...
358  Tenseur_sym ns_taij_tot (star.taij_auto + ns_taij_comp) ;
359  star.taij_tot = ns_taij_tot ;
360 
361  // Computation of tkij
367 
368  for (int i = 0 ; i<3 ; i++)
369  for (int j = i ; j<3 ; j++) {
370  star.tkij_tot.set(i,j) = 0.5*star.taij_tot(i,j)/star.nnn() ;
371  //star.tkij_auto.set(i,j) = 0.5*star.taij_tot(i,j)/star.nnn() ;
372  //star.tkij_comp.set(i,j) = 0.5*ns_taij_comp(i,j)/star.nnn() ;
373  hole.tkij_tot.set(i,j) = 0.5*hole.taij_tot(i,j)/hole.n_tot() ;
374  //hole.tkij_auto.set(i,j) = 0.5*hole.taij_auto(i,j)/hole.n_tot() ;
375  }
376 
377  for (int lig=0 ; lig<3 ; lig++)
378  for (int col=lig ; col<3 ; col++) {
379  star.tkij_auto.set(lig, col) = star.tkij_tot(lig, col)*
380  star.decouple ;
381  star.tkij_comp.set(lig, col) = star.tkij_tot(lig, col)*
382  (1-star.decouple) ;
383  hole.tkij_auto.set(lig, col) = hole.tkij_tot(lig, col)*
384  hole.decouple ;
385  }
389  }
390  else {
391 
394 
395  // On construit a_ij a partir du shift ...
396  // taij tot doit etre nul sur l'horizon.
397  hole.fait_taij_auto () ;
398  star.fait_taij_auto() ;
399 
400  // On trouve les trucs du compagnon
402  // Pas de membre pour NS
403  Tenseur_sym ns_taij_comp (star.mp, 2, CON, ref_triad) ;
404  ns_taij_comp.set_etat_qcq() ;
405 
406  Tenseur_sym copie_ns (star.taij_auto) ;
407  copie_ns.dec2_dzpuis() ;
408  Tenseur_sym copie_bh (hole.taij_auto) ;
409  copie_bh.dec2_dzpuis() ;
410 
411  // Les importations :
412  if (hole.taij_auto.get_etat() == ETATZERO)
413  ns_taij_comp.set_etat_zero() ;
414  else {
415  ns_taij_comp.set(0, 0).import_asymy(copie_bh(0, 0)) ;
416  ns_taij_comp.set(0, 1).import_symy(copie_bh(0, 1)) ;
417  ns_taij_comp.set(0, 2).import_asymy(copie_bh(0, 2)) ;
418  ns_taij_comp.set(1, 1).import_asymy(copie_bh(1, 1)) ;
419  ns_taij_comp.set(1, 2).import_symy(copie_bh(1, 2)) ;
420  ns_taij_comp.set(2, 2).import_asymy(copie_bh(2, 2)) ;
421  ns_taij_comp.set_triad(*copie_bh.get_triad()) ;
422  ns_taij_comp.change_triad(star.ref_triad) ;
423  }
424 
425  if (star.taij_auto.get_etat() == ETATZERO)
427  else {
428  hole.taij_comp.set(0, 0).import_asymy(copie_ns(0, 0)) ;
429  hole.taij_comp.set(0, 1).import_symy(copie_ns(0, 1)) ;
430  hole.taij_comp.set(0, 2).import_asymy(copie_ns(0, 2)) ;
431  hole.taij_comp.set(1, 1).import_asymy(copie_ns(1, 1)) ;
432  hole.taij_comp.set(1, 2).import_symy(copie_ns(1, 2)) ;
433  hole.taij_comp.set(2, 2).import_asymy(copie_ns(2, 2)) ;
434  hole.taij_comp.set_triad(*copie_ns.get_triad()) ;
436  }
437 
439  ns_taij_comp.set_std_base() ;
441  ns_taij_comp.inc2_dzpuis() ;
442 
443  // Et enfin les trucs totaux...
445  Tenseur_sym ns_taij_tot (star.taij_auto + ns_taij_comp) ;
446  star.taij_tot = ns_taij_tot ;
447 
448  int nz_ns = star.mp.get_mg()->get_nzone() ;
449  Cmp ntot_ns (star.get_nnn()()) ;
450 
451  Cmp ntot_bh (hole.n_tot()) ;
452  //des_coupe_z (ntot_bh, 0, -5, 0, -2.5, 2.5, "N") ;
453  ntot_bh = division_xpun (ntot_bh, 0) ;
454  ntot_bh.raccord(1) ;
455 
456  //des_coupe_z (ntot_bh, 0, -5, 0, -2.5, 2.5, "N/ (x+1)") ;
457 
458  // Boucle sur les composantes :
459  for (int lig = 0 ; lig<3 ; lig++)
460  for (int col = lig ; col<3 ; col++) {
461 
462  // Dans la grille du BH (pas de pb sauf pres horizon :
463  Cmp auxi_bh (hole.taij_tot(lig, col)/2.) ;
464  auxi_bh.dec2_dzpuis() ;
465  auxi_bh = division_xpun (auxi_bh, 0) ;
466 
467  //des_coupe_z (auxi_bh, 0, -10, 20, -7, 7, "Aij/ (x+1)") ;
468  auxi_bh = auxi_bh / ntot_bh ;
469  auxi_bh.raccord(1) ;
470 
471  // Pour la NS on doit faire attention a pas etre pres du trou
472  Cmp auxi_ns (star.mp) ;
473  auxi_ns.allocate_all() ;
474 
475  // copie :
476  Cmp copie_ns_bis (ns_taij_tot(lig, col)) ;
477  copie_ns_bis.dec2_dzpuis() ;
478 
479  // Double le rayon limite :
480  double lim_bh = hole.mp.get_alpha()[1] + hole.mp.get_beta()[1] ;
481 
482  Mtbl xabs_ns (star.mp.xa) ;
483  Mtbl yabs_ns (star.mp.ya) ;
484  Mtbl zabs_ns (star.mp.za) ;
485 
486  double xabs, yabs, zabs, air, theta, phi ;
487 
488  // On boucle sur les zones
489  for (int l=0 ; l<nz_ns ; l++) {
490 
491  int nr = star.mp.get_mg()->get_nr (l) ;
492 
493  if (l==nz_ns-1)
494  nr -- ;
495 
496  int np = star.mp.get_mg()->get_np (l) ;
497  int nt = star.mp.get_mg()->get_nt (l) ;
498 
499  for (int k=0 ; k<np ; k++)
500  for (int j=0 ; j<nt ; j++)
501  for (int i=0 ; i<nr ; i++) {
502 
503  xabs = xabs_ns (l, k, j, i) ;
504  yabs = yabs_ns (l, k, j, i) ;
505  zabs = zabs_ns (l, k, j, i) ;
506 
507  // les coordonnees du point vis a vis BH :
509  (xabs, yabs, zabs, air, theta, phi) ;
510 
511  if (air >= lim_bh)
512  // On est loin du trou 2 : pas de pb :
513  auxi_ns.set(l, k, j, i) =
514  copie_ns_bis(l, k, j, i) / ntot_ns (l, k, j, i)/2. ;
515  else
516  // On est pres du trou (faut pas tomber dedans !) :
517  auxi_ns.set(l, k, j, i) = auxi_bh.val_point (air, theta, phi) ;
518  }
519 
520  // Cas infini :
521  if (l==nz_ns-1)
522  for (int k=0 ; k<np ; k++)
523  for (int j=0 ; j<nt ; j++)
524  auxi_ns.set(nz_ns-1, k, j, nr) = 0 ;
525  }
526 
527 
528  star.tkij_tot.set(lig, col) = auxi_ns ;
529  hole.tkij_tot.set(lig, col) = auxi_bh ;
530  }
531 
535 
537 
538  //Cmp dessin_un (hole.get_tkij_tot()(0,0)) ;
539  //des_coupe_z (dessin_un, 0, -10, 20, -7, 7, "Kij tot BH") ;
540  //Cmp dessin_deux (ns_tkij_tot(0,0)) ;
541  //des_coupe_z (dessin_deux, 0, -10, 20, -7, 7, "Kij tot NS") ;
542 
546 
547  for (int lig=0 ; lig<3 ; lig++)
548  for (int col=lig ; col<3 ; col++) {
549  star.tkij_auto.set(lig, col) = star.tkij_tot(lig, col)*
550  star.decouple ;
551  star.tkij_comp.set(lig, col) = star.tkij_tot(lig, col)*
552  (1-star.decouple) ;
553  hole.tkij_auto.set(lig, col) = hole.tkij_tot(lig, col)*
554  hole.decouple ;
555  }
559 
560  }
561 
562  // On doit mettre a jour les champs akcar de NS :
564  star.akcar_auto.set() = 0 ;
565 
566  for (int i=0; i<3; i++)
567  for (int j=0; j<3; j++)
568  star.akcar_auto.set() += star.tkij_auto(i, j) % star.tkij_auto(i, j) ;
569 
572 
574  star.akcar_comp.set() = 0 ;
575 
576  for (int i=0; i<3; i++)
577  for (int j=0; j<3; j++)
578  star.akcar_comp.set() += star.tkij_auto(i, j) % star.tkij_comp(i, j) ;
579 
582  }
583 }
584 }
Coord xa
Absolute x coordinate.
Definition: map.h:748
const Base_vect & ref_triad
Reference triad ("absolute frame"), with respect to which the components of all the member Tenseur &#39;s...
Definition: etoile.h:831
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition: cmp.h:446
const double * get_alpha() const
Returns the pointer on the array alpha.
Definition: map_af.C:607
Tenseur n_tot
Total N .
Definition: bhole.h:288
void dec2_dzpuis()
dzpuis -= 2 ;
Definition: tenseur.C:1146
void set_triad(const Base_vect &new_triad)
Assigns a new vectorial basis (triad) of decomposition.
Definition: tenseur.C:690
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
Definition: grilles.h:479
void set_std_base()
Set the standard spectal basis of decomposition for each component.
Definition: tenseur.C:1186
Multi-domain array.
Definition: mtbl.h:118
Class intended to describe tensors with a symmetry on the two last indices *** DEPRECATED : use class...
Definition: tenseur.h:1256
Lorene prototypes.
Definition: app_hor.h:67
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:783
Tenseur nnn
Total lapse function.
Definition: etoile.h:512
double get_ori_x() const
Returns the x coordinate of the origin.
Definition: map.h:786
Tenseur n_comp
Part of N generated by the companion hole.
Definition: bhole.h:287
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
void fait_taij_auto()
Calculates the part of generated by shift_auto .
Definition: bhole.C:382
Tenseur_sym taij_auto
Part of generated by the hole.
Definition: bhole.h:299
Cmp cos(const Cmp &)
Cosine.
Definition: cmp_math.C:97
void inc2_dzpuis()
dzpuis += 2 ;
Definition: tenseur.C:1159
Tenseur shift_auto
Part of the shift vector generated principaly by the star.
Definition: etoile.h:892
Tenseur_sym taij_comp
Part of generated by the companion hole.
Definition: bhole.h:300
const Base_vect_cart ref_triad
Cartesian triad of the absolute reference frame.
Definition: bin_ns_bh.h:128
void fait_taij_auto()
Computes (LB)^{ij} auto.
Definition: bin_ns_bh_kij.C:99
Tenseur_sym taij_tot
Total extrinsic curvature tensor $ A^{ij} = 2 N K^{ij}$ generated by { shift_auto} and { shift_comp}...
Definition: et_bin_nsbh.h:140
Cmp & set()
Read/write for a scalar (see also operator=(const Cmp&) ).
Definition: tenseur.C:840
void change_triad(const Base_vect &new_triad)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
Definition: tenseur.C:684
Tbl norme(const Cmp &)
Sums of the absolute values of all the values of the Cmp in each domain.
Definition: cmp_math.C:484
Et_bin_nsbh star
The neutron star.
Definition: bin_ns_bh.h:131
const Base_vect * get_triad() const
Returns the vectorial basis (triad) on which the components are defined.
Definition: tenseur.h:707
Map_af & mp
Affine mapping.
Definition: bhole.h:273
Cmp decouple
Function used to construct the part of generated by the star from the total .
Definition: etoile.h:1003
const double * get_beta() const
Returns the pointer on the array beta.
Definition: map_af.C:611
Tenseur akcar_comp
Part of the scalar generated by shift_auto and shift_comp , i.e.
Definition: etoile.h:947
void fait_decouple()
Function used to compute the { decouple} functions for both the NS and the BH.
Map & mp
Mapping associated with the star.
Definition: etoile.h:432
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:465
Tbl max(const Cmp &)
Maximum values of a Cmp in each domain.
Definition: cmp_math.C:438
const Tenseur & get_nnn() const
Returns the total lapse function N.
Definition: etoile.h:730
int get_etat() const
Returns the logical state.
Definition: tenseur.h:710
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:183
Tenseur_sym tkij_auto
Auto .
Definition: bhole.h:307
Cmp pow(const Cmp &, int)
Power .
Definition: cmp_math.C:351
Tenseur contract(const Tenseur &, int id1, int id2)
Self contraction of two indices of a Tenseur .
void import_symy(const Cmp &ci)
Assignment to another Cmp defined on a different mapping.
void std_base_scal()
Sets the spectral bases of the Valeur va to the standard ones for a scalar.
Definition: cmp.C:647
Tenseur_sym tkij_tot
Total .
Definition: bhole.h:308
int get_nr(int l) const
Returns the number of points in the radial direction ( ) in domain no. l.
Definition: grilles.h:469
Tenseur a_car
Total conformal factor .
Definition: etoile.h:518
Coord ya
Absolute y coordinate.
Definition: map.h:749
void allocate_all()
Sets the logical state to ETATQCQ (ordinary state) and performs the memory allocation of all the elem...
Definition: cmp.C:326
Tbl & set(int l)
Read/write of the value in a given domain.
Definition: cmp.h:724
const Base_vect_cart & get_bvect_cart() const
Returns the Cartesian basis associated with the coordinates (x,y,z) of the mapping, i.e.
Definition: map.h:809
Tenseur_sym taij_auto
Part of the extrinsic curvature tensor $ A^{ij} = 2 N K^{ij}$ generated by { shift_auto}.
Definition: et_bin_nsbh.h:126
Cmp decouple
Function used to construct the part of generated by the hole from the total .
Definition: bhole.h:318
void raccord(int n)
Performs the matching of the nucleus with respect to the first shell.
Definition: cmp_raccord.C:173
Coord za
Absolute z coordinate.
Definition: map.h:750
void fait_tkij(int bound_nn=-1, double lim_nn=0)
Computation of the extrinsic curvature tensor for both { star} and { bhole}.
const Tenseur & get_shift_auto() const
Returns the part of generated by the hole.
Definition: bhole.h:440
Tenseur_sym tkij_tot
Total extrinsic curvature tensor $K^{ij}$ generated by { shift_auto} and { shift_comp}.
Definition: et_bin_nsbh.h:152
Bhole hole
The black hole.
Definition: bin_ns_bh.h:134
Cmp sin(const Cmp &)
Sine.
Definition: cmp_math.C:72
Tenseur n_auto
Part of N generated by the hole.
Definition: bhole.h:286
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
Definition: grilles.h:474
Tenseur_sym tkij_auto
Part of the extrinsic curvature tensor $K^{ij}$ generated by { shift_auto}.
Definition: et_bin_nsbh.h:146
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: tenseur.C:652
void import(const Cmp &ci)
Assignment to another Cmp defined on a different mapping.
Definition: cmp_import.C:76
const Tenseur & get_shift_auto() const
Returns the part of the shift vector generated principaly by the star.
Definition: etoile.h:1155
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:305
Tenseur akcar_auto
Part of the scalar generated by shift_auto , i.e.
Definition: etoile.h:941
void import_asymy(const Cmp &ci)
Assignment to another Cmp defined on a different mapping.
void set_etat_zero()
Sets the logical state to ETATZERO (zero state).
Definition: tenseur.C:661
Tenseur_sym tkij_comp
Part of the extrinsic curvature tensor generated by shift_comp .
Definition: etoile.h:935
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:735
Tensor handling *** DEPRECATED : use class Tensor instead ***.
Definition: tenseur.h:304
const Tenseur & gradient() const
Returns the gradient of *this (Cartesian coordinates)
Definition: tenseur.C:1558