LORENE
star_bhns_kinema.C
1/*
2 * Method of class Star_bhns to compute kinematic quantities
3 *
4 * (see file star_bhns.h for documentation).
5 *
6 */
7
8/*
9 * Copyright (c) 2005-2007 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 star_bhns_kinema_C[] = "$Header: /cvsroot/Lorene/C++/Source/Star_bhns/star_bhns_kinema.C,v 1.4 2014/10/13 08:53:41 j_novak Exp $" ;
29
30/*
31 * $Id: star_bhns_kinema.C,v 1.4 2014/10/13 08:53:41 j_novak Exp $
32 * $Log: star_bhns_kinema.C,v $
33 * Revision 1.4 2014/10/13 08:53:41 j_novak
34 * Lorene classes and functions now belong to the namespace Lorene.
35 *
36 * Revision 1.3 2014/10/06 15:13:16 j_novak
37 * Modified #include directives to use c++ syntax.
38 *
39 * Revision 1.2 2008/05/15 19:16:06 k_taniguchi
40 * Change of a parameter.
41 *
42 * Revision 1.1 2007/06/22 01:32:00 k_taniguchi
43 * *** empty log message ***
44 *
45 *
46 * $Header: /cvsroot/Lorene/C++/Source/Star_bhns/star_bhns_kinema.C,v 1.4 2014/10/13 08:53:41 j_novak Exp $
47 *
48 */
49
50// C++ headers
51//#include <>
52
53// C headers
54#include <cmath>
55
56// Lorene headers
57#include "star_bhns.h"
58#include "unites.h"
59
60namespace Lorene {
61void Star_bhns::kinema_bhns(bool kerrschild, const double& mass_bh,
62 const double& sepa, double omega,
63 double x_rot, double y_rot) {
64
65 // Fundamental constants and units
66 // -------------------------------
67 using namespace Unites ;
68
69 int nz = mp.get_mg()->get_nzone() ;
70 int nzm1 = nz - 1 ;
71
72 //----------------------
73 // Computation of B^i/N
74 //----------------------
75
76 // 1/ Computation of omega m^i
77
78 const Coord& xa = mp.xa ;
79 const Coord& ya = mp.ya ;
80
81 // bsn.change_triad(mp.get_bvect_cart()) ;
82
83 if (fabs(mp.get_rot_phi()) < 1.e-10) {
84
85 bsn.set(1) = - omega * (ya - y_rot) ;
86 bsn.set(2) = omega * (xa - x_rot) ;
87 bsn.set(3) = 0. ;
88
89 }
90 else {
91
92 bsn.set(1) = omega * (ya - y_rot) ;
93 bsn.set(2) = - omega * (xa - x_rot) ;
94 bsn.set(3) = 0. ;
95
96 }
97
100
101 // 2/ Addition of shift_tot and division by lapse
102
103 for (int i=1; i<=3; i++) {
104 bsn.set(i) = confo_tot * ( bsn(i) + shift_tot(i) ) / lapconf_tot ;
105 }
108
109 //-----------------------------------------------------
110 // Lorentz factor between the co-orbiting ---> gam0
111 // observer and the Eulerian one
112 // See Eq (23) and (24) from Gourgoulhon et al. (2001)
113 //-----------------------------------------------------
114
115 Scalar bsn2(mp) ;
116 bsn2 = bsn(1)%bsn(1) + bsn(2)%bsn(2) + bsn(3)%bsn(3) ;
117 bsn2.std_spectral_base() ;
118
119 if (kerrschild) {
120
121 double mass = ggrav * mass_bh ;
122
123 Scalar xx(mp) ;
124 xx = mp.x ;
125 xx.std_spectral_base() ;
126 Scalar yy(mp) ;
127 yy = mp.y ;
128 yy.std_spectral_base() ;
129 Scalar zz(mp) ;
130 zz = mp.z ;
131 zz.std_spectral_base() ;
132
133 double yns = mp.get_ori_y() ;
134
135 Scalar rbh(mp) ;
136 rbh = sqrt( (xx+sepa)*(xx+sepa) + (yy+yns)*(yy+yns) + zz*zz ) ;
137 rbh.std_spectral_base() ;
138
139 Vector ll(mp, CON, mp.get_bvect_cart()) ;
140 ll.set_etat_qcq() ;
141 ll.set(1) = (xx+sepa) / rbh ;
142 ll.set(2) = (yy+yns) / rbh ;
143 ll.set(3) = zz / rbh ;
144 ll.std_spectral_base() ;
145
146 Scalar msr(mp) ;
147 msr = mass / rbh ;
148 msr.std_spectral_base() ;
149
150 Scalar llbsn(mp) ;
151 llbsn = ll(1)%bsn(1) + ll(2)%bsn(2) + ll(3)%bsn(3) ;
152 llbsn.std_spectral_base() ;
153
154 Scalar tmp1(mp) ;
155 tmp1 = 2. * msr % llbsn % llbsn ;
156 tmp1.std_spectral_base() ;
157
158 gam0 = 1. / sqrt(1. - psi4*(bsn2+tmp1)) ;
160
161 }
162 else { // Isotropic coordinates with the maximal slicing
163
164 gam0 = 1. / sqrt(1. - psi4%bsn2) ;
166
167 }
168
169 //-----------------------
170 // Centrifugal potential
171 //-----------------------
172
173 pot_centri = - log( gam0 ) ;
176
177 // The derived quantities are obsolete
178 // -----------------------------------
179
180 del_deriv() ;
181
182}
183}
Active physical coordinates and mapping derivatives.
Definition coord.h:90
Time evolution with partial storage (*** under development ***).
Definition evolution.h:371
const Base_vect_cart & get_bvect_cart() const
Returns the Cartesian basis associated with the coordinates (x,y,z) of the mapping,...
Definition map.h:791
Coord y
y coordinate centered on the grid
Definition map.h:727
Coord ya
Absolute y coordinate.
Definition map.h:731
double get_ori_y() const
Returns the y coordinate of the origin.
Definition map.h:770
Coord z
z coordinate centered on the grid
Definition map.h:728
Coord x
x coordinate centered on the grid
Definition map.h:726
double get_rot_phi() const
Returns the angle between the x –axis and X –axis.
Definition map.h:775
Coord xa
Absolute x coordinate.
Definition map.h:730
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition map.h:765
Tensor field of valence 0 (or component of a tensorial field).
Definition scalar.h:387
virtual void std_spectral_base()
Sets the spectral bases of the Valeur va to the standard ones for a scalar field.
Definition scalar.C:784
Vector shift_tot
Total shift vector.
Definition star_bhns.h:144
Vector bsn
3-vector shift, divided by N , of the rotating coordinates, .
Definition star_bhns.h:99
virtual void del_deriv() const
Deletes all the derived quantities.
Definition star_bhns.C:341
void kinema_bhns(bool kerrschild, const double &mass_bh, const double &sepa, double omega, double x_rot, double y_rot)
Computes the quantities bsn and pot_centri .
Scalar confo_tot
Total conformal factor.
Definition star_bhns.h:163
Scalar psi4
Fourth power of the total conformal factor.
Definition star_bhns.h:176
Scalar gam0
Lorentz factor between the co-orbiting observer and the Eulerian one.
Definition star_bhns.h:107
Scalar lapconf_tot
Total lapconf function.
Definition star_bhns.h:119
Scalar pot_centri
Centrifugal potential.
Definition star_bhns.h:110
Map & mp
Mapping associated with the star.
Definition star.h:180
Tensor field of valence 1.
Definition vector.h:188
virtual void std_spectral_base()
Sets the standard spectal bases of decomposition for each component.
Definition vector.C:316
Scalar & set(int)
Read/write access to a component.
Definition vector.C:296
Cmp sqrt(const Cmp &)
Square root.
Definition cmp_math.C:220
Cmp log(const Cmp &)
Neperian logarithm.
Definition cmp_math.C:296
void annule_domain(int l)
Sets the Tensor to zero in a given domain.
Definition tensor.C:666
virtual void set_etat_qcq()
Sets the logical state of all components to ETATQCQ (ordinary state).
Definition tensor.C:481
Lorene prototypes.
Definition app_hor.h:64
Standard units of space, time and mass.