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