LORENE
hole_bhns_bc.C
1 /*
2  * Methods of class Hole_bhns to compute the inner boundary condition
3  * at the excised surface
4  *
5  * (see file hile_bhns.h for documentation).
6  *
7  */
8 
9 /*
10  * Copyright (c) 2005-2007 Keisuke Taniguchi
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 version 2
16  * as published by the Free Software Foundation.
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  * $Id: hole_bhns_bc.C,v 1.5 2016/12/05 16:17:55 j_novak Exp $
33  * $Log: hole_bhns_bc.C,v $
34  * Revision 1.5 2016/12/05 16:17:55 j_novak
35  * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
36  *
37  * Revision 1.4 2014/10/13 08:53:00 j_novak
38  * Lorene classes and functions now belong to the namespace Lorene.
39  *
40  * Revision 1.3 2014/10/06 15:13:10 j_novak
41  * Modified #include directives to use c++ syntax.
42  *
43  * Revision 1.2 2008/05/15 19:04:10 k_taniguchi
44  * Change of some parameters.
45  *
46  * Revision 1.1 2007/06/22 01:23:56 k_taniguchi
47  * *** empty log message ***
48  *
49  *
50  * $Header: /cvsroot/Lorene/C++/Source/Hole_bhns/hole_bhns_bc.C,v 1.5 2016/12/05 16:17:55 j_novak Exp $
51  *
52  */
53 
54 // C++ headers
55 //#include <>
56 
57 // C headers
58 #include <cmath>
59 
60 // Lorene headers
61 #include "hole_bhns.h"
62 #include "valeur.h"
63 #include "grilles.h"
64 #include "unites.h"
65 
66  //----------------------------------//
67  // Inner boundary condition //
68  //----------------------------------//
69 
70 namespace Lorene {
72 
73  // Fundamental constants and units
74  // -------------------------------
75  using namespace Unites ;
76 
77  const Mg3d* mg = mp.get_mg() ;
78  const Mg3d* mg_angu = mg->get_angu() ;
79  Valeur bc(mg_angu) ;
80 
81  int nt = mg->get_nt(0) ;
82  int np = mg->get_np(0) ;
83 
84  Scalar tmp(mp) ;
85 
86  // double cc ; // C/M^2
87 
88  if (bc_lapconf_nd) {
89 
90  Scalar st(mp) ;
91  st = mp.sint ;
92  st.std_spectral_base() ;
93  Scalar ct(mp) ;
94  ct = mp.cost ;
95  ct.std_spectral_base() ;
96  Scalar sp(mp) ;
97  sp = mp.sinp ;
98  sp.std_spectral_base() ;
99  Scalar cp(mp) ;
100  cp = mp.cosp ;
101  cp.std_spectral_base() ;
102 
103  Scalar rr(mp) ;
104  rr = mp.r ;
105  rr.std_spectral_base() ;
106 
107  if (bc_lapconf_fs) { // dlapconf/dr = 0
108 
109  if (kerrschild) {
110  cout << "!!!!! WARNING: Not yet available !!!!!" << endl ;
111  abort() ;
112  }
113  else { // Isotropic coordinates with the maximal slicing
114  tmp = - d_lapconf_comp(1) % st % cp
115  - d_lapconf_comp(2) % st % sp - d_lapconf_comp(3) % ct ;
116  }
117 
118  }
119  else { // dlapconf/dr = 0.5*lapconf/rr
120 
121  Scalar tmp1(mp) ;
122  tmp1 = 0.5 * (lapconf_auto_rs + lapconf_comp) / rr ;
123  tmp1.std_spectral_base() ;
124  tmp1.inc_dzpuis(2) ; // dzpuis : 0 -> 2
125 
126  if (kerrschild) { // dlapconf/dr = 0.5*lapconf/rr
127  cout << "!!!!! WARNING: Not yet available !!!!!" << endl ;
128  abort() ;
129  }
130  else { // Isotropic coordinates with the maximal slicing
131  // dlapconf/dr = 0.5*lapconf/rr
132 
133  tmp = - d_lapconf_comp(1) % st % cp
134  - d_lapconf_comp(2) % st % sp - d_lapconf_comp(3) % ct
135  + tmp1 ;
136  }
137 
138  }
139  }
140  else {
141 
142  if (bc_lapconf_fs) { // The poisson solver in LORENE assumes
143  // the asymptotic behavior of the function -> 0
144 
145  if (kerrschild) {
146  cout << "!!!!! WARNING: Not yet available !!!!!" << endl ;
147  abort() ;
148  // lapconf_auto -> 0.5 <-> lapconf_auto_rs -> -0.5
149  }
150  else { // Isotropic coordinates with the maximal slicing
151  cout << "!!!!! WARNING: Not yet prepared !!!!!" << endl ;
152  abort() ;
153  // tmp = -lapconf_comp + 0.5 ; // lapconf = 0.5
154 
155  }
156 
157  }
158  else {
159 
160  if (kerrschild) {
161  cout << "!!!!! WARNING: Not yet available !!!!!" << endl ;
162  abort() ;
163  }
164  else { // Isotropic coordinates with the maximal slicing
165  cout << "!!!!! WARNING: Not yet prepared !!!!!" << endl ;
166  abort() ;
167  // tmp = -lapconf_comp + 0.5 ;
168 
169  }
170 
171  }
172  }
173 
174  bc = 1. ;
175  for (int j=0; j<nt; j++) {
176  for (int k=0; k<np; k++) {
177  bc.set(0,k,j,0) = tmp.val_grid_point(1,k,j,0) ;
178  }
179  }
180 
181  bc.std_base_scal() ;
182  return bc ;
183 
184 }
185 
186 const Valeur Hole_bhns::bc_shift_x(double ome_orb, double y_rot) const {
187 
188  // Fundamental constants and units
189  // -------------------------------
190  using namespace Unites ;
191 
192  const Mg3d* mg = mp.get_mg() ;
193  const Mg3d* mg_angu = mg->get_angu() ;
194  Valeur bc(mg_angu) ;
195 
196  Base_val** bases = mp.get_mg()->std_base_vect_cart() ;
197 
198  Scalar rr(mp) ;
199  rr = mp.r ;
200  rr.std_spectral_base() ;
201  Scalar st(mp) ;
202  st = mp.sint ;
203  st.std_spectral_base() ;
204  Scalar cp(mp) ;
205  cp = mp.cosp ;
206  cp.std_spectral_base() ;
207  Scalar yy(mp) ;
208  yy = mp.y ;
209  yy.std_spectral_base() ;
210 
211  double mass = ggrav * mass_bh ;
212  double ori_y_bh = mp.get_ori_y() ;
213 
214  Scalar tmp(mp) ;
215 
216  if (kerrschild) {
217 
218  // Exact solution of an isolated BH is extracted
219 
220  cout << "!!!!! WARNING: Not yet available !!!!!" << endl ;
221  abort() ;
222 
223  }
224  else { // Isotropic coordinates with the maximal slicing
225 
226  double cc ;
227 
228  // Sets C/M^2 for each case of the lapse boundary condition
229  // --------------------------------------------------------
230 
231  if (bc_lapconf_nd) { // Neumann boundary condition
232  if (bc_lapconf_fs) { // First condition
233  // d(\alpha \psi)/dr = 0
234  // ---------------------
235  cc = 2. * (sqrt(13.) - 1.) / 3. ;
236  }
237  else { // Second condition
238  // d(\alpha \psi)/dr = (\alpha \psi)/(2 rah)
239  // -----------------------------------------
240  cc = 4. / 3. ;
241  }
242  }
243  else { // Dirichlet boundary condition
244  if (bc_lapconf_fs) { // First condition
245  // (\alpha \psi) = 1/2
246  // -------------------
247  cout << "!!!!! WARNING: Not yet prepared !!!!!" << endl ;
248  abort() ;
249  }
250  else { // Second condition
251  // (\alpha \psi) = 1/sqrt(2.) \psi_KS
252  // ----------------------------------
253  cout << "!!!!! WARNING: Not yet prepared !!!!!" << endl ;
254  abort() ;
255  // cc = 2. * sqrt(2.) ;
256  }
257  }
258 
259  Scalar r_are(mp) ;
261  r_are.std_spectral_base() ;
262 
263  /*
264  tmp = ((lapse_tot / confo_tot / confo_tot)
265  - mass*mass*cc/rr/rr/pow(r_are,3.)) * st * cp
266  - shift_comp(1)
267  + (ome_orb - omega_spin) * yy + ome_orb * (ori_y_bh - y_rot) ;
268  */
269  tmp = ((lapconf_tot / pow(confo_tot,3.)) - (0.25*cc/r_are)) * st * cp
270  - shift_comp(1)
271  + (ome_orb - omega_spin) * yy + ome_orb * (ori_y_bh - y_rot) ;
272 
273  }
274 
275  int nt = mg->get_nt(0) ;
276  int np = mg->get_np(0) ;
277 
278  bc = 1. ;
279  for (int j=0; j<nt; j++) {
280  for (int k=0; k<np; k++) {
281  bc.set(0,k,j,0) = tmp.val_grid_point(1,k,j,0) ;
282  }
283  }
284 
285  bc.base = *bases[0] ;
286 
287  for (int i=0; i<3; i++)
288  delete bases[i] ;
289 
290  delete [] bases ;
291 
292  return bc ;
293 
294 }
295 
296 const Valeur Hole_bhns::bc_shift_y(double ome_orb, double x_rot) const {
297 
298  // Fundamental constants and units
299  // -------------------------------
300  using namespace Unites ;
301 
302  const Mg3d* mg = mp.get_mg() ;
303  const Mg3d* mg_angu = mg->get_angu() ;
304  Valeur bc(mg_angu) ;
305 
306  Base_val** bases = mp.get_mg()->std_base_vect_cart() ;
307 
308  Scalar rr(mp) ;
309  rr = mp.r ;
310  rr.std_spectral_base() ;
311  Scalar st(mp) ;
312  st = mp.sint ;
313  st.std_spectral_base() ;
314  Scalar sp(mp) ;
315  sp = mp.sinp ;
316  sp.std_spectral_base() ;
317  Scalar xx(mp) ;
318  xx = mp.x ;
319  xx.std_spectral_base() ;
320 
321  double mass = ggrav * mass_bh ;
322  double ori_x_bh = mp.get_ori_x() ;
323 
324  Scalar tmp(mp) ;
325 
326  if (kerrschild) {
327 
328  // Exact solution of an isolated BH is extracted
329 
330  cout << "!!!!! WARNING: Not yet available !!!!!" << endl ;
331  abort() ;
332 
333  }
334  else { // Isotropic coordinates with the maximal slicing
335 
336  double cc ;
337 
338  // Sets C/M^2 for each case of the lapse boundary condition
339  // --------------------------------------------------------
340 
341  if (bc_lapconf_nd) { // Neumann boundary condition
342  if (bc_lapconf_fs) { // First condition
343  // d(\alpha \psi)/dr = 0
344  // ---------------------
345  cc = 2. * (sqrt(13.) - 1.) / 3. ;
346  }
347  else { // Second condition
348  // d(\alpha \psi)/dr = (\alpha \psi)/(2 rah)
349  // -----------------------------------------
350  cc = 4. / 3. ;
351  }
352  }
353  else { // Dirichlet boundary condition
354  if (bc_lapconf_fs) { // First condition
355  // (\alpha \psi) = 1/2
356  // -------------------
357  cout << "!!!!! WARNING: Not yet prepared !!!!!" << endl ;
358  abort() ;
359  }
360  else { // Second condition
361  // (\alpha \psi) = 1/sqrt(2.) \psi_KS
362  // ----------------------------------
363  cout << "!!!!! WARNING: Not yet prepared !!!!!" << endl ;
364  abort() ;
365  // cc = 2. * sqrt(2.) ;
366  }
367  }
368 
369  Scalar r_are(mp) ;
371  r_are.std_spectral_base() ;
372 
373  /*
374  tmp = ((lapse_tot / confo_tot / confo_tot)
375  - mass*mass*cc/rr/rr/pow(r_are,3.)) * st * sp
376  - shift_comp(2)
377  - (ome_orb - omega_spin) * xx - ome_orb * (ori_x_bh - x_rot) ;
378  */
379  tmp = ((lapconf_tot / pow(confo_tot,3.)) - (0.25*cc/r_are)) * st * sp
380  - shift_comp(2)
381  - (ome_orb - omega_spin) * xx - ome_orb * (ori_x_bh - x_rot) ;
382 
383  }
384 
385  int nt = mg->get_nt(0) ;
386  int np = mg->get_np(0) ;
387 
388  bc = 1. ;
389  for (int j=0; j<nt; j++) {
390  for (int k=0; k<np; k++) {
391  bc.set(0,k,j,0) = tmp.val_grid_point(1,k,j,0) ;
392  }
393  }
394 
395  bc.base = *bases[1] ;
396 
397  for (int i=0; i<3; i++)
398  delete bases[i] ;
399 
400  delete [] bases ;
401 
402  return bc ;
403 
404 }
405 
407 
408  // Fundamental constants and units
409  // -------------------------------
410  using namespace Unites ;
411 
412  const Mg3d* mg = mp.get_mg() ;
413  const Mg3d* mg_angu = mg->get_angu() ;
414  Valeur bc(mg_angu) ;
415 
416  Base_val** bases = mp.get_mg()->std_base_vect_cart() ;
417 
418  Scalar rr(mp) ;
419  rr = mp.r ;
420  rr.std_spectral_base() ;
421  Scalar ct(mp) ;
422  ct = mp.cost ;
423  ct.std_spectral_base() ;
424 
425  double mass = ggrav * mass_bh ;
426 
427  Scalar tmp(mp) ;
428 
429  if (kerrschild) {
430 
431  // Exact solution of an isolated BH is extracted
432 
433  cout << "!!!!! WARNING: Not yet available !!!!!" << endl ;
434  abort() ;
435 
436  }
437  else { // Isotropic coordinates with the maximal slicing
438 
439  double cc ;
440 
441  // Sets C/M^2 for each case of the lapse boundary condition
442  // --------------------------------------------------------
443 
444  if (bc_lapconf_nd) { // Neumann boundary condition
445  if (bc_lapconf_fs) { // First condition
446  // d(\alpha \psi)/dr = 0
447  // ---------------------
448  cc = 2. * (sqrt(13.) - 1.) / 3. ;
449  }
450  else { // Second condition
451  // d(\alpha \psi)/dr = (\alpha \psi)/(2 rah)
452  // -----------------------------------------
453  cc = 4. / 3. ;
454  }
455  }
456  else { // Dirichlet boundary condition
457  if (bc_lapconf_fs) { // First condition
458  // (\alpha \psi) = 1/2
459  // -------------------
460  cout << "!!!!! WARNING: Not yet prepared !!!!!" << endl ;
461  abort() ;
462  }
463  else { // Second condition
464  // (\alpha \psi) = 1/sqrt(2.) \psi_KS
465  // ----------------------------------
466  cout << "!!!!! WARNING: Not yet prepared !!!!!" << endl ;
467  abort() ;
468  // cc = 2. * sqrt(2.) ;
469  }
470  }
471 
472  Scalar r_are(mp) ;
474  r_are.std_spectral_base() ;
475 
476  /*
477  tmp = ((lapse_tot / confo_tot / confo_tot)
478  - mass*mass*cc/rr/rr/pow(r_are,3.)) * ct - shift_comp(3) ;
479  */
480  tmp = ((lapconf_tot / pow(confo_tot,3.)) - (0.25*cc/r_are)) * ct
481  - shift_comp(3) ;
482 
483  }
484 
485  int nt = mg->get_nt(0) ;
486  int np = mg->get_np(0) ;
487 
488  bc = 1. ;
489  for (int j=0; j<nt; j++) {
490  for (int k=0; k<np; k++) {
491  bc.set(0,k,j,0) = tmp.val_grid_point(1,k,j,0) ;
492  }
493  }
494 
495  bc.base = *bases[2] ;
496 
497  for (int i=0; i<3; i++)
498  delete bases[i] ;
499 
500  delete [] bases ;
501 
502  return bc ;
503 
504 }
505 
506 const Valeur Hole_bhns::bc_confo(double ome_orb, double x_rot,
507  double y_rot) const {
508 
509  // Fundamental constants and units
510  // -------------------------------
511  using namespace Unites ;
512 
513  const Mg3d* mg = mp.get_mg() ;
514  const Mg3d* mg_angu = mg->get_angu() ;
515  Valeur bc(mg_angu) ;
516 
517  Scalar rr(mp) ;
518  rr = mp.r ;
519  rr.std_spectral_base() ;
520  Scalar st(mp) ;
521  st = mp.sint ;
522  st.std_spectral_base() ;
523  Scalar ct(mp) ;
524  ct = mp.cost ;
525  ct.std_spectral_base() ;
526  Scalar sp(mp) ;
527  sp = mp.sinp ;
528  sp.std_spectral_base() ;
529  Scalar cp(mp) ;
530  cp = mp.cosp ;
531  cp.std_spectral_base() ;
532 
533  Vector ll(mp, CON, mp.get_bvect_cart()) ;
534  ll.set_etat_qcq() ;
535  ll.set(1) = st % cp ;
536  ll.set(2) = st % sp ;
537  ll.set(3) = ct ;
538  ll.std_spectral_base() ;
539 
540  Scalar divshift(mp) ; // dzpuis = 2
541  divshift = shift_auto_rs(1).deriv(1) + shift_auto_rs(2).deriv(2)
542  + shift_auto_rs(3).deriv(3) + d_shift_comp(1,1) + d_shift_comp(2,2)
543  + d_shift_comp(3,3) ;
544  divshift.std_spectral_base() ;
545 
546  Scalar llshift(mp) ; // dzpuis = 0
547  llshift = ll(1) % (shift_auto_rs(1) + shift_comp(1))
548  + ll(2) % (shift_auto_rs(2) + shift_comp(2))
549  + ll(3) % (shift_auto_rs(3) + shift_comp(3)) ;
550  llshift.std_spectral_base() ;
551 
552  Scalar llshift_auto_rs(mp) ; // dzpuis = 0
553  llshift_auto_rs = ll(1)%shift_auto_rs(1) + ll(2)%shift_auto_rs(2)
554  + ll(3)%shift_auto_rs(3) ;
555  llshift_auto_rs.std_spectral_base() ;
556 
557  Scalar lldllsh = llshift_auto_rs.dsdr()
558  + ll(1) * ( ll(1)%d_shift_comp(1,1) + ll(2)%d_shift_comp(1,2)
559  + ll(3)%d_shift_comp(1,3) )
560  + ll(2) * ( ll(1)%d_shift_comp(2,1) + ll(2)%d_shift_comp(2,2)
561  + ll(3)%d_shift_comp(2,3) )
562  + ll(3) * ( ll(1)%d_shift_comp(3,1) + ll(2)%d_shift_comp(3,2)
563  + ll(3)%d_shift_comp(3,3) ) ; // dzpuis = 2
564  lldllsh.std_spectral_base() ;
565 
566  Scalar tmp2 = divshift ;
567  Scalar tmp3 = -3.*lldllsh ;
568 
569  tmp2.dec_dzpuis(2) ;
570  tmp3.dec_dzpuis(2) ;
571 
572  Scalar tmp(mp) ;
573 
574  double mass = ggrav * mass_bh ;
575 
576  if (kerrschild) {
577 
578  cout << "!!!!! WARNING: Not yet available !!!!!" << endl ;
579  abort() ;
580 
581  }
582  else { // Isotropic coordinates with the maximal slicing
583 
584  double cc ;
585 
586  // Sets C/M^2 for each case of the lapse boundary condition
587  // --------------------------------------------------------
588 
589  if (bc_lapconf_nd) { // Neumann boundary condition
590  if (bc_lapconf_fs) { // First condition
591  // d(\alpha \psi)/dr = 0
592  // ---------------------
593  cc = 2. * (sqrt(13.) - 1.) / 3. ;
594  }
595  else { // Second condition
596  // d(\alpha \psi)/dr = (\alpha \psi)/(2 rah)
597  // -----------------------------------------
598  cc = 4. / 3. ;
599  }
600  }
601  else { // Dirichlet boundary condition
602  if (bc_lapconf_fs) { // First condition
603  // (\alpha \psi) = 1/2
604  // -------------------
605  cout << "!!!!! WARNING: Not yet prepared !!!!!" << endl ;
606  abort() ;
607  }
608  else { // Second condition
609  // (\alpha \psi) = 1/sqrt(2.) \psi_KS
610  // ----------------------------------
611  cout << "!!!!! WARNING: Not yet prepared !!!!!" << endl ;
612  abort() ;
613  // cc = 2. * sqrt(2.) ;
614  }
615  }
616 
617  Scalar r_are(mp) ;
619  r_are.std_spectral_base() ;
620 
621  Scalar tmp1 = - 0.5 * (confo_auto_rs + confo_comp) / rr ;
622  Scalar tmp7 = - ll(1)%d_confo_comp(1) - ll(2)%d_confo_comp(2)
623  - ll(3)%d_confo_comp(3) ;
624  tmp7.std_spectral_base() ;
625  tmp7.dec_dzpuis(2) ; // dzpuis : 2 -> 0
626 
627  /*
628  Scalar tmp8 = 0.5 * sqrt(1. - 2.*mass/r_are/rr
629  + cc*cc*pow(mass/r_are/rr,4.))
630  * (pow(confo_tot,3.)*mass*mass*cc/lapse_tot/pow(r_are*rr,3.)
631  - sqrt(r_are) / rr) ;
632  */
633  Scalar tmp8 = 0.125*cc*(0.25*cc*pow(confo_tot,4.)/r_are/lapconf_tot
634  - sqrt(r_are)) / rr ;
635  tmp8.std_spectral_base() ;
636 
637  tmp = tmp7 + tmp1
638  + pow(confo_tot,4.) * (tmp2 + tmp3) / 12. / lapconf_tot
639  + tmp8 ;
640 
641  }
642 
643  int nt = mg->get_nt(0) ;
644  int np = mg->get_np(0) ;
645 
646  bc = 1. ;
647  for (int j=0; j<nt; j++) {
648  for (int k=0; k<np; k++) {
649  bc.set(0,k,j,0) = tmp.val_grid_point(1,k,j,0) ;
650  }
651  }
652 
653  bc.std_base_scal() ;
654  return bc ;
655 
656 }
657 }
Map & mp
Mapping associated with the black hole.
Definition: blackhole.h:80
Base_val ** std_base_vect_cart() const
Returns the standard spectral bases for the Cartesian components of a vector.
virtual void set_etat_qcq()
Sets the logical state of all components to ETATQCQ (ordinary state).
Definition: tensor.C:490
const Valeur bc_shift_z() const
Boundary condition on the apparent horizon of the black hole for the shift vector of the direction: ...
Definition: hole_bhns_bc.C:406
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
Definition: grilles.h:479
Cmp sqrt(const Cmp &)
Square root.
Definition: cmp_math.C:223
double get_ori_y() const
Returns the y coordinate of the origin.
Definition: map.h:788
Lorene prototypes.
Definition: app_hor.h:67
double mass_bh
Gravitational mass of BH.
Definition: blackhole.h:88
Standard units of space, time and mass.
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:783
Tensor field of valence 0 (or component of a tensorial field).
Definition: scalar.h:393
bool kerrschild
true for a Kerr-Schild background, false for a conformally flat background
Definition: blackhole.h:85
double get_ori_x() const
Returns the x coordinate of the origin.
Definition: map.h:786
virtual void std_spectral_base()
Sets the spectral bases of the Valeur va to the standard ones for a scalar field. ...
Definition: scalar.C:790
Values and coefficients of a (real-value) function.
Definition: valeur.h:297
Tensor d_shift_comp
Derivative of the shift vector generated by the companion star.
Definition: hole_bhns.h:154
Tensor field of valence 1.
Definition: vector.h:188
void std_base_scal()
Sets the bases for spectral expansions (member base ) to the standard ones for a scalar.
Definition: valeur.C:827
double val_grid_point(int l, int k, int j, int i) const
Returns the value of the field at a specified grid point.
Definition: scalar.h:643
Scalar lapconf_comp
Lapconf function generated by the companion star.
Definition: hole_bhns.h:98
const Valeur bc_shift_x(double ome_orb, double y_rot) const
Boundary condition on the apparent horizon of the black hole for the shift vector of the direction: ...
Definition: hole_bhns_bc.C:186
Coord sint
Definition: map.h:739
const Valeur bc_shift_y(double ome_orb, double x_rot) const
Boundary condition on the apparent horizon of the black hole for the shift vector of the direction: ...
Definition: hole_bhns_bc.C:296
Scalar confo_auto_rs
Part of the conformal factor from the numerical computation.
Definition: hole_bhns.h:157
Vector d_lapconf_comp
Derivative of the lapconf function generated by the companion star.
Definition: hole_bhns.h:123
virtual void inc_dzpuis(int inc=1)
Increases by inc units the value of dzpuis and changes accordingly the values of the Scalar in the co...
const Mg3d * get_angu() const
Returns the pointer on the associated angular grid.
Definition: mg3d.C:604
Scalar confo_comp
Conformal factor generated by the companion star.
Definition: hole_bhns.h:166
Base_val base
Bases on which the spectral expansion is performed.
Definition: valeur.h:315
Vector shift_auto_rs
Part of the shift vector from the numerical computation.
Definition: hole_bhns.h:126
Coord sinp
Definition: map.h:741
Vector shift_comp
Shift vector generated by the companion star.
Definition: hole_bhns.h:135
Cmp pow(const Cmp &, int)
Power .
Definition: cmp_math.C:351
Scalar confo_tot
Total conformal factor.
Definition: hole_bhns.h:169
Multi-domain grid.
Definition: grilles.h:279
Bases of the spectral expansions.
Definition: base_val.h:325
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
Scalar & set(const Itbl &ind)
Returns the value of a component (read/write version).
Definition: tensor.C:663
Coord y
y coordinate centered on the grid
Definition: map.h:745
double omega_spin
Spin angular velocity of the black hole.
Definition: hole_bhns.h:86
Vector d_confo_comp
Derivative of the conformal factor generated by the companion star.
Definition: hole_bhns.h:185
Coord cosp
Definition: map.h:742
Scalar lapconf_auto_rs
Part of the lapconf function from the numerical computation.
Definition: hole_bhns.h:89
const Scalar & dsdr() const
Returns of *this .
Definition: scalar_deriv.C:113
Coord x
x coordinate centered on the grid
Definition: map.h:744
bool bc_lapconf_fs
true for the first type BC for the lapconf function, false for the second type BH ...
Definition: hole_bhns.h:78
bool bc_lapconf_nd
true for the Neumann type BC for the lapconf function, false for the Dirichlet type BH ...
Definition: hole_bhns.h:73
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
Definition: grilles.h:474
const Scalar r_coord(bool neumann, bool first) const
Expresses the areal radial coordinate by that in spatially isotropic coordinates. ...
virtual void dec_dzpuis(int dec=1)
Decreases by dec units the value of dzpuis and changes accordingly the values of the Scalar in the co...
const Valeur bc_lapconf() const
Boundary condition on the apparent horizon of the black hole for the lapconf function: 2-D Valeur...
Definition: hole_bhns_bc.C:71
const Valeur bc_confo() const
Boundary condition on the apparent horizon of the black hole for the conformal factor: 2-D Valeur...
Definition: blackhole_bc.C:413
Scalar lapconf_tot
Total lapconf function.
Definition: hole_bhns.h:101
Coord r
r coordinate centered on the grid
Definition: map.h:736
Tbl & set(int l)
Read/write of the value in a given domain (configuration space).
Definition: valeur.h:373
Coord cost
Definition: map.h:740