LORENE
bhole.C
1 /*
2  * Methods of class Bhole
3  *
4  * (see file bhole.h for documentation)
5  *
6  */
7 
8 /*
9  * Copyright (c) 2000-2001 Philippe Grandclement
10  *
11  * This file is part of LORENE.
12  *
13  * LORENE is free software; you can redistribute it and/or modify
14  * it under the terms of the GNU General Public License as published by
15  * the Free Software Foundation; either version 2 of the License, or
16  * (at your option) any later version.
17  *
18  * LORENE is distributed in the hope that it will be useful,
19  * but WITHOUT ANY WARRANTY; without even the implied warranty of
20  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
21  * GNU General Public License for more details.
22  *
23  * You should have received a copy of the GNU General Public License
24  * along with LORENE; if not, write to the Free Software
25  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
26  *
27  */
28 
29 
30 
31 
32 /*
33  * $Id: bhole.C,v 1.16 2016/12/05 16:17:45 j_novak Exp $
34  * $Log: bhole.C,v $
35  * Revision 1.16 2016/12/05 16:17:45 j_novak
36  * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
37  *
38  * Revision 1.15 2014/10/13 08:52:39 j_novak
39  * Lorene classes and functions now belong to the namespace Lorene.
40  *
41  * Revision 1.14 2014/10/06 15:12:57 j_novak
42  * Modified #include directives to use c++ syntax.
43  *
44  * Revision 1.13 2007/04/26 13:16:21 f_limousin
45  * Correction of an error in the computation of grad_n_tot and grad_psi_tot
46  *
47  * Revision 1.12 2007/04/24 20:14:04 f_limousin
48  * Implementation of Dirichlet and Neumann BC for the lapse
49  *
50  * Revision 1.11 2006/09/25 10:01:48 p_grandclement
51  * Addition of N-dimensional Tbl
52  *
53  * Revision 1.10 2006/06/01 12:47:51 p_grandclement
54  * update of the Bin_ns_bh project
55  *
56  * Revision 1.9 2006/04/27 09:12:31 p_grandclement
57  * First try at irrotational black holes
58  *
59  * Revision 1.8 2005/08/29 15:10:13 p_grandclement
60  * Addition of things needed :
61  * 1) For BBH with different masses
62  * 2) Provisory files for the mixted binaries (Bh and NS) : THIS IS NOT
63  * WORKING YET !!!
64  *
65  * Revision 1.7 2003/12/16 05:27:07 k_taniguchi
66  * Change some arguments.
67  *
68  * Revision 1.6 2003/11/25 07:10:05 k_taniguchi
69  * Change some arguments from the class Etoile_bin to Et_bin_nsbh.
70  *
71  * Revision 1.5 2003/02/13 16:40:25 p_grandclement
72  * Addition of various things for the Bin_ns_bh project, non of them being
73  * completely tested
74  *
75  * Revision 1.4 2003/01/31 16:57:12 p_grandclement
76  * addition of the member Cmp decouple used to compute the K_ij auto, once
77  * the K_ij total is known
78  *
79  * Revision 1.3 2002/10/16 14:36:31 j_novak
80  * Reorganization of #include instructions of standard C++, in order to
81  * use experimental version 3 of gcc.
82  *
83  * Revision 1.2 2001/12/04 21:27:52 e_gourgoulhon
84  *
85  * All writing/reading to a binary file are now performed according to
86  * the big endian convention, whatever the system is big endian or
87  * small endian, thanks to the functions fwrite_be and fread_be
88  *
89  * Revision 1.1.1.1 2001/11/20 15:19:27 e_gourgoulhon
90  * LORENE
91  *
92  * Revision 2.23 2001/06/28 15:08:27 eric
93  * Ajout de la fonction area().
94  *
95  * Revision 2.22 2001/05/16 14:22:14 phil
96  * modif init_bhole_seul
97  *
98  * Revision 2.21 2001/05/16 11:33:08 phil
99  * ajout de init_bhole_seul
100  *
101  * Revision 2.20 2001/05/07 12:45:29 phil
102  * *** empty log message ***
103  *
104  * Revision 2.19 2001/05/07 09:30:34 phil
105  * *** empty log message ***
106  *
107  * Revision 2.18 2001/05/07 09:11:36 phil
108  * *** empty log message ***
109  *
110  * Revision 2.17 2001/04/26 12:04:27 phil
111  * *** empty log message ***
112  *
113  * Revision 2.16 2001/04/02 12:17:21 phil
114  * *** empty log message ***
115  *
116  * Revision 2.15 2001/02/28 13:22:43 phil
117  * vire kk_auto
118  *
119  * Revision 2.14 2001/01/29 14:30:37 phil
120  * ajout type_rotation
121  *
122  * Revision 2.13 2001/01/10 09:31:31 phil
123  * vire fait_kk_auto
124  *
125  * Revision 2.12 2001/01/09 13:38:26 phil
126  * ajout memebre kk_auto
127  *
128  * Revision 2.11 2000/12/21 10:48:59 phil
129  * retour a version 2.9
130  *
131  * Revision 2.9 2000/12/18 16:40:45 phil
132  * *** empty log message ***
133  *
134  * Revision 2.8 2000/12/14 10:44:48 phil
135  * ATTENTION : PASSAGE DE PHI A PSI
136  *
137  * Revision 2.7 2000/12/04 14:29:47 phil
138  * ajout de grad_n_tot
139  *
140  * Revision 2.6 2000/12/01 16:39:25 phil
141  * correction impoetation grad_phi_comp
142  *
143  * Revision 2.5 2000/12/01 16:21:25 phil
144  * modification init_bhole
145  *
146  * Revision 2.4 2000/12/01 16:17:42 phil
147  * correction base import de grad_phi_comp
148  *
149  * Revision 2.3 2000/12/01 16:14:33 phil
150  * *** empty log message ***
151  *
152  * Revision 2.2 2000/12/01 16:12:55 phil
153  * *** empty log message ***
154  *
155  * Revision 2.1 2000/10/23 09:22:07 phil
156  * rearrangement dans constructeurs.
157  *
158  * Revision 2.0 2000/10/20 09:18:47 phil
159  * *** empty log message ***
160  *
161  *
162  * $Header: /cvsroot/Lorene/C++/Source/Bhole/bhole.C,v 1.16 2016/12/05 16:17:45 j_novak Exp $
163  *
164  */
165 
166 //standard
167 #include <cstdlib>
168 #include <cmath>
169 
170 // Lorene
171 #include "tenseur.h"
172 #include "bhole.h"
173 #include "proto.h"
174 #include "utilitaires.h"
175 #include "et_bin_nsbh.h"
176 #include "graphique.h"
177 
178 // Constructeur standard
179 namespace Lorene {
180 Bhole::Bhole (Map_af& mpi) : mp(mpi),
181  rayon ((mp.get_alpha())[0]),
182  omega(0), omega_local(0), rot_state(COROT), boost (new double[3]), regul(0),
183  n_auto(mpi), n_comp(mpi), n_tot(mpi),
184  psi_auto(mpi), psi_comp(mpi), psi_tot(mpi),
185  grad_n_tot (mpi, 1, COV, mpi.get_bvect_cart()),
186  grad_psi_tot (mpi, 1, COV, mpi.get_bvect_cart()),
187  shift_auto(mpi, 1, CON, mpi.get_bvect_cart()),
188  taij_auto(mpi, 2, CON, mpi.get_bvect_cart()),
189  taij_comp(mpi, 2, CON, mpi.get_bvect_cart()),
190  taij_tot (mpi, 2, CON, mpi.get_bvect_cart()),
191  tkij_auto (mpi, 2, CON, mpi.get_bvect_cart()),
192  tkij_tot (mpi, 2, CON, mpi.get_bvect_cart()),
193  decouple (mpi) {
194  for (int i=0 ; i<3 ; i++)
195  boost[i] = 0 ;
196 }
197 
198 
199 // Constructeur par recopie
200 Bhole::Bhole (const Bhole& source) :
201  mp(source.mp), rayon(source.rayon), omega (source.omega), omega_local(source.omega_local), rot_state(source.rot_state),
202  boost (new double [3]), regul(source.regul), n_auto(source.n_auto),
203  n_comp(source.n_comp), n_tot(source.n_tot), psi_auto(source.psi_auto),
204  psi_comp(source.psi_comp), psi_tot(source.psi_tot),
205  grad_n_tot (source.grad_n_tot), grad_psi_tot (source.grad_psi_tot),
206  shift_auto(source.shift_auto),
207  taij_auto (source.taij_auto),
208  taij_comp(source.taij_comp), taij_tot(source.taij_tot),
209  tkij_auto(source.tkij_auto), tkij_tot(source.tkij_tot), decouple(source.decouple) {
210 
211  for (int i=0 ; i<3 ; i++)
212  boost[i] = source.boost[i] ;
213 }
214 
215 // Destructeur
217  delete [] boost ;
218 }
219 
220 // Affectation
221 void Bhole::operator= (const Bhole& source) {
222 
223  assert (&mp == &source.mp) ;
224 
225  rayon = source.rayon ;
226  omega = source.omega ;
227  omega_local = source.omega_local ;
228  rot_state = source.rot_state ;
229  for (int i=0 ; i<3 ; i++)
230  boost[i] = source.boost[i] ;
231  regul = regul ;
232 
233  n_auto = source.n_auto ;
234  n_comp = source.n_comp ;
235  n_tot = source.n_tot ;
236 
237  psi_auto = source.psi_auto ;
238  psi_comp = source.psi_comp ;
239  psi_tot = source.psi_tot ;
240 
241  grad_n_tot = source.grad_n_tot ;
242  grad_psi_tot = source.grad_psi_tot ;
243 
244  shift_auto = source.shift_auto ;
245 
246  taij_auto = source.taij_auto ;
247  taij_comp = source.taij_comp ;
248  taij_tot = source.taij_tot ;
249 
250  tkij_auto = source.tkij_auto ;
251  tkij_tot = source.tkij_tot ;
252  decouple = source.decouple ;
253 }
254 
255 // Importe le lapse du compagnon (Bhole case)
256 
257 void Bhole::fait_n_comp (const Bhole& comp) {
258 
259  n_comp.set_etat_qcq() ;
260  n_comp.set().import_symy(comp.n_auto()) ;
261  n_comp.set_std_base() ;
262 
263  n_tot = n_comp + n_auto ;
264  n_tot.set_std_base() ;
265 
266  Tenseur auxi (comp.n_auto.gradient()) ;
267  auxi.dec2_dzpuis() ;
268 
269  Tenseur grad_comp (mp, 1, COV, *auxi.get_triad()) ;
270  grad_comp.set_etat_qcq() ;
271  grad_comp.set(0).import_symy(auxi(0)) ;
272  grad_comp.set(1).import_asymy(auxi(1)) ;
273  grad_comp.set(2).import_symy(auxi(2)) ;
274  grad_comp.set_std_base() ;
275  grad_comp.inc2_dzpuis() ;
276  grad_comp.change_triad(mp.get_bvect_cart()) ;
277 
278  grad_n_tot = n_auto.gradient() + grad_comp ;
279 }
280 
281 // Importe le facteur conforme du compagnon (Bhole case)
282 
283 void Bhole::fait_psi_comp (const Bhole& comp) {
284 
286  psi_comp.set().import_symy(comp.psi_auto()) ;
288 
291 
292 
293  Tenseur auxi (comp.psi_auto.gradient()) ;
294  auxi.dec2_dzpuis() ;
295 
296  Tenseur grad_comp (mp, 1, COV, *auxi.get_triad()) ;
297  grad_comp.set_etat_qcq() ;
298  grad_comp.set(0).import_symy(auxi(0)) ;
299  grad_comp.set(1).import_asymy(auxi(1)) ;
300  grad_comp.set(2).import_symy(auxi(2)) ;
301  grad_comp.set_std_base() ;
302  grad_comp.inc2_dzpuis() ;
303  grad_comp.change_triad(mp.get_bvect_cart()) ;
304 
305  grad_psi_tot = psi_auto.gradient() + grad_comp ;
306 }
307 
308 
309 // Importe le lapse du compagnon (NS case)
310 void Bhole::fait_n_comp (const Et_bin_nsbh& comp) {
311 
312  n_comp.set_etat_qcq() ;
313  if (regul == 0)
314  n_comp.set().import(comp.get_n_auto()()) ;
315  else
316  n_comp.set().import_symy(comp.get_n_auto()()) ;
317 
318  n_comp.set_std_base() ;
319 
320  n_tot = n_comp + n_auto ;
321  n_tot.set_std_base() ;
322 
323  Tenseur grad_comp (mp, 1, COV, mp.get_bvect_cart()) ;
324  Tenseur auxi (comp.get_d_n_auto()) ;
325  auxi.dec2_dzpuis() ;
326  grad_comp.set_etat_qcq() ;
327  if (regul == 0){
328  grad_comp.set(0).import(auxi(0)) ;
329  grad_comp.set(1).import(auxi(1)) ;
330  grad_comp.set(2).import(auxi(2)) ;
331  }
332  else{
333  grad_comp.set(0).import_symy(auxi(0)) ;
334  grad_comp.set(1).import_asymy(auxi(1)) ;
335  grad_comp.set(2).import_symy(auxi(2)) ;
336  }
337  grad_comp.set_std_base() ;
338  grad_comp.inc2_dzpuis() ;
339  grad_comp.set_triad( *(auxi.get_triad()) ) ;
340  grad_comp.change_triad(mp.get_bvect_cart()) ;
341 
342  grad_n_tot = n_auto.gradient() + grad_comp ;
343 }
344 
345 // Importe le facteur conforme du compagnon (Bhole case)
346 
347 void Bhole::fait_psi_comp (const Et_bin_nsbh& comp) {
348 
350  if (regul == 0)
351  psi_comp.set().import(comp.get_confpsi_auto()()) ;
352  else
355 
358 
359  Tenseur grad_comp (mp, 1, COV, mp.get_bvect_cart()) ;
360  Tenseur auxi (comp.get_d_confpsi_auto()) ;
361  auxi.dec2_dzpuis() ;
362  grad_comp.set_etat_qcq() ;
363  if (regul == 0){
364  grad_comp.set(0).import(auxi(0)) ;
365  grad_comp.set(1).import(auxi(1)) ;
366  grad_comp.set(2).import(auxi(2)) ;
367  }
368  else{
369  grad_comp.set(0).import_symy(auxi(0)) ;
370  grad_comp.set(1).import_asymy(auxi(1)) ;
371  grad_comp.set(2).import_symy(auxi(2)) ;
372  }
373  grad_comp.set_std_base() ;
374  grad_comp.inc2_dzpuis() ;
375  grad_comp.set_triad( *(auxi.get_triad()) ) ;
376  grad_comp.change_triad(mp.get_bvect_cart()) ;
377 
378  grad_psi_tot = psi_auto.gradient() + grad_comp ;
379 }
380 
381 // Calcul Taij auto (nul sur H que si la regularisation a ete faite)
383 
384  if (shift_auto.get_etat() == ETATZERO)
386  else {
387  Tenseur grad (shift_auto.gradient()) ;
388  Tenseur trace (mp) ;
389  trace = grad(0, 0)+grad(1, 1)+grad(2, 2) ;
390 
393  for (int i=0 ; i<3 ; i++) {
394  for (int j=i+1 ; j<3 ; j++)
395  taij_auto.set(i, j) = grad(i, j)+grad(j, i) ;
396  taij_auto.set(i, i) = 2*grad(i, i) -2./3.*trace() ;
397  }
398 
399  for (int i=0 ; i<3 ; i++)
400  for (int j=0 ; j<3 ; j++)
401  taij_auto.set(i, j).raccord(1) ;
402  }
403 }
404 
405 // Regularise le shift, untilisant le compagnon.
406 void Bhole::regularise_shift (const Tenseur& shift_comp) {
407  regul = regle (shift_auto, shift_comp, omega, omega_local) ;
408 }
409 
410 //Initialise a Schwartz
411 // Compagnon a 0
413 
414  Cmp auxi(mp) ;
415 
416  auxi = 1./2.-rayon/mp.r ;
417  auxi.annule(0);
418  auxi.set_dzpuis(0) ;
419  n_auto = auxi;
420  n_auto.set_std_base() ;
421  n_auto.set().raccord(1) ;
422  n_comp =0.5;
423  n_comp.set_std_base() ;
424  n_tot = n_comp+n_auto ;
425 
426  auxi = 0.5+rayon/mp.r ;
427  auxi.annule(0);
428  auxi.set_dzpuis(0) ;
429  psi_auto = auxi;
431  psi_auto.set().raccord(1) ;
432  psi_comp = 0.5;
435 
438 
446 }
447 
449 
450  Cmp auxi(mp) ;
451 
452  auxi = (1-rayon/mp.r)/(1+rayon/mp.r) ;
453  auxi.annule(0);
454  auxi.set_val_inf(1) ;
455  auxi.set_dzpuis(0) ;
456  n_auto = auxi;
457  n_auto.set_std_base() ;
458  n_auto.set().raccord(1) ;
460  n_tot.set_etat_zero() ;
461 
462  auxi = 1+rayon/mp.r ;
463  auxi.annule(0);
464  auxi.set_val_inf(1) ;
465  auxi.set_dzpuis(0) ;
466  psi_auto = auxi;
468  psi_auto.set().raccord(1) ;
471 
474 
482 }
483 
484 // Sauvegarde dans fichier
485 void Bhole::sauve (FILE* fich) const {
486 
487  fwrite_be (&omega, sizeof(double), 1, fich) ;
488  fwrite_be (&omega_local, sizeof(double), 1, fich) ;
489  fwrite_be (&rot_state, sizeof(int), 1, fich) ;
490  fwrite_be (boost, sizeof(double), 3, fich) ;
491  fwrite_be (&regul, sizeof(double), 1, fich) ;
492 
493  n_auto.sauve(fich) ;
494  psi_auto.sauve(fich) ;
495  shift_auto.sauve(fich) ;
496 }
497 
498 //Constructeur par fichier
499 Bhole::Bhole(Map_af& mpi, FILE* fich, bool old) :
500  mp(mpi), rayon ((mp.get_alpha())[0]),
501  boost (new double[3]), n_auto(mpi), n_comp(mpi), n_tot(mpi),
502  psi_auto(mpi), psi_comp(mpi), psi_tot(mpi),
503  grad_n_tot (mpi, 1, COV, mpi.get_bvect_cart()),
504  grad_psi_tot (mpi, 1, COV, mpi.get_bvect_cart()),
505  shift_auto(mpi, 1, CON, mpi.get_bvect_cart()),
506  taij_auto(mpi, 2, CON, mpi.get_bvect_cart()),
507  taij_comp(mpi, 2, CON, mpi.get_bvect_cart()),
508  taij_tot (mpi, 2, CON, mpi.get_bvect_cart()),
509  tkij_auto (mpi, 2, CON, mpi.get_bvect_cart()) ,
510  tkij_tot (mpi, 2, CON, mpi.get_bvect_cart()) ,
511  decouple (mpi) {
512 
513  fread_be(&omega, sizeof(double), 1, fich) ;
514 
515  if (!old) {
516  fread_be(&omega_local, sizeof(double), 1, fich) ;
517  fread_be(&rot_state, sizeof(int), 1, fich) ;
518  }
519  fread_be(boost, sizeof(double), 3, fich) ;
520  fread_be(&regul, sizeof(double), 1, fich) ;
521 
522  Tenseur n_auto_file (mp, fich) ;
523  n_auto = n_auto_file ;
524 
525  Tenseur psi_auto_file (mp, fich) ;
526  psi_auto = psi_auto_file ;
527 
528  Tenseur shift_auto_file (mp, mp.get_bvect_cart(), fich) ;
529  shift_auto = shift_auto_file ;
530 
533 
537 
541 }
542 
543 
544  //---------------------------------//
545  // Computation of throat area //
546  //---------------------------------//
547 
548 double Bhole::area() const {
549 
550  Cmp integrant( pow(psi_tot(), 4) ) ; // Psi^4
551 
552  integrant.std_base_scal() ;
553  integrant.raccord(1) ;
554 
555  return mp.integrale_surface(integrant, rayon) ;
556 
557 }
558 
559  //---------------------------------//
560  // Moment local //
561  //---------------------------------//
562 
563 double Bhole::local_momentum() const {
564 
565  Cmp xx (mp) ;
566  xx = mp.x ;
567  xx.std_base_scal() ;
568 
569  Cmp yy (mp) ;
570  yy = mp.y ;
571  yy.std_base_scal() ;
572 
573  Tenseur vecteur (mp, 1, CON, mp.get_bvect_cart()) ;
574  vecteur.set_etat_qcq() ;
575  for (int i=0 ; i<3 ; i++)
576  vecteur.set(i) = (-yy*tkij_tot(0, i)+xx*tkij_tot(1, i)) ;
577  vecteur.set_std_base() ;
578  vecteur.annule(mp.get_mg()->get_nzone()-1) ;
579  vecteur.change_triad (mp.get_bvect_spher()) ;
580 
581  Cmp integrant (pow(psi_tot(), 6)*vecteur(0)) ;
582  integrant.std_base_scal() ;
583  double moment = mp.integrale_surface(integrant, rayon)/8/M_PI ;
584 
585  return moment ;
586 }
587 
588 }
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition: cmp.h:446
double integrale_surface(const Cmp &ci, double rayon) const
Performs the surface integration of ci on the sphere of radius rayon .
Tenseur n_tot
Total N .
Definition: bhole.h:288
void dec2_dzpuis()
dzpuis -= 2 ;
Definition: tenseur.C:1146
Tenseur grad_n_tot
Total gradient of N .
Definition: bhole.h:294
void annule(int l)
Sets the Cmp to zero in a given domain.
Definition: cmp.C:351
void set_std_base()
Set the standard spectal basis of decomposition for each component.
Definition: tenseur.C:1186
const Base_vect_spher & get_bvect_spher() const
Returns the orthonormal vectorial basis associated with the coordinates of the mapping.
Definition: map.h:801
~Bhole()
Destructor.
Definition: bhole.C:216
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
Bhole(Map_af &mapping)
Standard constructor.
Definition: bhole.C:180
void init_bhole_seul()
Initiates for a single the black hole.
Definition: bhole.C:448
double omega_local
local angular velocity
Definition: bhole.h:276
Tenseur psi_auto
Part of generated by the hole.
Definition: bhole.h:290
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
void fait_n_comp(const Bhole &comp)
Imports the part of N due to the companion hole comp .
Definition: bhole.C:257
Black hole.
Definition: bhole.h:268
Tenseur shift_auto
Part of generated by the hole.
Definition: bhole.h:297
Tenseur_sym taij_comp
Part of generated by the companion hole.
Definition: bhole.h:300
void set_etat_zero()
Sets the logical state to ETATZERO (zero).
Definition: cmp.C:292
Cmp & set()
Read/write for a scalar (see also operator=(const Cmp&) ).
Definition: tenseur.C:840
void sauve(FILE *) const
Save in a file.
Definition: tenseur.C:1341
void regularise_shift(const Tenseur &comp)
Corrects shift_auto in such a way that the total is equal to zero in the horizon, which should ensure the regularity of , using the stored values of the boost and the angular velocity.
Definition: bhole.C:406
Tenseur psi_comp
Part of generated by the companion hole.
Definition: bhole.h:291
double area() const
Computes the area of the throat.
Definition: bhole.C:548
Map_af & mp
Affine mapping.
Definition: bhole.h:273
const Tenseur & get_n_auto() const
Returns the part of the lapse { N} generated principaly by the star.
Definition: et_bin_nsbh.h:239
void set_val_inf(double val)
Sets the value of the Cmp to val at infinity.
Definition: cmp_manip.C:129
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:465
int get_etat() const
Returns the logical state.
Definition: tenseur.h:710
const Tenseur & get_d_confpsi_auto() const
Returns the gradient of { confpsi_auto} (Cartesian components with respect to { ref_triad}) ...
Definition: et_bin_nsbh.h:272
int fwrite_be(const int *aa, int size, int nb, FILE *fich)
Writes integer(s) into a binary file according to the big endian convention.
Definition: fwrite_be.C:73
void sauve(FILE *fich) const
Write on a file.
Definition: bhole.C:485
double rayon
Radius of the horizon in LORENE&#39;s units.
Definition: bhole.h:274
Cmp pow(const Cmp &, int)
Power .
Definition: cmp_math.C:351
Tenseur_sym tkij_auto
Auto .
Definition: bhole.h:307
void import_symy(const Cmp &ci)
Assignment to another Cmp defined on a different mapping.
double * boost
Lapse on the horizon.
Definition: bhole.h:283
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 fread_be(int *aa, int size, int nb, FILE *fich)
Reads integer(s) from a binary file according to the big endian convention.
Definition: fread_be.C:72
Tenseur grad_psi_tot
Total gradient of .
Definition: bhole.h:295
void init_bhole()
Sets the values of the fields to :
Definition: bhole.C:412
void operator=(const Bhole &)
Affectation.
Definition: bhole.C:221
Affine radial mapping.
Definition: map.h:2048
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
Class for a star in a NS-BH binary system.
Definition: et_bin_nsbh.h:79
const Tenseur & get_confpsi_auto() const
Returns the part of the conformal factor $$ generated principaly by the star.
Definition: et_bin_nsbh.h:262
int rot_state
State of rotation.
Definition: bhole.h:277
Coord y
y coordinate centered on the grid
Definition: map.h:745
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
Tenseur psi_tot
Total .
Definition: bhole.h:292
Coord x
x coordinate centered on the grid
Definition: map.h:744
double omega
Angular velocity in LORENE&#39;s units.
Definition: bhole.h:275
void set_dzpuis(int)
Set a value to dzpuis.
Definition: cmp.C:657
Tenseur n_auto
Part of N generated by the hole.
Definition: bhole.h:286
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: tenseur.C:652
void fait_psi_comp(const Bhole &comp)
Imports the part of due to the companion hole comp .
Definition: bhole.C:283
void import(const Cmp &ci)
Assignment to another Cmp defined on a different mapping.
Definition: cmp_import.C:76
void set_etat_zero()
Sets the logical state to ETATZERO (zero state).
Definition: tenseur.C:661
Tensor handling *** DEPRECATED : use class Tensor instead ***.
Definition: tenseur.h:304
double regul
Intensity of the correction on the shift vector.
Definition: bhole.h:284
const Tenseur & gradient() const
Returns the gradient of *this (Cartesian coordinates)
Definition: tenseur.C:1558
Coord r
r coordinate centered on the grid
Definition: map.h:736
const Tenseur & get_d_n_auto() const
Returns the gradient of { n_auto} (Cartesian components with respect to { ref_triad}) ...
Definition: et_bin_nsbh.h:249