LORENE
star_bin_kinema.C
1 /*
2  * Method Star_bin::kinematics
3  *
4  * (see file star.h for documentation)
5  *
6  */
7 
8 /*
9  * Copyright (c) 2004 Francois Limousin
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 as published by
15  * the Free Software Foundation; either version 2 of the License, or
16  * (at your option) any later version.
17  *
18  * LORENE is distributed in the hope that it will be useful,
19  * but WITHOUT ANY WARRANTY; without even the implied warranty of
20  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
21  * GNU General Public License for more details.
22  *
23  * You should have received a copy of the GNU General Public License
24  * along with LORENE; if not, write to the Free Software
25  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
26  *
27  */
28 
29 
30 
31 
32 /*
33  * $Id: star_bin_kinema.C,v 1.11 2016/12/05 16:18:14 j_novak Exp $
34  * $Log: star_bin_kinema.C,v $
35  * Revision 1.11 2016/12/05 16:18:14 j_novak
36  * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
37  *
38  * Revision 1.10 2014/10/13 08:53:38 j_novak
39  * Lorene classes and functions now belong to the namespace Lorene.
40  *
41  * Revision 1.9 2014/10/06 15:13:17 j_novak
42  * Modified #include directives to use c++ syntax.
43  *
44  * Revision 1.8 2006/04/11 14:24:44 f_limousin
45  * New version of the code : improvement of the computation of some
46  * critical sources, estimation of the dirac gauge, helical symmetry...
47  *
48  * Revision 1.7 2005/09/13 19:38:31 f_limousin
49  * Reintroduction of the resolution of the equations in cartesian coordinates.
50  *
51  * Revision 1.6 2005/02/17 17:33:54 f_limousin
52  * Change the name of some quantities to be consistent with other classes
53  * (for instance nnn is changed to nn, shift to beta, beta to lnq...)
54  *
55  * Revision 1.5 2004/06/22 12:51:59 f_limousin
56  * Simplify the computation of gam_pot to improve the convergence of the code.
57  *
58  * Revision 1.4 2004/05/25 14:21:26 f_limousin
59  * New method to compute pot_centri to improve the convergence of the code.
60  *
61  * Revision 1.3 2004/02/27 09:56:10 f_limousin
62  * Correction of an error on the computation of bsn.
63  *
64  * Revision 1.2 2004/01/20 15:18:45 f_limousin
65  * First version
66  *
67  *
68  * $Header: /cvsroot/Lorene/C++/Source/Star/star_bin_kinema.C,v 1.11 2016/12/05 16:18:14 j_novak Exp $
69  *
70  */
71 
72 #include <cmath>
73 
74 // Headers Lorene
75 #include "star.h"
76 
77 namespace Lorene {
78 void Star_bin::kinematics(double omega, double x_axe) {
79 
80  int nz = mp.get_mg()->get_nzone() ;
81  int nzm1 = nz - 1 ;
82 
83  // --------------------
84  // Computation of B^i/N
85  // --------------------
86 
87  // 1/ Computation of omega m^i
88 
89  const Coord& xa = mp.xa ;
90  const Coord& ya = mp.ya ;
91 
93 
94  if (fabs(mp.get_rot_phi()) < 1e-10){
95  bsn.set(1) = - omega * ya ;
96  bsn.set(2) = omega * (xa - x_axe) ;
97  bsn.set(3) = 0 ;
98  }
99  else {
100  bsn.set(1) = omega * ya ;
101  bsn.set(2) = - omega * (xa - x_axe) ;
102  bsn.set(3) = 0 ;
103  }
104 
106 
107  bsn.annule(nzm1, nzm1) ; // set to zero in the ZEC
108 
109  // 2/ Addition of beta and division by lapse
110  // See Eq (47) from Gourgoulhon et al. (2001)
111  // New convention : l = Nn + B ==> B = \beta + \Omega d\phi
112 
113  bsn = ( bsn + beta ) / nn ;
114 
115  bsn.annule(nzm1, nzm1) ; // set to zero in the ZEC
116 
117  //-------------------------
118  // Centrifugal potential
119  //-------------------------
120 
121  // Lorentz factor between the co-orbiting observer and the Eulerian one
122  // See Eq (23) from Gourgoulhon et al. (2001)
123 
124  Sym_tensor gamma_cov (gamma.cov()) ;
125  gamma_cov.change_triad(mp.get_bvect_cart()) ;
126 
127  Scalar gam0 = 1 / sqrt(1 - contract(gamma_cov, 0, 1, bsn * bsn, 0, 1)) ;
128  pot_centri = - log( gam0 ) ;
129 
130  pot_centri.annule(nzm1, nzm1) ; // set to zero in the external domain
131  pot_centri.std_spectral_base() ; // set the bases for spectral expansions
132 
133  // The derived quantities are obsolete
134  // -----------------------------------
135 
136  del_deriv() ;
137 
138 }
139 }
Coord xa
Absolute x coordinate.
Definition: map.h:748
Cmp log(const Cmp &)
Neperian logarithm.
Definition: cmp_math.C:299
Map & mp
Mapping associated with the star.
Definition: star.h:180
Cmp sqrt(const Cmp &)
Square root.
Definition: cmp_math.C:223
virtual void annule(int l_min, int l_max)
Sets the Scalar to zero in several domains.
Definition: scalar.C:397
Metric gamma
3-metric
Definition: star.h:235
Lorene prototypes.
Definition: app_hor.h:67
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:783
Tensor field of valence 0 (or component of a tensorial field).
Definition: scalar.h:393
virtual void std_spectral_base()
Sets the spectral bases of the Valeur va to the standard ones for a scalar field. ...
Definition: scalar.C:790
virtual void change_triad(const Base_vect &)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
double get_rot_phi() const
Returns the angle between the x –axis and X –axis.
Definition: map.h:793
Vector beta
Shift vector.
Definition: star.h:228
virtual void std_spectral_base()
Sets the standard spectal bases of decomposition for each component.
Definition: vector.C:322
Scalar pot_centri
Centrifugal potential.
Definition: star.h:521
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:465
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: star_bin.C:372
Tenseur contract(const Tenseur &, int id1, int id2)
Self contraction of two indices of a Tenseur .
Active physical coordinates and mapping derivatives.
Definition: coord.h:90
virtual const Sym_tensor & cov() const
Read-only access to the covariant representation.
Definition: metric.C:283
Coord ya
Absolute y coordinate.
Definition: map.h:749
const Base_vect_cart & get_bvect_cart() const
Returns the Cartesian basis associated with the coordinates (x,y,z) of the mapping, i.e.
Definition: map.h:809
Scalar nn
Lapse function N .
Definition: star.h:225
void kinematics(double omega, double x_axe)
Computes the quantities bsn and pot_centri.
virtual void change_triad(const Base_vect &new_triad)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
Scalar & set(int)
Read/write access to a component.
Definition: vector.C:302
virtual void annule(int l_min, int l_max)
Sets the Tensor to zero in several domains.
Definition: tensor.C:680
Vector bsn
3-vector shift, divided by N, of the rotating coordinates, .
Definition: star.h:518
Class intended to describe valence-2 symmetric tensors.
Definition: sym_tensor.h:226