LORENE
et_bin_nsbh_angul.C
1/*
2 * Copyright (c) 2005 Philippe Grandclément
3 *
4 * This file is part of LORENE.
5 *
6 * LORENE is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation; either version 2 of the License, or
9 * (at your option) any later version.
10 *
11 * LORENE is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with LORENE; if not, write to the Free Software
18 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
19 *
20 */
21
22
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 $" ;
24
25/*
26 * $Id: et_bin_nsbh_angul.C,v 1.5 2014/10/13 08:52:56 j_novak Exp $
27 * $Log: et_bin_nsbh_angul.C,v $
28 * Revision 1.5 2014/10/13 08:52:56 j_novak
29 * Lorene classes and functions now belong to the namespace Lorene.
30 *
31 * Revision 1.4 2014/10/06 15:13:08 j_novak
32 * Modified #include directives to use c++ syntax.
33 *
34 * Revision 1.3 2005/08/31 09:53:19 m_saijo
35 * Delete unnecessary words around headers
36 *
37 * Revision 1.2 2005/08/31 09:13:45 p_grandclement
38 * add math.h
39 *
40 * Revision 1.1 2005/08/29 15:10:16 p_grandclement
41 * Addition of things needed :
42 * 1) For BBH with different masses
43 * 2) Provisory files for the mixted binaries (Bh and NS) : THIS IS NOT
44 * WORKING YET !!!
45 *
46 *
47 * $Header: /cvsroot/Lorene/C++/Source/Etoile/et_bin_nsbh_angul.C,v 1.5 2014/10/13 08:52:56 j_novak Exp $
48 *
49 */
50
51// C Headers
52#include <cmath>
53
54// Headers Lorene
55#include "et_bin_nsbh.h"
56#include "graphique.h"
57
58namespace Lorene {
59double Et_bin_nsbh::compute_angul() const {
60
61 // On récupère les trucs qui vont servir :
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) ;
65
66 Cmp xabs (mp) ;
67 xabs = mp.xa ;
68
69 Cmp yabs (mp) ;
70 yabs = mp.ya ;
71
72 // Derivee_square
73 Cmp G_square_cmp (-pow(confpsi(), 4)/nnn()/nnn()*(xabs*xabs+yabs*yabs)) ;
74 G_square_cmp.std_base_scal() ;
75 double dG_square = G_square_cmp.dsdx()(0,0,0,0) ;
76
77 // Derive single
78 Cmp G_single_cmp (-2*pow(confpsi(), 4)/nnn()/nnn()*(shift(1)*xabs-shift(0)*yabs)) ;
79 G_single_cmp.std_base_scal() ;
80 double dG_single = G_single_cmp.dsdx()(0,0,0,0) ;
81
82 // Derive const
83 Cmp G_const_cmp = 1-pow(confpsi(), 4)/nnn()/nnn()*(shift(0)*shift(0) + shift(1)*shift(1) + shift(2)*shift(2)) ;
84 G_const_cmp.std_base_scal() ;
85 double dG_const = G_const_cmp.dsdx()(0,0,0,0) ;
86
87 // coefficients de G
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) ;
91
92 // Les coefficients du polynome :
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;
96
97 double determinant = b_coef*b_coef - 4*a_coef*c_coef ;
98
99 bool soluce = true ;
100 double res ;
101
102 if (determinant <0) {
103 soluce = false ;
104 }
105
106 else {
107 double sol_un = (-b_coef - sqrt(determinant))/2/a_coef ;
108 double sol_deux = (-b_coef + sqrt(determinant))/2/a_coef ;
109
110 bool signe_un = (sol_un >0) ? true : false ;
111 bool signe_deux = (sol_deux >0) ? true : false ;
112
113
114 if (signe_un == signe_deux) {
115 soluce = false ;
116 }
117 else {
118 res = (signe_un) ? sol_un : sol_deux ;
119 }
120}
121
122 if (soluce == false)
123 return 0 ;
124 else
125 return res ;
126}
127}
Tenseur confpsi
Total conformal factor $\Psi$.
Tenseur loggam
Logarithm of the Lorentz factor between the fluid and the co-orbiting observer.
Definition etoile.h:849
Tenseur nnn
Total lapse function.
Definition etoile.h:509
Map & mp
Mapping associated with the star.
Definition etoile.h:429
Tenseur shift
Total shift vector.
Definition etoile.h:512
Coord ya
Absolute y coordinate.
Definition map.h:731
Coord xa
Absolute x coordinate.
Definition map.h:730
Cmp sqrt(const Cmp &)
Square root.
Definition cmp_math.C:220
Cmp pow(const Cmp &, int)
Power .
Definition cmp_math.C:348
Lorene prototypes.
Definition app_hor.h:64