LORENE
et_rot_hydro.C
1/*
2 * Function Etoile_rot::hydro_euler
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_rot_hydro_C[] = "$Header: /cvsroot/Lorene/C++/Source/Etoile/et_rot_hydro.C,v 1.4 2014/10/13 08:52:58 j_novak Exp $" ;
31
32/*
33 * $Id: et_rot_hydro.C,v 1.4 2014/10/13 08:52:58 j_novak Exp $
34 * $Log: et_rot_hydro.C,v $
35 * Revision 1.4 2014/10/13 08:52:58 j_novak
36 * Lorene classes and functions now belong to the namespace Lorene.
37 *
38 * Revision 1.3 2014/10/06 15:13:09 j_novak
39 * Modified #include directives to use c++ syntax.
40 *
41 * Revision 1.2 2006/09/14 07:37:47 j_novak
42 * Removal of a test on u_euler.
43 *
44 * Revision 1.1.1.1 2001/11/20 15:19:28 e_gourgoulhon
45 * LORENE
46 *
47 * Revision 2.6 2000/10/12 15:36:06 eric
48 * Le test sur u_euler est desormais effectue sur la totalite de u_euler
49 * (ie comprenant le shift).
50 *
51 * Revision 2.5 2000/10/06 15:08:13 eric
52 * Calcul 3-D de u_euler.
53 * uuu s'en deduit.
54 *
55 * Revision 2.4 2000/08/31 15:38:34 eric
56 * Appel du nouvel operateur Cmp::mult_rsint pour le calcul de uuu.
57 *
58 * Revision 2.3 2000/08/25 12:28:21 eric
59 * *** empty log message ***
60 *
61 * Revision 2.2 2000/08/18 14:01:49 eric
62 * *** empty log message ***
63 *
64 * Revision 2.1 2000/08/17 12:39:56 eric
65 * *** empty log message ***
66 *
67 * Revision 2.0 2000/07/21 16:31:00 eric
68 * *** empty log message ***
69 *
70 *
71 * $Header: /cvsroot/Lorene/C++/Source/Etoile/et_rot_hydro.C,v 1.4 2014/10/13 08:52:58 j_novak Exp $
72 *
73 */
74
75// Headers C
76#include <cstdlib>
77
78// Headers Lorene
79#include "etoile.h"
80#include "utilitaires.h"
81
82namespace Lorene {
84
85 int nz = mp.get_mg()->get_nzone() ;
86 int nzm1 = nz - 1 ;
87
88 // Computation of u_euler
89 // ----------------------
90
91 const Coord& x = mp.x ;
92 const Coord& y = mp.y ;
93
95
96 u_euler.set(0) = - omega * y ; // Cartesian components of solid rotation
97 u_euler.set(1) = omega * x ;
98 u_euler.set(2) = 0 ;
99 u_euler.annule(nzm1) ;
100
101 u_euler.set_triad( mp.get_bvect_cart() ) ; // Triad = Cartesian triad
102
103 u_euler.set_std_base() ; // sets the standard bases for spectral expansions
104
105 u_euler = ( u_euler - shift ) / nnn ;
106
107 u_euler.set_std_base() ; // sets the standard bases for spectral expansions
108
109 if ( (u_euler(0).get_etat() == ETATZERO) &&
110 (u_euler(1).get_etat() == ETATZERO) &&
111 (u_euler(2).get_etat() == ETATZERO) ) {
112
113 u_euler = 0 ;
114 }
115
116
117 // Computation of uuu (norme of u_euler)
118 // ------------------
119
120 // The scalar product is performed on the spherical components:
121
122 Tenseur us = u_euler ;
124
125 Cmp uuu2 = a_car() * ( us(0) * us(0) + us(1) * us(1) )
126 + b_car() * us(2) * us(2) ;
127
128 uuu = sqrt( uuu2 ) ;
129
130 if (uuu.get_etat() == ETATQCQ) {
131 ((uuu.set()).va).set_base( us(2).va.base ) ; // Same basis as
132 } // (Omega -N^phi) r sin(theta)
133
134
135 // Lorentz factor
136 // --------------
137
138 Tenseur u2(mp) ;
139 u2 = unsurc2 * uuu2 ;
140
141 Tenseur gam2 = 1 / (1 - u2) ;
142
143 gam_euler = sqrt(gam2) ;
144
145 gam_euler.set_std_base() ; // sets the standard spectral bases for
146 // a scalar field
147
148 // Energy density E with respect to the Eulerian observer
149 //------------------------------------
150
151 ener_euler = gam2 * ( ener + press ) - press ;
152
154
155 // Trace of the stress tensor with respect to the Eulerian observer
156 //------------------------------------
157
158 s_euler = 3 * press + ( ener_euler + press ) * u2 ;
159
161
162 // The derived quantities are obsolete
163 // -----------------------------------
164
165 del_deriv() ;
166
167
168}
169}
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition cmp.h:446
Active physical coordinates and mapping derivatives.
Definition coord.h:90
Tenseur uuu
Norm of u_euler.
Definition etoile.h:1518
virtual void del_deriv() const
Deletes all the derived quantities.
Definition etoile_rot.C:335
double omega
Rotation angular velocity ([f_unit] )
Definition etoile.h:1501
virtual void hydro_euler()
Computes the hydrodynamical quantities relative to the Eulerian observer from those in the fluid fram...
Tenseur b_car
Square of the metric factor B.
Definition etoile.h:1507
Tenseur nnn
Total lapse function.
Definition etoile.h:509
Tenseur u_euler
Fluid 3-velocity with respect to the Eulerian observer.
Definition etoile.h:474
Tenseur gam_euler
Lorentz factor between the fluid and Eulerian observers.
Definition etoile.h:471
Map & mp
Mapping associated with the star.
Definition etoile.h:429
Tenseur ener
Total energy density in the fluid frame.
Definition etoile.h:460
Tenseur press
Fluid pressure.
Definition etoile.h:461
Tenseur ener_euler
Total energy density in the Eulerian frame.
Definition etoile.h:465
Tenseur shift
Total shift vector.
Definition etoile.h:512
Tenseur s_euler
Trace of the stress tensor in the Eulerian frame.
Definition etoile.h:468
Tenseur a_car
Total conformal factor .
Definition etoile.h:515
double unsurc2
: unsurc2=1 for a relativistic star, 0 for a Newtonian one.
Definition etoile.h:442
const Base_vect_cart & get_bvect_cart() const
Returns the Cartesian basis associated with the coordinates (x,y,z) of the mapping,...
Definition map.h:791
Coord y
y coordinate centered on the grid
Definition map.h:727
const Base_vect_spher & get_bvect_spher() const
Returns the orthonormal vectorial basis associated with the coordinates of the mapping.
Definition map.h:783
Coord x
x coordinate centered on the grid
Definition map.h:726
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
void change_triad(const Base_vect &new_triad)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
Definition tenseur.C:668
void set_triad(const Base_vect &new_triad)
Assigns a new vectorial basis (triad) of decomposition.
Definition tenseur.C:674
int get_etat() const
Returns the logical state.
Definition tenseur.h:707
Cmp sqrt(const Cmp &)
Square root.
Definition cmp_math.C:220
Lorene prototypes.
Definition app_hor.h:64