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