LORENE
et_rot_diff_hydro.C
1 /*
2  * Function Et_rot_diff::hydro_euler
3  *
4  * (see file et_rot_diff.h for documentation)
5  *
6  */
7 
8 /*
9  * Copyright (c) 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_rot_diff_hydro.C,v 1.4 2016/12/05 16:17:54 j_novak Exp $
34  * $Log: et_rot_diff_hydro.C,v $
35  * Revision 1.4 2016/12/05 16:17:54 j_novak
36  * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
37  *
38  * Revision 1.3 2014/10/13 08:52:57 j_novak
39  * Lorene classes and functions now belong to the namespace Lorene.
40  *
41  * Revision 1.2 2014/10/06 15:13:09 j_novak
42  * Modified #include directives to use c++ syntax.
43  *
44  * Revision 1.1.1.1 2001/11/20 15:19:28 e_gourgoulhon
45  * LORENE
46  *
47  * Revision 1.1 2001/10/19 08:18:36 eric
48  * Initial revision
49  *
50  *
51  * $Header: /cvsroot/Lorene/C++/Source/Etoile/et_rot_diff_hydro.C,v 1.4 2016/12/05 16:17:54 j_novak Exp $
52  *
53  */
54 
55 
56 // Headers C
57 #include <cstdlib>
58 
59 // Headers Lorene
60 #include "et_rot_diff.h"
61 #include "utilitaires.h"
62 
63 namespace Lorene {
65 
66  int nz = mp.get_mg()->get_nzone() ;
67  int nzm1 = nz - 1 ;
68 
69  // Computation of u_euler
70  // ----------------------
71 
72  Cmp x(mp) ;
73  Cmp y(mp) ;
74  x = mp.x ;
75  y = mp.y ;
76 
78 
79  // Cartesian components of differential rotation:
80 
81  u_euler.set(0) = - omega_field() * y ;
82  u_euler.set(1) = omega_field() * x ;
83  u_euler.set(2) = 0 ;
84  u_euler.annule(nzm1) ;
85 
86  u_euler.set_triad( mp.get_bvect_cart() ) ; // Triad = Cartesian triad
87 
88  u_euler.set_std_base() ; // sets the standard bases for spectral expansions
89 
90  u_euler = ( u_euler - shift ) / nnn ;
91 
92  u_euler.set_std_base() ; // sets the standard bases for spectral expansions
93 
94 //## Test
95  Tenseur utest(mp, 1, CON, mp.get_bvect_spher()) ;
96  utest.set_etat_qcq() ;
97 
98  utest.set(0) = 0 ; // Spherical components of differential rotation
99  utest.set(1) = 0 ;
100  utest.set(2) = ( omega_field() - nphi() ) / nnn();
101 
102  utest.set(2).annule(nzm1) ;
103  utest.set(2).std_base_scal() ;
104  utest.set(2).mult_rsint() ; // Multiplication by r sin(theta)
105 
106  utest.set_triad( mp.get_bvect_spher() ) ;
107 
108  utest.change_triad( mp.get_bvect_cart() ) ;
109 
110  for (int i=0; i<3; i++) {
111  Valeur& uu = u_euler.set(i).va ;
112  Valeur& ut = utest.set(i).va ;
113 
114  if (uu.get_etat() != ETATZERO) {
115  uu.coef() ;
116 
117  if (ut.get_etat() == ETATZERO) {
118  ut.set_etat_cf_qcq() ;
119  *(ut.c_cf) = 0 ;
120  ut.c_cf->base = uu.c_cf->base ;
121  }
122  else {
123  ut.coef() ;
124  }
125 
126  Mtbl_cf diff = *(uu.c_cf) - *(ut.c_cf) ;
127  cout << "Et_rot_diff::hydro_euler: test u_euler(" << i << ") : "
128  << max( abs(diff) )(0) << endl ;
129 
130  }
131  }
132 //##
133 
134  if ( (u_euler(0).get_etat() == ETATZERO) &&
135  (u_euler(1).get_etat() == ETATZERO) &&
136  (u_euler(2).get_etat() == ETATZERO) ) {
137 
138  u_euler = 0 ;
139  }
140 
141 
142  // Computation of uuu (norme of u_euler)
143  // ------------------
144 
145  // The scalar product is performed on the spherical components:
146 
147  Tenseur us = u_euler ;
148  us.change_triad( mp.get_bvect_spher() ) ;
149 
150  Cmp uuu2 = a_car() * ( us(0) * us(0) + us(1) * us(1) )
151  + b_car() * us(2) * us(2) ;
152 
153  uuu = sqrt( uuu2 ) ;
154 
155  if (uuu.get_etat() == ETATQCQ) {
156  ((uuu.set()).va).set_base( us(2).va.base ) ; // Same basis as
157  } // (Omega -N^phi) r sin(theta)
158 
159 
160  // Lorentz factor
161  // --------------
162 
163  Tenseur u2(mp) ;
164  u2 = unsurc2 * uuu2 ;
165 
166  Tenseur gam2 = 1 / (1 - u2) ;
167 
168  gam_euler = sqrt(gam2) ;
169 
170  gam_euler.set_std_base() ; // sets the standard spectral bases for
171  // a scalar field
172 
173  // Energy density E with respect to the Eulerian observer
174  //------------------------------------
175 
176  ener_euler = gam2 * ( ener + press ) - press ;
177 
179 
180  // Trace of the stress tensor with respect to the Eulerian observer
181  //------------------------------------
182 
183  s_euler = 3 * press + ( ener_euler + press ) * u2 ;
184 
185  s_euler.set_std_base() ;
186 
187  // The derived quantities are obsolete
188  // -----------------------------------
189 
190  del_deriv() ;
191 
192 
193 }
194 }
void annule(int l)
Sets the Tenseur to zero in a given domain.
Definition: tenseur.C:916
Mtbl_cf * c_cf
Coefficients of the spectral expansion of the function.
Definition: valeur.h:312
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition: cmp.h:446
void set_triad(const Base_vect &new_triad)
Assigns a new vectorial basis (triad) of decomposition.
Definition: tenseur.C:690
void set_etat_cf_qcq()
Sets the logical state to ETATQCQ (ordinary state) for values in the configuration space (Mtbl_cf c_c...
Definition: valeur.C:715
void coef() const
Computes the coeffcients of *this.
Definition: valeur_coef.C:151
Cmp sqrt(const Cmp &)
Square root.
Definition: cmp_math.C:223
void annule(int l)
Sets the Cmp to zero in a given domain.
Definition: cmp.C:351
void set_std_base()
Set the standard spectal basis of decomposition for each component.
Definition: tenseur.C:1186
const Base_vect_spher & get_bvect_spher() const
Returns the orthonormal vectorial basis associated with the coordinates of the mapping.
Definition: map.h:801
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
virtual void hydro_euler()
Computes the hydrodynamical quantities relative to the Eulerian observer from those in the fluid fram...
Tenseur nnn
Total lapse function.
Definition: etoile.h:512
Tenseur nphi
Metric coefficient .
Definition: etoile.h:1513
Tenseur s_euler
Trace of the stress tensor in the Eulerian frame.
Definition: etoile.h:471
Tenseur b_car
Square of the metric factor B.
Definition: etoile.h:1510
double unsurc2
: unsurc2=1 for a relativistic star, 0 for a Newtonian one.
Definition: etoile.h:445
Values and coefficients of a (real-value) function.
Definition: valeur.h:297
Tenseur press
Fluid pressure.
Definition: etoile.h:464
Tenseur shift
Total shift vector.
Definition: etoile.h:515
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: etoile_rot.C:344
int get_etat() const
Returns the logical state.
Definition: valeur.h:760
Cmp & set()
Read/write for a scalar (see also operator=(const Cmp&) ).
Definition: tenseur.C:840
Tenseur u_euler
Fluid 3-velocity with respect to the Eulerian observer.
Definition: etoile.h:477
void change_triad(const Base_vect &new_triad)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
Definition: tenseur.C:684
Tenseur gam_euler
Lorentz factor between the fluid and Eulerian observers.
Definition: etoile.h:474
Map & mp
Mapping associated with the star.
Definition: etoile.h:432
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:465
Tbl max(const Cmp &)
Maximum values of a Cmp in each domain.
Definition: cmp_math.C:438
int get_etat() const
Returns the logical state.
Definition: tenseur.h:710
Tenseur uuu
Norm of u_euler.
Definition: etoile.h:1521
Tenseur a_car
Total conformal factor .
Definition: etoile.h:518
Tenseur ener
Total energy density in the fluid frame.
Definition: etoile.h:463
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
Coord y
y coordinate centered on the grid
Definition: map.h:745
Cmp abs(const Cmp &)
Absolute value.
Definition: cmp_math.C:413
Coefficients storage for the multi-domain spectral method.
Definition: mtbl_cf.h:196
Coord x
x coordinate centered on the grid
Definition: map.h:744
Base_val base
Bases of the spectral expansions.
Definition: mtbl_cf.h:210
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: tenseur.C:652
Tenseur ener_euler
Total energy density in the Eulerian frame.
Definition: etoile.h:468
Valeur va
The numerical value of the Cmp.
Definition: cmp.h:464
Tensor handling *** DEPRECATED : use class Tensor instead ***.
Definition: tenseur.h:304
Tenseur omega_field
Field .
Definition: et_rot_diff.h:108
Tbl & set(int l)
Read/write of the value in a given domain (configuration space).
Definition: valeur.h:373