LORENE
star_bin_hydro.C
1 /*
2  * Methods of the class Star_bin for computing hydro quantities
3  *
4  * (see file star.h for documentation)
5  */
6 
7 /*
8  * Copyright (c) 2004 Francois Limousin
9  *
10  * This file is part of LORENE.
11  *
12  * LORENE is free software; you can redistribute it and/or modify
13  * it under the terms of the GNU General Public License as published by
14  * the Free Software Foundation; either version 2 of the License, or
15  * (at your option) any later version.
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 char star_bin_hydro_C[] = "$Header: /cvsroot/Lorene/C++/Source/Star/star_bin_hydro.C,v 1.8 2014/10/13 08:53:38 j_novak Exp $" ;
30 
31 /*
32  * $Id: star_bin_hydro.C,v 1.8 2014/10/13 08:53:38 j_novak Exp $
33  * $Log: star_bin_hydro.C,v $
34  * Revision 1.8 2014/10/13 08:53:38 j_novak
35  * Lorene classes and functions now belong to the namespace Lorene.
36  *
37  * Revision 1.7 2005/09/13 19:38:31 f_limousin
38  * Reintroduction of the resolution of the equations in cartesian coordinates.
39  *
40  * Revision 1.6 2005/02/11 18:15:16 f_limousin
41  * Minor modification (to improve the convergence of the code).
42  *
43  * Revision 1.5 2004/06/22 12:50:43 f_limousin
44  * Change qq, qq_auto and qq_comp to beta, beta_auto and beta_comp.
45  *
46  * Revision 1.4 2004/04/08 16:34:39 f_limousin
47  * Changes for irrotationnal binaries.
48  *
49  * Revision 1.3 2004/02/27 09:54:48 f_limousin
50  * Many minor changes.
51  *
52  * Revision 1.2 2004/01/20 15:18:31 f_limousin
53  * First version
54  *
55  *
56  * $Header: /cvsroot/Lorene/C++/Source/Star/star_bin_hydro.C,v 1.8 2014/10/13 08:53:38 j_novak Exp $
57  *
58  */
59 
60 // Headers C
61 
62 // Headers Lorene
63 #include "star.h"
64 
65 namespace Lorene {
67 
68  cout << "loggam 1" << norme(loggam) << endl ;
69 
70  int nz = mp.get_mg()->get_nzone() ;
71  int nzm1 = nz - 1 ;
72 
73  Sym_tensor gamma_cov (gamma.cov()) ;
74  Sym_tensor gamma_con (gamma.con()) ;
75  gamma_cov.change_triad(mp.get_bvect_cart()) ;
76  gamma_con.change_triad(mp.get_bvect_cart()) ;
77 
78  //----------------------------------
79  // Specific relativistic enthalpy ---> hhh
80  //----------------------------------
81 
82  Scalar hhh = exp(ent) ; // = 1 at the Newtonian limit
83  hhh.std_spectral_base() ;
84 
85  //---------------------------------------------------
86  // Lorentz factor between the co-orbiting ---> gam00
87  // observer and the Eulerian one
88  // See Eq (23) and (24) from Gourgoulhon et al. (2001)
89  //---------------------------------------------------
90 
91  Scalar gam0 = 1 / sqrt( 1 - contract(gamma_cov, 0, 1, bsn * bsn, 0, 1)) ;
92  gam0.std_spectral_base() ;
93 
94  //------------------------------------------
95  // Lorentz factor and 3-velocity of the fluid
96  // with respect to the Eulerian observer
97  //------------------------------------------
98 
99  if (irrotational) {
100 
102 
103  // See Eq (32) from Gourgoulhon et al. (2001)
104  gam_euler = sqrt( 1 + contract(gamma_con, 0, 1, d_psi * d_psi, 0, 1)
105  / (hhh%hhh) ) ;
106 
108 
109  u_euler = contract(gamma_con, 0, d_psi, 0)/( hhh % gam_euler ) ;
111 
112  }
113  else {
114  // Rigid rotation
115  // --------------
116 
117  gam_euler = gam0 ;
119  u_euler = bsn ;
120 
121  }
122 
123  //------------------------------------
124  // Energy density E with respect to the Eulerian observer
125  // See Eq (53) from Gourgoulhon et al. (2001)
126  //------------------------------------
127 
128  ener_euler = gam_euler % gam_euler % ( ener + press ) - press ;
129 
130  //------------------------------------
131  // Trace of the stress tensor with respect to the Eulerian observer
132  // See Eq (54) from Gourgoulhon et al. (2001)
133  //------------------------------------
134 
135  s_euler = 3 * press + ( ener_euler + press ) %
136  contract(gamma_cov, 0, 1, u_euler * u_euler, 0 ,1) ;
137 
138  //-------------------------------------------
139  // Spatial part of the stress-energy tensor with respect
140  // to the Eulerian observer.
141  //-------------------------------------------
142 
143  for(int i=1; i<=3; i++){
144  for(int j=1; j<=3; j++){
145  stress_euler.set(i,j) = (ener_euler + press )*u_euler(i)
146  *u_euler(j) + press * gamma_con(i,j) ;
147  }
148  }
149 
150  //-------------------------------------------
151  // Lorentz factor between the fluid and ---> gam
152  // co-orbiting observers
153  // See Eq (58) from Gourgoulhon et al. (2001)
154  //--------------------------------------------
155 
156  if (irrotational) {
157 
158  Scalar tmp = ( 1 - contract(gamma_cov, 0, 1, bsn * u_euler, 0, 1) ) ;
159  tmp.std_spectral_base() ;
160  Scalar gam = gam0 % gam_euler % tmp ;
161 
162  //-------------------------------------------
163  // Spatial projection of the fluid 3-velocity
164  // with respect to the co-orbiting observer
165  //--------------------------------------------
166 
167  wit_w = - gam_euler / gam * u_euler + gam0 * bsn ;
168 
169  wit_w.std_spectral_base() ; // set the bases for spectral expansions
170 
171  wit_w.annule_domain(nzm1) ; // zero in the ZEC
172 
173 
174  //-------------------------------------------
175  // Logarithm of the Lorentz factor between
176  // the fluid and co-orbiting observers
177  //--------------------------------------------
178 
179  loggam = log( gam ) ;
180 
181  loggam.std_spectral_base() ; // set the bases for spectral expansions
182 
183  //-------------------------------------------------
184  // Velocity fields set to zero in external domains
185  //-------------------------------------------------
186 
187  loggam.annule_domain(nzm1) ; // zero in the ZEC only
188 
189  wit_w.annule_domain(nzm1) ; // zero outside the star
190 
191  u_euler.annule_domain(nzm1) ; // zero outside the star
192 
193  loggam.set_dzpuis(0) ;
194  }
195  else {
196 
197  loggam = 0 ;
198  wit_w.set_etat_zero() ;
199  }
200 
201  // The derived quantities are obsolete
202  // -----------------------------------
203 
204  del_deriv() ;
205 
206 }
207 }
void annule_domain(int l)
Sets the Tensor to zero in a given domain.
Definition: tensor.C:666
bool irrotational
true for an irrotational star, false for a corotating one
Definition: star.h:491
Cmp log(const Cmp &)
Neperian logarithm.
Definition: cmp_math.C:296
virtual const Sym_tensor & con() const
Read-only access to the contravariant representation.
Definition: metric.C:290
Cmp exp(const Cmp &)
Exponential.
Definition: cmp_math.C:270
Vector wit_w
Spatial projection of the fluid 3-velocity with respect to the co-orbiting observer.
Definition: star.h:507
Map & mp
Mapping associated with the star.
Definition: star.h:180
Cmp sqrt(const Cmp &)
Square root.
Definition: cmp_math.C:220
Metric gamma
3-metric
Definition: star.h:235
Vector d_psi
Gradient of (in the irrotational case) (Spherical components with respect to the mapping of the star...
Definition: star.h:501
Lorene prototypes.
Definition: app_hor.h:64
Sym_tensor stress_euler
Spatial part of the stress-energy tensor with respect to the Eulerian observer.
Definition: star.h:212
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
Scalar loggam
Logarithm of the Lorentz factor between the fluid and the co-orbiting observer.
Definition: star.h:512
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
Scalar ent
Log-enthalpy.
Definition: star.h:190
Scalar s_euler
Trace of the stress scalar in the Eulerian frame.
Definition: star.h:201
void set_dzpuis(int)
Modifies the dzpuis flag.
Definition: scalar.C:808
virtual void std_spectral_base()
Sets the standard spectal bases of decomposition for each component.
Definition: vector.C:316
Scalar gam_euler
Lorentz factor between the fluid and Eulerian observers.
Definition: star.h:204
Scalar ener
Total energy density in the fluid frame.
Definition: star.h:193
Tbl norme(const Cmp &)
Sums of the absolute values of all the values of the Cmp in each domain.
Definition: cmp_math.C:481
Scalar press
Fluid pressure.
Definition: star.h:194
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:448
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: star_bin.C:369
Vector u_euler
Fluid 3-velocity with respect to the Eulerian observer.
Definition: star.h:207
Tenseur contract(const Tenseur &, int id1, int id2)
Self contraction of two indices of a Tenseur .
virtual void hydro_euler()
Computes the hydrodynamical quantities relative to the Eulerian observer from those in the fluid fram...
virtual const Sym_tensor & cov() const
Read-only access to the covariant representation.
Definition: metric.C:280
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:791
Scalar & set(const Itbl &ind)
Returns the value of a component (read/write version).
Definition: tensor.C:654
virtual void change_triad(const Base_vect &new_triad)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
virtual void set_etat_zero()
Sets the logical state of all components to ETATZERO (zero state).
Definition: tensor.C:497
Vector bsn
3-vector shift, divided by N, of the rotating coordinates, .
Definition: star.h:518
Scalar ener_euler
Total energy density in the Eulerian frame.
Definition: star.h:198
Class intended to describe valence-2 symmetric tensors.
Definition: sym_tensor.h:223