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 
28 
29 
30 /*
31  * $Id: bin_bhns_extr_anashift.C,v 1.4 2016/12/05 16:17:46 j_novak Exp $
32  * $Log: bin_bhns_extr_anashift.C,v $
33  * Revision 1.4 2016/12/05 16:17:46 j_novak
34  * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
35  *
36  * Revision 1.3 2014/10/13 08:52:41 j_novak
37  * Lorene classes and functions now belong to the namespace Lorene.
38  *
39  * Revision 1.2 2014/10/06 15:13:00 j_novak
40  * Modified #include directives to use c++ syntax.
41  *
42  * Revision 1.1 2004/11/30 20:45:29 k_taniguchi
43  * *** empty log message ***
44  *
45  *
46  * $Header: /cvsroot/Lorene/C++/Source/Bin_bhns_extr/bin_bhns_extr_anashift.C,v 1.4 2016/12/05 16:17:46 j_novak Exp $
47  *
48  */
49 
50 // C headers
51 #include <cmath>
52 
53 // Lorene headers
54 #include "bin_bhns_extr.h"
55 #include "unites.h"
56 
57 namespace Lorene {
59 
60  using namespace Unites ;
61 
62  // BH-NS binary systems should be relativistic
63  // -------------------------------------------
64 
65  if ( !star.is_relativistic() ) {
66 
67  cout << "BH-NS binary systems should be relativistic !!!" << endl ;
68  abort() ;
69  }
70 
71  // Radius of the neutron star
72  double a0 = star.ray_eq() ;
73 
74  // G M Omega R
75  double www = ggrav * star.mass_g() * omega * separ ;
76  // Approximates the mass ratio -> 0
77 
78  const Map& mp = star.get_mp() ;
79  Tenseur tmp(mp) ;
80  Tenseur tmp_ext(mp) ;
81  int nzet = star.get_nzet() ;
82  int nzm1 = mp.get_mg()->get_nzone() - 1 ;
83 
84  //-------------------
85  // Irrotational case
86  //-------------------
87  // Since this formula is only an initial guess, we use it
88  // also for the corotating case.
89 
90  // Computation of w_shift
91  // ----------------------
93 
94  // X component
95  // -----------
96  star.set_w_shift().set(0) = 0. ;
97 
98  // Y component
99  // -----------
100 
101  // For the incompressible case :
102  tmp.set_etat_qcq() ;
103  tmp.set() = 6. * www / a0 * ( 1. - (mp.r)*(mp.r) / (3.*a0*a0) ) ;
104  tmp.set().annule(nzet, nzm1) ;
105  tmp.set_std_base() ;
106 
107  tmp_ext.set_etat_qcq() ;
108  tmp_ext.set() = 4. * www / mp.r ;
109  tmp_ext.set().annule(0, nzet-1) ;
110  tmp_ext.set_std_base() ;
111 
112  star.set_w_shift().set(1) = tmp() + tmp_ext() ;
113 
114  // Z component
115  // -----------
116  star.set_w_shift().set(2) = 0. ;
117 
118  // Sets the standard spectral bases for Cartesian components
120 
121  // Computation of khi_shift
122  //-------------------------
123 
124  tmp.set() = 2. * www / a0 * (mp.y)
125  * ( 1. - 3.*(mp.r)*(mp.r) / (5.*a0*a0) ) ;
126  tmp.set().annule(nzet, nzm1) ;
127  tmp.set_std_base() ;
128 
129  tmp_ext.set() = 0.8 * www * a0 * a0 * (mp.sint) * (mp.sinp)
130  / ((mp.r)*(mp.r)) ;
131  tmp_ext.set().annule(0, nzet-1) ;
132  tmp_ext.set_std_base() ;
133 
134  star.set_khi_shift() = tmp + tmp_ext ;
135 
136  // Sets the standard spectral bases for a scalar field
138 
139 }
140 }
const Map & get_mp() const
Returns the mapping.
Definition: etoile.h:662
void annule(int l)
Sets the Cmp to zero in a given domain.
Definition: cmp.C:351
void set_std_base()
Set the standard spectal basis of decomposition for each component.
Definition: tenseur.C:1186
Lorene prototypes.
Definition: app_hor.h:67
Standard units of space, time and mass.
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:783
double ray_eq() const
Coordinate radius at , [r_unit].
Base class for coordinate mappings.
Definition: map.h:688
int get_nzet() const
Returns the number of domains occupied by the star.
Definition: etoile.h:665
virtual double mass_g() const
Gravitational mass.
Et_bin_bhns_extr star
Neutron star.
Definition: bin_bhns_extr.h:66
Coord sint
Definition: map.h:739
Cmp & set()
Read/write for a scalar (see also operator=(const Cmp&) ).
Definition: tenseur.C:840
Coord sinp
Definition: map.h:741
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:465
Tenseur & set_w_shift()
Read/write of w_shift.
Definition: etoile_bin.C:541
void analytical_shift()
Sets some analytical template for the shift vector (via the members { w_shift} and { khi_shift} of { ...
Coord y
y coordinate centered on the grid
Definition: map.h:745
double separ
Absolute orbital separation between two centers of BH and NS.
Definition: bin_bhns_extr.h:74
double omega
Angular velocity with respect to an asymptotically inertial observer.
Definition: bin_bhns_extr.h:71
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: tenseur.C:652
Tenseur & set_khi_shift()
Read/write of khi_shift.
Definition: etoile_bin.C:548
bool is_relativistic() const
Returns true for a relativistic star, false for a Newtonian one.
Definition: etoile.h:670
Tensor handling *** DEPRECATED : use class Tensor instead ***.
Definition: tenseur.h:304
Coord r
r coordinate centered on the grid
Definition: map.h:736