LORENE
et_bin_kinema.C
1 /*
2  * Method Etoile_bin::kinematics
3  *
4  * (see file etoile.h for documentation)
5  *
6  */
7 
8 /*
9  * Copyright (c) 2000-2001 Eric Gourgoulhon
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: et_bin_kinema.C,v 1.7 2016/12/05 16:17:53 j_novak Exp $
34  * $Log: et_bin_kinema.C,v $
35  * Revision 1.7 2016/12/05 16:17:53 j_novak
36  * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
37  *
38  * Revision 1.6 2014/10/13 08:52:56 j_novak
39  * Lorene classes and functions now belong to the namespace Lorene.
40  *
41  * Revision 1.5 2005/08/29 15:10:16 p_grandclement
42  * Addition of things needed :
43  * 1) For BBH with different masses
44  * 2) Provisory files for the mixted binaries (Bh and NS) : THIS IS NOT
45  * WORKING YET !!!
46  *
47  * Revision 1.4 2003/03/03 19:18:12 f_limousin
48  * Suppression of bsn.set_triad(ref_triad).
49  *
50  * Revision 1.3 2003/01/17 13:35:48 f_limousin
51  * Add comments and replace A^2*flat_scalar_prod by sprod
52  *
53  * Revision 1.2 2002/12/10 15:56:43 k_taniguchi
54  * Change the multiplication "*" to "%".
55  *
56  * Revision 1.1.1.1 2001/11/20 15:19:28 e_gourgoulhon
57  * LORENE
58  *
59  * Revision 2.4 2000/02/17 18:51:22 eric
60  * pot_centri is set to zero only in the external compactified domain.
61  *
62  * Revision 2.3 2000/02/12 18:37:17 eric
63  * Appel de set_std_base sur les quantites calculees.
64  *
65  * Revision 2.2 2000/02/08 19:29:27 eric
66  * Calcul sur les tenseurs.
67  *
68  * Revision 2.1 2000/02/04 17:14:27 eric
69  * *** empty log message ***
70  *
71  * Revision 2.0 2000/02/04 16:37:51 eric
72  * *** empty log message ***
73  *
74  *
75  * $Header: /cvsroot/Lorene/C++/Source/Etoile/et_bin_kinema.C,v 1.7 2016/12/05 16:17:53 j_novak Exp $
76  *
77  */
78 
79 // Headers Lorene
80 #include "etoile.h"
81 #include "graphique.h"
82 
83 namespace Lorene {
84 void Etoile_bin::kinematics(double omega, double x_axe) {
85 
86  int nz = mp.get_mg()->get_nzone() ;
87  int nzm1 = nz - 1 ;
88 
89  // --------------------
90  // Computation of B^i/N
91  // --------------------
92 
93  // 1/ Computation of - omega m^i
94 
95  const Coord& xa = mp.xa ;
96  const Coord& ya = mp.ya ;
97 
98  bsn.set_etat_qcq() ;
99 
100  bsn.set(0) = omega * ya ;
101  bsn.set(1) = - omega * (xa - x_axe) ;
102  bsn.set(2) = 0 ;
103 
104  bsn.annule(nzm1, nzm1) ; // set to zero in the ZEC
105 
106  // 2/ Addition of shift and division by lapse
107  // See Eq (47) from Gourgoulhon et al. (2001)
108 
109  bsn = ( bsn + shift ) / nnn ;
110 
111  bsn.annule(nzm1, nzm1) ; // set to zero in the ZEC
112  bsn.set_std_base() ; // set the bases for spectral expansions
113 
114  //-------------------------
115  // Centrifugal potentatial
116  //-------------------------
117 
118  if (relativistic) {
119 
120  // Lorentz factor between the co-orbiting observer and the Eulerian one
121  // See Eq (23) from Gourgoulhon et al. (2001)
122  Tenseur gam0 = 1 / sqrt( 1-sprod(bsn, bsn) ) ;
123 
124  pot_centri = - log( gam0 ) ;
125 
126  }
127  else {
128 
130 
131 
132  // See Eq (40) from Gourgoulhon et al. (2001)
133  pot_centri.set() = - 0.5 * omega * omega * (
134  (xa - x_axe) * (xa - x_axe) + ya * ya ) ;
135 
136  }
137 
138  pot_centri.annule(nzm1, nzm1) ; // set to zero in the external domain
139  pot_centri.set_std_base() ; // set the bases for spectral expansions
140 
141  // The derived quantities are obsolete
142  // -----------------------------------
143 
144  del_deriv() ;
145 
146 
147 }
148 }
Coord xa
Absolute x coordinate.
Definition: map.h:748
Cmp log(const Cmp &)
Neperian logarithm.
Definition: cmp_math.C:299
void annule(int l)
Sets the Tenseur to zero in a given domain.
Definition: tenseur.C:916
Cmp sqrt(const Cmp &)
Square root.
Definition: cmp_math.C:223
void set_std_base()
Set the standard spectal basis of decomposition for each component.
Definition: tenseur.C:1186
Tenseur pot_centri
Centrifugal potential.
Definition: etoile.h:956
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
Tenseur nnn
Total lapse function.
Definition: etoile.h:512
Tenseur shift
Total shift vector.
Definition: etoile.h:515
Cmp & set()
Read/write for a scalar (see also operator=(const Cmp&) ).
Definition: tenseur.C:840
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: etoile_bin.C:450
Tenseur bsn
3-vector shift, divided by N , of the rotating coordinates, .
Definition: etoile.h:953
Map & mp
Mapping associated with the star.
Definition: etoile.h:432
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:465
Active physical coordinates and mapping derivatives.
Definition: coord.h:90
Coord ya
Absolute y coordinate.
Definition: map.h:749
bool relativistic
Indicator of relativity: true for a relativistic star, false for a Newtonian one. ...
Definition: etoile.h:440
virtual void kinematics(double omega, double x_axe)
Computes the quantities bsn and pot_centri .
Definition: et_bin_kinema.C:84
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: tenseur.C:652
virtual Tenseur sprod(const Tenseur &t1, const Tenseur &t2) const
Performs the scalar product of two tensors by contracting the last index of t1 with the first index o...
Definition: etoile_bin.C:751
Tensor handling *** DEPRECATED : use class Tensor instead ***.
Definition: tenseur.h:304