LORENE
bin_bhns_extr_anashift.C
1/*
2 * Method of class Bin_bhns_extr to set some analytical form
3 *
4 * (see file bin_bhns_extr.h for documentation).
5 *
6 */
7
8/*
9 * Copyright (c) 2004 Keisuke Taniguchi
10 *
11 * This file is part of LORENE.
12 *
13 * LORENE is free software; you can redistribute it and/or modify
14 * it under the terms of the GNU General Public License version 2
15 * as published by the Free Software Foundation.
16 *
17 * LORENE is distributed in the hope that it will be useful,
18 * but WITHOUT ANY WARRANTY; without even the implied warranty of
19 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
20 * GNU General Public License for more details.
21 *
22 * You should have received a copy of the GNU General Public License
23 * along with LORENE; if not, write to the Free Software
24 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
25 *
26 */
27
28char bin_bhns_extr_anashift_C[] = "$Header: /cvsroot/Lorene/C++/Source/Bin_bhns_extr/bin_bhns_extr_anashift.C,v 1.3 2014/10/13 08:52:41 j_novak Exp $" ;
29
30/*
31 * $Id: bin_bhns_extr_anashift.C,v 1.3 2014/10/13 08:52:41 j_novak Exp $
32 * $Log: bin_bhns_extr_anashift.C,v $
33 * Revision 1.3 2014/10/13 08:52:41 j_novak
34 * Lorene classes and functions now belong to the namespace Lorene.
35 *
36 * Revision 1.2 2014/10/06 15:13:00 j_novak
37 * Modified #include directives to use c++ syntax.
38 *
39 * Revision 1.1 2004/11/30 20:45:29 k_taniguchi
40 * *** empty log message ***
41 *
42 *
43 * $Header: /cvsroot/Lorene/C++/Source/Bin_bhns_extr/bin_bhns_extr_anashift.C,v 1.3 2014/10/13 08:52:41 j_novak Exp $
44 *
45 */
46
47// C headers
48#include <cmath>
49
50// Lorene headers
51#include "bin_bhns_extr.h"
52#include "unites.h"
53
54namespace Lorene {
56
57 using namespace Unites ;
58
59 // BH-NS binary systems should be relativistic
60 // -------------------------------------------
61
62 if ( !star.is_relativistic() ) {
63
64 cout << "BH-NS binary systems should be relativistic !!!" << endl ;
65 abort() ;
66 }
67
68 // Radius of the neutron star
69 double a0 = star.ray_eq() ;
70
71 // G M Omega R
72 double www = ggrav * star.mass_g() * omega * separ ;
73 // Approximates the mass ratio -> 0
74
75 const Map& mp = star.get_mp() ;
76 Tenseur tmp(mp) ;
77 Tenseur tmp_ext(mp) ;
78 int nzet = star.get_nzet() ;
79 int nzm1 = mp.get_mg()->get_nzone() - 1 ;
80
81 //-------------------
82 // Irrotational case
83 //-------------------
84 // Since this formula is only an initial guess, we use it
85 // also for the corotating case.
86
87 // Computation of w_shift
88 // ----------------------
90
91 // X component
92 // -----------
93 star.set_w_shift().set(0) = 0. ;
94
95 // Y component
96 // -----------
97
98 // For the incompressible case :
99 tmp.set_etat_qcq() ;
100 tmp.set() = 6. * www / a0 * ( 1. - (mp.r)*(mp.r) / (3.*a0*a0) ) ;
101 tmp.set().annule(nzet, nzm1) ;
102 tmp.set_std_base() ;
103
104 tmp_ext.set_etat_qcq() ;
105 tmp_ext.set() = 4. * www / mp.r ;
106 tmp_ext.set().annule(0, nzet-1) ;
107 tmp_ext.set_std_base() ;
108
109 star.set_w_shift().set(1) = tmp() + tmp_ext() ;
110
111 // Z component
112 // -----------
113 star.set_w_shift().set(2) = 0. ;
114
115 // Sets the standard spectral bases for Cartesian components
117
118 // Computation of khi_shift
119 //-------------------------
120
121 tmp.set() = 2. * www / a0 * (mp.y)
122 * ( 1. - 3.*(mp.r)*(mp.r) / (5.*a0*a0) ) ;
123 tmp.set().annule(nzet, nzm1) ;
124 tmp.set_std_base() ;
125
126 tmp_ext.set() = 0.8 * www * a0 * a0 * (mp.sint) * (mp.sinp)
127 / ((mp.r)*(mp.r)) ;
128 tmp_ext.set().annule(0, nzet-1) ;
129 tmp_ext.set_std_base() ;
130
131 star.set_khi_shift() = tmp + tmp_ext ;
132
133 // Sets the standard spectral bases for a scalar field
135
136}
137}
Et_bin_bhns_extr star
Neutron star.
double separ
Absolute orbital separation between two centers of BH and NS.
void analytical_shift()
Sets some analytical template for the shift vector (via the members {\tt w_shift} and {\tt khi_shift}...
double omega
Angular velocity with respect to an asymptotically inertial observer.
void annule(int l)
Sets the Cmp to zero in a given domain.
Definition cmp.C:348
Tenseur & set_w_shift()
Read/write of w_shift.
Definition etoile_bin.C:538
virtual double mass_g() const
Gravitational mass.
Tenseur & set_khi_shift()
Read/write of khi_shift.
Definition etoile_bin.C:545
int get_nzet() const
Returns the number of domains occupied by the star.
Definition etoile.h:662
double ray_eq() const
Coordinate radius at , [r_unit].
const Map & get_mp() const
Returns the mapping.
Definition etoile.h:659
bool is_relativistic() const
Returns true for a relativistic star, false for a Newtonian one.
Definition etoile.h:667
Base class for coordinate mappings.
Definition map.h:670
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
Coord sinp
Definition map.h:723
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition map.h:765
int get_nzone() const
Returns the number of domains.
Definition grilles.h:448
Tensor handling *** DEPRECATED : use class Tensor instead ***.
Definition tenseur.h:301
Cmp & set()
Read/write for a scalar (see also operator=(const Cmp&) ).
Definition tenseur.C:824
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition tenseur.C:636
void set_std_base()
Set the standard spectal basis of decomposition for each component.
Definition tenseur.C:1170
Lorene prototypes.
Definition app_hor.h:64
Standard units of space, time and mass.