23char et_bin_nsbh_angul_C[] =
"$Header: /cvsroot/Lorene/C++/Source/Etoile/et_bin_nsbh_angul.C,v 1.5 2014/10/13 08:52:56 j_novak Exp $" ;
55#include "et_bin_nsbh.h"
59double Et_bin_nsbh::compute_angul()
const {
62 Cmp dx_mu (
nnn().dsdx() /
nnn());
63 Cmp dx_loggam (
loggam().dsdx()) ;
64 double part_dx = dx_mu(0,0,0,0) + dx_loggam(0,0,0,0) ;
74 G_square_cmp.std_base_scal() ;
75 double dG_square = G_square_cmp.dsdx()(0,0,0,0) ;
79 G_single_cmp.std_base_scal() ;
80 double dG_single = G_single_cmp.dsdx()(0,0,0,0) ;
84 G_const_cmp.std_base_scal() ;
85 double dG_const = G_const_cmp.dsdx()(0,0,0,0) ;
88 double G_square = G_square_cmp (0,0,0,0) ;
89 double G_single = G_single_cmp (0,0,0,0) ;
90 double G_const = G_const_cmp (0,0,0,0) ;
93 double a_coef = dG_square + 2*G_square*part_dx ;
94 double b_coef = dG_single + 2*G_single*part_dx ;
95 double c_coef = dG_const + 2*G_const *part_dx;
97 double determinant = b_coef*b_coef - 4*a_coef*c_coef ;
102 if (determinant <0) {
107 double sol_un = (-b_coef -
sqrt(determinant))/2/a_coef ;
108 double sol_deux = (-b_coef +
sqrt(determinant))/2/a_coef ;
110 bool signe_un = (sol_un >0) ?
true : false ;
111 bool signe_deux = (sol_deux >0) ?
true : false ;
114 if (signe_un == signe_deux) {
118 res = (signe_un) ? sol_un : sol_deux ;
Tenseur confpsi
Total conformal factor $\Psi$.
Tenseur loggam
Logarithm of the Lorentz factor between the fluid and the co-orbiting observer.
Tenseur nnn
Total lapse function.
Map & mp
Mapping associated with the star.
Tenseur shift
Total shift vector.
Coord ya
Absolute y coordinate.
Coord xa
Absolute x coordinate.
Cmp sqrt(const Cmp &)
Square root.
Cmp pow(const Cmp &, int)
Power .