LORENE
bin_bhns_extr_orbit.C
1/*
2 * Method of class Bin_bhns_extr to compute the orbital angular velocity
3 * {\tt omega}
4 *
5 * (see file bin_bhns_extr.h for documentation).
6 *
7 */
8
9/*
10 * Copyright (c) 2004-2005 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: bin_bhns_extr_orbit.C,v 1.6 2016/12/05 16:17:46 j_novak Exp $
33 * $Log: bin_bhns_extr_orbit.C,v $
34 * Revision 1.6 2016/12/05 16:17:46 j_novak
35 * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
36 *
37 * Revision 1.5 2014/10/24 14:10:23 j_novak
38 * Minor change to prevent weird error from g++-4.8...
39 *
40 * Revision 1.4 2014/10/13 08:52:42 j_novak
41 * Lorene classes and functions now belong to the namespace Lorene.
42 *
43 * Revision 1.3 2014/10/06 15:13:00 j_novak
44 * Modified #include directives to use c++ syntax.
45 *
46 * Revision 1.2 2005/02/28 23:07:47 k_taniguchi
47 * Addition of the code for the conformally flat case
48 *
49 * Revision 1.1 2004/11/30 20:46:57 k_taniguchi
50 * *** empty log message ***
51 *
52 *
53 * $Header: /cvsroot/Lorene/C++/Source/Bin_bhns_extr/bin_bhns_extr_orbit.C,v 1.6 2016/12/05 16:17:46 j_novak Exp $
54 *
55 */
56
57// C headers
58#include <cmath>
59
60// Lorene headers
61#include "bin_bhns_extr.h"
62#include "eos.h"
63#include "param.h"
64#include "utilitaires.h"
65#include "unites.h"
66
67namespace Lorene {
68double fonc_bhns_orbit_ks(double, const Param& ) ;
69double fonc_bhns_orbit_cf(double, const Param& ) ;
70
71//************************************************************************
72
73void Bin_bhns_extr::orbit_omega_ks(double fact_omeg_min,
74 double fact_omeg_max) {
75
76 using namespace Unites ;
77
78 //------------------------------------------------------------------
79 // Evaluation of various quantities at the center of a neutron star
80 //------------------------------------------------------------------
81
82 double dnulg, asn2, dasn2, nx, dnx, ny, dny, npn, dnpn ;
83 double msr ;
84
85 const Map& mp = star.get_mp() ;
86
87 const Cmp& logn_auto_regu = star.get_logn_auto_regu()() ;
88 const Cmp& loggam = star.get_loggam()() ;
89 const Cmp& nnn = star.get_nnn()() ;
90 const Cmp& a_car = star.get_a_car()() ;
91 const Tenseur& shift = star.get_shift() ;
92 const Tenseur& d_logn_auto_div = star.get_d_logn_auto_div() ;
93
94 Tenseur dln_auto_div = d_logn_auto_div ;
95
96 if ( *(dln_auto_div.get_triad()) != ref_triad ) {
97
98 // Change the basis from spherical coordinate to Cartesian one
99 dln_auto_div.change_triad( mp.get_bvect_cart() ) ;
100
101 // Change the basis from mapping coordinate to absolute one
102 dln_auto_div.change_triad( ref_triad ) ;
103
104 }
105
106 const Tenseur& d_logn_comp = star.get_d_logn_comp() ;
107
108 Tenseur dln_comp = d_logn_comp ;
109
110 if ( *(dln_comp.get_triad()) != ref_triad ) {
111
112 // Change the basis from spherical coordinate to Cartesian one
113 dln_comp.change_triad( mp.get_bvect_cart() ) ;
114
115 // Change the basis from mapping coordinate to absolute one
116 dln_comp.change_triad( ref_triad ) ;
117
118 }
119
120 //-------------------------------
121 // Parameters with respect to BH
122 //-------------------------------
123
124 msr = ggrav * mass_bh / separ ;
125
126 //-----------------------------------------------------------------
127 // Calculation of d/dX( nu + ln(Gamma) ) at the center of the star
128 // ---> dnulg
129 //-----------------------------------------------------------------
130
131 // Factor for the coordinate transformation x --> X :
132 double factx ;
133 if (fabs(mp.get_rot_phi()) < 1.e-14) {
134 factx = 1. ;
135 }
136 else {
137 if (fabs(mp.get_rot_phi() - M_PI) < 1.e-14) {
138 factx = - 1. ;
139 }
140 else {
141 cout << "Bin_bhns_extr::orbit_omega : unknown value of rot_phi !"
142 << endl ;
143 abort() ;
144 }
145 }
146
147 Cmp tmp = logn_auto_regu + loggam ;
148
149 // ... gradient
150 dnulg = dln_auto_div(0)(0, 0, 0, 0) + dln_comp(0)(0, 0, 0, 0)
151 + factx * tmp.dsdx()(0, 0, 0, 0) ;
152
153 //------------------------------------------------------------
154 // Calculation of A^2/N^2 at the center of the star ---> asn2
155 //------------------------------------------------------------
156 double nc = nnn(0, 0, 0, 0) ;
157 double a2c = a_car(0, 0, 0, 0) ;
158 asn2 = a2c / (nc * nc) ;
159
160 if ( star.is_relativistic() ) {
161
162 //------------------------------------------------------------------
163 // Calculation of d/dX(A^2/N^2) at the center of the star ---> dasn
164 //------------------------------------------------------------------
165
166 double da2c = factx * a_car.dsdx()(0, 0, 0, 0) ;
167 double dnc = factx * nnn.dsdx()(0, 0, 0, 0) ;
168
169 dasn2 = ( da2c - 2. * a2c / nc * dnc ) / (nc*nc) ;
170
171 //------------------------------------------------------
172 // Calculation of N^X at the center of the star ---> nx
173 //------------------------------------------------------
174 nx = shift(0)(0, 0, 0, 0) ;
175
176 //-----------------------------------------------------------
177 // Calculation of dN^X/dX at the center of the star ---> dnx
178 //-----------------------------------------------------------
179 dnx = factx * shift(0).dsdx()(0, 0, 0, 0) ;
180
181 //------------------------------------------------------
182 // Calculation of N^Y at the center of the star ---> ny
183 //------------------------------------------------------
184 ny = shift(1)(0, 0, 0, 0) ;
185
186 //-----------------------------------------------------------
187 // Calculation of dN^Y/dX at the center of the star ---> dny
188 //-----------------------------------------------------------
189 dny = factx * shift(1).dsdx()(0, 0, 0, 0) ;
190
191 //--------------------------------------------
192 // Calculation of (N^X)^2 + (N^Y)^2 + (N^Z)^2
193 // at the center of the star ---> npn
194 //--------------------------------------------
195 tmp = flat_scalar_prod(shift, shift)() ; // For the flat part
196 npn = tmp(0, 0, 0, 0) ;
197
198 //----------------------------------------------------
199 // Calculation of d/dX( (N^X)^2 + (N^Y)^2 + (N^Z)^2 )
200 // at the center of the star ---> dnpn
201 //----------------------------------------------------
202 dnpn = factx * tmp.dsdx()(0, 0, 0, 0) ;
203
204 } // Finish of the relativistic case
205 else {
206 cout << "Bin_bhns_extr::orbit_omega : "
207 << "It should be the relativistic calculation !" << endl ;
208 abort() ;
209 }
210
211 cout << "Bin_bhns_extr::orbit_omega: coord. ori. d(nu+log(Gam))/dX : "
212 << dnulg << endl ;
213 cout << "Bin_bhns_extr::orbit_omega: coord. ori. A^2/N^2 : "
214 << asn2 << endl ;
215 cout << "Bin_bhns_extr::orbit_omega: coord. ori. d(A^2/N^2)/dX : "
216 << dasn2 << endl ;
217 cout << "Bin_bhns_extr::orbit_omega: coord. ori. N^X : " << nx << endl ;
218 cout << "Bin_bhns_extr::orbit_omega: coord. ori. dN^X/dX : "
219 << dnx << endl ;
220 cout << "Bin_bhns_extr::orbit_omega: coord. ori. N^Y : " << ny << endl ;
221 cout << "Bin_bhns_extr::orbit_omega: coord. ori. dN^Y/dX : "
222 << dny << endl ;
223 cout << "Bin_bhns_extr::orbit_omega: coord. ori. N.N : " << npn << endl ;
224 cout << "Bin_bhns_extr::orbit_omega: coord. ori. d(N.N)/dX : "
225 << dnpn << endl ;
226
227 //------------------------------------------------------
228 // Start of calculation of the orbital angular velocity
229 //------------------------------------------------------
230 int relat = ( star.is_relativistic() ) ? 1 : 0 ;
231
232 double ori_x = (star.get_mp()).get_ori_x() ;
233 Param parf ;
234 parf.add_int(relat) ;
235 parf.add_double( ori_x, 0) ;
236 parf.add_double( dnulg, 1) ;
237 parf.add_double( asn2, 2) ;
238 parf.add_double( dasn2, 3) ;
239 parf.add_double( nx, 4) ;
240 parf.add_double( dnx, 5) ;
241 parf.add_double( ny, 6) ;
242 parf.add_double( dny, 7) ;
243 parf.add_double( npn, 8) ;
244 parf.add_double( dnpn, 9) ;
245 parf.add_double( msr, 10) ;
246
247 double omega1 = fact_omeg_min * omega ;
248 double omega2 = fact_omeg_max * omega ;
249
250 cout << "Bin_bhns_extr::orbit_omega: omega1, omega2 [rad/s] : "
251 << omega1 * f_unit << " " << omega2 * f_unit << endl ;
252
253 // Search for the various zeros in the interval [omega1,omega2]
254 // ------------------------------------------------------------
255 int nsub = 50 ; // total number of subdivisions of the interval
256 Tbl* azer = 0x0 ;
257 Tbl* bzer = 0x0 ;
258 zero_list(fonc_bhns_orbit_ks, parf, omega1, omega2, nsub,
259 azer, bzer) ;
260
261 // Search for the zero closest to the previous value of omega
262 // ----------------------------------------------------------
263 double omeg_min, omeg_max ;
264 int nzer = azer->get_taille() ; // number of zeros found by zero_list
265 cout << "Bin_bhns_extr::orbit_omega : " << nzer <<
266 "zero(s) found in the interval [omega1, omega2]." << endl ;
267 cout << "omega, omega1, omega2 : " << omega << " " << omega1
268 << " " << omega2 << endl ;
269 cout << "azer : " << *azer << endl ;
270 cout << "bzer : " << *bzer << endl ;
271
272 if (nzer == 0) {
273 cout << "Bin_bhns_extr::orbit_omega: WARNING : "
274 << "no zero detected in the interval" << endl
275 << " [" << omega1 * f_unit << ", "
276 << omega2 * f_unit << "] rad/s !" << endl ;
277 omeg_min = omega1 ;
278 omeg_max = omega2 ;
279 }
280 else {
281 double dist_min = fabs(omega2 - omega1) ;
282 int i_dist_min = -1 ;
283 for (int i=0; i<nzer; i++) {
284 // Distance of previous value of omega from the center of the
285 // interval [azer(i), bzer(i)]
286 double dist = fabs( omega - 0.5 * ( (*azer)(i) + (*bzer)(i) ) ) ;
287
288 if (dist < dist_min) {
289 dist_min = dist ;
290 i_dist_min = i ;
291 }
292 }
293 omeg_min = (*azer)(i_dist_min) ;
294 omeg_max = (*bzer)(i_dist_min) ;
295 }
296
297 delete azer ; // Tbl allocated by zero_list
298 delete bzer ; //
299
300 cout << "Bin_bhns_extr:orbit_omega : "
301 << "interval selected for the search of the zero : "
302 << endl << " [" << omeg_min << ", " << omeg_max << "] = ["
303 << omeg_min * f_unit << ", " << omeg_max * f_unit << "] rad/s "
304 << endl ;
305
306 // Computation of the zero in the selected interval by the secant method
307 // ---------------------------------------------------------------------
308
309 int nitermax = 200 ;
310 int niter ;
311 double precis = 1.e-13 ;
312 omega = zerosec_b(fonc_bhns_orbit_ks, parf, omeg_min, omeg_max,
313 precis, nitermax, niter) ;
314
315 cout << "Bin_bhns_extr::orbit_omega : "
316 << "Number of iterations in zerosec for omega : "
317 << niter << endl ;
318
319 cout << "Bin_bhns_extr::orbit_omega : omega [rad/s] : "
320 << omega * f_unit << endl ;
321
322}
323
324
325void Bin_bhns_extr::orbit_omega_cf(double fact_omeg_min,
326 double fact_omeg_max) {
327
328 using namespace Unites ;
329
330 //------------------------------------------------------------------
331 // Evaluation of various quantities at the center of a neutron star
332 //------------------------------------------------------------------
333
334 double dnulg, asn2, dasn2, ny, dny, npn, dnpn ;
335
336 const Map& mp = star.get_mp() ;
337
338 const Cmp& logn_auto_regu = star.get_logn_auto_regu()() ;
339 const Cmp& logn_comp = star.get_logn_comp()() ;
340 const Cmp& loggam = star.get_loggam()() ;
341 const Cmp& nnn = star.get_nnn()() ;
342 const Cmp& a_car = star.get_a_car()() ;
343 const Tenseur& shift = star.get_shift() ;
344 const Tenseur& d_logn_auto_div = star.get_d_logn_auto_div() ;
345
346 Tenseur dln_auto_div = d_logn_auto_div ;
347
348 if ( *(dln_auto_div.get_triad()) != ref_triad ) {
349
350 // Change the basis from spherical coordinate to Cartesian one
351 dln_auto_div.change_triad( mp.get_bvect_cart() ) ;
352
353 // Change the basis from mapping coordinate to absolute one
354 dln_auto_div.change_triad( ref_triad ) ;
355
356 }
357
358 //-----------------------------------------------------------------
359 // Calculation of d/dX( nu + ln(Gamma) ) at the center of the star
360 // ---> dnulg
361 //-----------------------------------------------------------------
362
363 // Factor for the coordinate transformation x --> X :
364 double factx ;
365 if (fabs(mp.get_rot_phi()) < 1.e-14) {
366 factx = 1. ;
367 }
368 else {
369 if (fabs(mp.get_rot_phi() - M_PI) < 1.e-14) {
370 factx = - 1. ;
371 }
372 else {
373 cout <<
374 "Bin_bhns_extr::orbit_omega_cfc : unknown value of rot_phi !"
375 << endl ;
376 abort() ;
377 }
378 }
379
380 Cmp tmp = logn_auto_regu + logn_comp + loggam ;
381
382 // ... gradient
383 dnulg = dln_auto_div(0)(0, 0, 0, 0)
384 + factx * tmp.dsdx()(0, 0, 0, 0) ;
385
386 //------------------------------------------------------------
387 // Calculation of A^2/N^2 at the center of the star ---> asn2
388 //------------------------------------------------------------
389 double nc = nnn(0, 0, 0, 0) ;
390 double a2c = a_car(0, 0, 0, 0) ;
391 asn2 = a2c / (nc * nc) ;
392
393 if ( star.is_relativistic() ) {
394
395 //------------------------------------------------------------------
396 // Calculation of d/dX(A^2/N^2) at the center of the star ---> dasn
397 //------------------------------------------------------------------
398
399 double da2c = factx * a_car.dsdx()(0, 0, 0, 0) ;
400 double dnc = factx * nnn.dsdx()(0, 0, 0, 0) ;
401
402 dasn2 = ( da2c - 2. * a2c / nc * dnc ) / (nc*nc) ;
403
404 //------------------------------------------------------
405 // Calculation of N^Y at the center of the star ---> ny
406 //------------------------------------------------------
407 ny = shift(1)(0, 0, 0, 0) ;
408
409 //-----------------------------------------------------------
410 // Calculation of dN^Y/dX at the center of the star ---> dny
411 //-----------------------------------------------------------
412 dny = factx * shift(1).dsdx()(0, 0, 0, 0) ;
413
414 //--------------------------------------------
415 // Calculation of (N^X)^2 + (N^Y)^2 + (N^Z)^2
416 // at the center of the star ---> npn
417 //--------------------------------------------
418 tmp = flat_scalar_prod(shift, shift)() ; // For the flat part
419 npn = tmp(0, 0, 0, 0) ;
420
421 //----------------------------------------------------
422 // Calculation of d/dX( (N^X)^2 + (N^Y)^2 + (N^Z)^2 )
423 // at the center of the star ---> dnpn
424 //----------------------------------------------------
425 dnpn = factx * tmp.dsdx()(0, 0, 0, 0) ;
426
427 } // Finish of the relativistic case
428 else {
429 cout << "Bin_bhns_extr::orbit_omega_cfc : "
430 << "It should be the relativistic calculation !" << endl ;
431 abort() ;
432 }
433
434 cout << "Bin_bhns_extr::orbit_omega_cfc: coord. ori. d(nu+log(Gam))/dX : "
435 << dnulg << endl ;
436 cout << "Bin_bhns_extr::orbit_omega_cfc: coord. ori. A^2/N^2 : "
437 << asn2 << endl ;
438 cout << "Bin_bhns_extr::orbit_omega_cfc: coord. ori. d(A^2/N^2)/dX : "
439 << dasn2 << endl ;
440 cout << "Bin_bhns_extr::orbit_omega_cfc: coord. ori. N^Y : "
441 << ny << endl ;
442 cout << "Bin_bhns_extr::orbit_omega_cfc: coord. ori. dN^Y/dX : "
443 << dny << endl ;
444 cout << "Bin_bhns_extr::orbit_omega_cfc: coord. ori. N.N : "
445 << npn << endl ;
446 cout << "Bin_bhns_extr::orbit_omega_cfc: coord. ori. d(N.N)/dX : "
447 << dnpn << endl ;
448
449 //------------------------------------------------------
450 // Start of calculation of the orbital angular velocity
451 //------------------------------------------------------
452 int relat = ( star.is_relativistic() ) ? 1 : 0 ;
453 double ori_x = (star.get_mp()).get_ori_x() ;
454
455 Param parf ;
456 parf.add_int(relat) ;
457 parf.add_double( ori_x, 0) ;
458 parf.add_double( dnulg, 1) ;
459 parf.add_double( asn2, 2) ;
460 parf.add_double( dasn2, 3) ;
461 parf.add_double( ny, 4) ;
462 parf.add_double( dny, 5) ;
463 parf.add_double( npn, 6) ;
464 parf.add_double( dnpn, 7) ;
465
466 double omega1 = fact_omeg_min * omega ;
467 double omega2 = fact_omeg_max * omega ;
468
469 cout << "Bin_bhns_extr::orbit_omega_cfc: omega1, omega2 [rad/s] : "
470 << omega1 * f_unit << " " << omega2 * f_unit << endl ;
471
472 // Search for the various zeros in the interval [omega1,omega2]
473 // ------------------------------------------------------------
474 int nsub = 50 ; // total number of subdivisions of the interval
475 Tbl* azer = 0x0 ;
476 Tbl* bzer = 0x0 ;
477 zero_list(fonc_bhns_orbit_cf, parf, omega1, omega2, nsub,
478 azer, bzer) ;
479
480 // Search for the zero closest to the previous value of omega
481 // ----------------------------------------------------------
482 double omeg_min, omeg_max ;
483 int nzer = azer->get_taille() ; // number of zeros found by zero_list
484 cout << "Bin_bhns_extr::orbit_omega_cfc : " << nzer <<
485 "zero(s) found in the interval [omega1, omega2]." << endl ;
486 cout << "omega, omega1, omega2 : " << omega << " " << omega1
487 << " " << omega2 << endl ;
488 cout << "azer : " << *azer << endl ;
489 cout << "bzer : " << *bzer << endl ;
490
491 if (nzer == 0) {
492 cout << "Bin_bhns_extr::orbit_omega_cfc: WARNING : "
493 << "no zero detected in the interval" << endl
494 << " [" << omega1 * f_unit << ", "
495 << omega2 * f_unit << "] rad/s !" << endl ;
496 omeg_min = omega1 ;
497 omeg_max = omega2 ;
498 }
499 else {
500 double dist_min = fabs(omega2 - omega1) ;
501 int i_dist_min = -1 ;
502 for (int i=0; i<nzer; i++) {
503 // Distance of previous value of omega from the center of the
504 // interval [azer(i), bzer(i)]
505 double dist = fabs( omega - 0.5 * ( (*azer)(i) + (*bzer)(i) ) ) ;
506
507 if (dist < dist_min) {
508 dist_min = dist ;
509 i_dist_min = i ;
510 }
511 }
512 omeg_min = (*azer)(i_dist_min) ;
513 omeg_max = (*bzer)(i_dist_min) ;
514 }
515
516 delete azer ; // Tbl allocated by zero_list
517 delete bzer ; //
518
519 cout << "Bin_bhns_extr:orbit_omega_cfc : "
520 << "interval selected for the search of the zero : "
521 << endl << " [" << omeg_min << ", " << omeg_max << "] = ["
522 << omeg_min * f_unit << ", " << omeg_max * f_unit << "] rad/s "
523 << endl ;
524
525 // Computation of the zero in the selected interval by the secant method
526 // ---------------------------------------------------------------------
527
528 int nitermax = 200 ;
529 int niter ;
530 double precis = 1.e-13 ;
531 omega = zerosec_b(fonc_bhns_orbit_cf, parf, omeg_min, omeg_max,
532 precis, nitermax, niter) ;
533
534 cout << "Bin_bhns_extr::orbit_omega_cfc : "
535 << "Number of iterations in zerosec for omega : "
536 << niter << endl ;
537
538 cout << "Bin_bhns_extr::orbit_omega_cfc : omega [rad/s] : "
539 << omega * f_unit << endl ;
540
541}
542
543
544//***********************************************************
545// Function used for search of the orbital angular velocity
546//***********************************************************
547
548double fonc_bhns_orbit_ks(double om, const Param& parf) {
549
550 int relat = parf.get_int() ;
551
552 double xc = parf.get_double(0) ;
553 double dnulg = parf.get_double(1) ;
554 double asn2 = parf.get_double(2) ;
555 double dasn2 = parf.get_double(3) ;
556 double nx = parf.get_double(4) ;
557 double dnx = parf.get_double(5) ;
558 double ny = parf.get_double(6) ;
559 double dny = parf.get_double(7) ;
560 double npn = parf.get_double(8) ;
561 double dnpn = parf.get_double(9) ;
562 double msr = parf.get_double(10) ;
563
564 double om2 = om*om ;
565
566 double dphi_cent ;
567
568 if (relat == 1) {
569 double bpb = om2 * xc*xc - 2.*om * ny * xc + npn + 2.*msr * nx*nx;
570
571 dphi_cent = ( asn2* ( om* (ny + xc*dny) - om2*xc - 0.5*dnpn
572 - 2.*msr * nx * dnx + msr * nx*nx / xc )
573 - 0.5*bpb* dasn2 ) / ( 1 - asn2 * bpb ) ;
574 }
575 else {
576 cout << "Bin_bhns_extr::orbit_omega_ks : "
577 << "It should be the relativistic calculation !" << endl ;
578 abort() ;
579 }
580
581 return dnulg + dphi_cent ;
582
583}
584
585double fonc_bhns_orbit_cf(double om, const Param& parf) {
586
587 int relat = parf.get_int() ;
588
589 double xc = parf.get_double(0) ;
590 double dnulg = parf.get_double(1) ;
591 double asn2 = parf.get_double(2) ;
592 double dasn2 = parf.get_double(3) ;
593 double ny = parf.get_double(4) ;
594 double dny = parf.get_double(5) ;
595 double npn = parf.get_double(6) ;
596 double dnpn = parf.get_double(7) ;
597
598 double om2 = om*om ;
599
600 double dphi_cent ;
601
602 if (relat == 1) {
603 double bpb = om2 * xc*xc - 2.*om * ny * xc + npn ;
604
605 dphi_cent = ( asn2* ( om* (ny + xc*dny) - om2*xc - 0.5*dnpn )
606 - 0.5*bpb* dasn2 ) / ( 1 - asn2 * bpb ) ;
607 }
608 else {
609 cout << "Bin_bhns_extr::orbit_omega_cf : "
610 << "It should be the relativistic calculation !" << endl ;
611 abort() ;
612 }
613
614 return dnulg + dphi_cent ;
615
616}
617}
void orbit_omega_ks(double fact_omeg_min, double fact_omeg_max)
Computes the orbital angular velocity {\tt omega} in the Kerr-Schild background metric.
const Base_vect_cart ref_triad
Cartesian triad of the absolute reference frame.
Et_bin_bhns_extr star
Neutron star.
double separ
Absolute orbital separation between two centers of BH and NS.
void orbit_omega_cf(double fact_omeg_min, double fact_omeg_max)
Computes the orbital angular velocity {\tt omega} in the conformally flat background metric.
double mass_bh
Gravitational mass of BH.
double omega
Angular velocity with respect to an asymptotically inertial observer.
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition cmp.h:446
const Cmp & dsdx() const
Returns of *this , where .
Definition cmp_deriv.C:151
Parameter storage.
Definition param.h:125
void add_double(const double &x, int position=0)
Adds the the address of a new double to the list.
Definition param.C:318
const int & get_int(int position=0) const
Returns the reference of a int stored in the list.
Definition param.C:295
const double & get_double(int position=0) const
Returns the reference of a double stored in the list.
Definition param.C:364
void add_int(const int &n, int position=0)
Adds the address of a new int to the list.
Definition param.C:249
Basic array class.
Definition tbl.h:161
int get_taille() const
Gives the total size (ie dim.taille).
Definition tbl.h:397
Tensor handling *** DEPRECATED : use class Tensor instead ***.
Definition tenseur.h:304
const Base_vect * get_triad() const
Returns the vectorial basis (triad) on which the components are defined.
Definition tenseur.h:707
void change_triad(const Base_vect &new_triad)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
Definition tenseur.C:674
void zero_list(double(*f)(double, const Param &), const Param &par, double xmin, double xmax, int nsub, Tbl *&az, Tbl *&bz)
Locates approximatively all the zeros of a function in a given interval.
Definition zero_list.C:60
double zerosec_b(double(*f)(double, const Param &), const Param &par, double a, double b, double precis, int nitermax, int &niter)
Finding the zero a function on a bounded domain.
Tenseur flat_scalar_prod(const Tenseur &t1, const Tenseur &t2)
Scalar product of two Tenseur when the metric is : performs the contraction of the last index of t1 w...
Lorene prototypes.
Definition app_hor.h:67
Map(const Mg3d &)
Constructor from a multi-domain 3D grid.
Definition map.C:142
Standard units of space, time and mass.