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
30char et_rot_diff_hydro_C[] = "$Header: /cvsroot/Lorene/C++/Source/Etoile/et_rot_diff_hydro.C,v 1.3 2014/10/13 08:52:57 j_novak Exp $" ;
31
32/*
33 * $Id: et_rot_diff_hydro.C,v 1.3 2014/10/13 08:52:57 j_novak Exp $
34 * $Log: et_rot_diff_hydro.C,v $
35 * Revision 1.3 2014/10/13 08:52:57 j_novak
36 * Lorene classes and functions now belong to the namespace Lorene.
37 *
38 * Revision 1.2 2014/10/06 15:13:09 j_novak
39 * Modified #include directives to use c++ syntax.
40 *
41 * Revision 1.1.1.1 2001/11/20 15:19:28 e_gourgoulhon
42 * LORENE
43 *
44 * Revision 1.1 2001/10/19 08:18:36 eric
45 * Initial revision
46 *
47 *
48 * $Header: /cvsroot/Lorene/C++/Source/Etoile/et_rot_diff_hydro.C,v 1.3 2014/10/13 08:52:57 j_novak Exp $
49 *
50 */
51
52
53// Headers C
54#include <cstdlib>
55
56// Headers Lorene
57#include "et_rot_diff.h"
58#include "utilitaires.h"
59
60namespace Lorene {
62
63 int nz = mp.get_mg()->get_nzone() ;
64 int nzm1 = nz - 1 ;
65
66 // Computation of u_euler
67 // ----------------------
68
69 Cmp x(mp) ;
70 Cmp y(mp) ;
71 x = mp.x ;
72 y = mp.y ;
73
75
76 // Cartesian components of differential rotation:
77
78 u_euler.set(0) = - omega_field() * y ;
79 u_euler.set(1) = omega_field() * x ;
80 u_euler.set(2) = 0 ;
81 u_euler.annule(nzm1) ;
82
83 u_euler.set_triad( mp.get_bvect_cart() ) ; // Triad = Cartesian triad
84
85 u_euler.set_std_base() ; // sets the standard bases for spectral expansions
86
87 u_euler = ( u_euler - shift ) / nnn ;
88
89 u_euler.set_std_base() ; // sets the standard bases for spectral expansions
90
91//## Test
92 Tenseur utest(mp, 1, CON, mp.get_bvect_spher()) ;
93 utest.set_etat_qcq() ;
94
95 utest.set(0) = 0 ; // Spherical components of differential rotation
96 utest.set(1) = 0 ;
97 utest.set(2) = ( omega_field() - nphi() ) / nnn();
98
99 utest.set(2).annule(nzm1) ;
100 utest.set(2).std_base_scal() ;
101 utest.set(2).mult_rsint() ; // Multiplication by r sin(theta)
102
103 utest.set_triad( mp.get_bvect_spher() ) ;
104
105 utest.change_triad( mp.get_bvect_cart() ) ;
106
107 for (int i=0; i<3; i++) {
108 Valeur& uu = u_euler.set(i).va ;
109 Valeur& ut = utest.set(i).va ;
110
111 if (uu.get_etat() != ETATZERO) {
112 uu.coef() ;
113
114 if (ut.get_etat() == ETATZERO) {
115 ut.set_etat_cf_qcq() ;
116 *(ut.c_cf) = 0 ;
117 ut.c_cf->base = uu.c_cf->base ;
118 }
119 else {
120 ut.coef() ;
121 }
122
123 Mtbl_cf diff = *(uu.c_cf) - *(ut.c_cf) ;
124 cout << "Et_rot_diff::hydro_euler: test u_euler(" << i << ") : "
125 << max( abs(diff) )(0) << endl ;
126
127 }
128 }
129//##
130
131 if ( (u_euler(0).get_etat() == ETATZERO) &&
132 (u_euler(1).get_etat() == ETATZERO) &&
133 (u_euler(2).get_etat() == ETATZERO) ) {
134
135 u_euler = 0 ;
136 }
137
138
139 // Computation of uuu (norme of u_euler)
140 // ------------------
141
142 // The scalar product is performed on the spherical components:
143
144 Tenseur us = u_euler ;
146
147 Cmp uuu2 = a_car() * ( us(0) * us(0) + us(1) * us(1) )
148 + b_car() * us(2) * us(2) ;
149
150 uuu = sqrt( uuu2 ) ;
151
152 if (uuu.get_etat() == ETATQCQ) {
153 ((uuu.set()).va).set_base( us(2).va.base ) ; // Same basis as
154 } // (Omega -N^phi) r sin(theta)
155
156
157 // Lorentz factor
158 // --------------
159
160 Tenseur u2(mp) ;
161 u2 = unsurc2 * uuu2 ;
162
163 Tenseur gam2 = 1 / (1 - u2) ;
164
165 gam_euler = sqrt(gam2) ;
166
167 gam_euler.set_std_base() ; // sets the standard spectral bases for
168 // a scalar field
169
170 // Energy density E with respect to the Eulerian observer
171 //------------------------------------
172
173 ener_euler = gam2 * ( ener + press ) - press ;
174
176
177 // Trace of the stress tensor with respect to the Eulerian observer
178 //------------------------------------
179
180 s_euler = 3 * press + ( ener_euler + press ) * u2 ;
181
183
184 // The derived quantities are obsolete
185 // -----------------------------------
186
187 del_deriv() ;
188
189
190}
191}
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition cmp.h:446
void mult_rsint()
Multiplication by .
Valeur va
The numerical value of the Cmp
Definition cmp.h:464
void std_base_scal()
Sets the spectral bases of the Valeur va to the standard ones for a scalar.
Definition cmp.C:644
void annule(int l)
Sets the Cmp to zero in a given domain.
Definition cmp.C:348
virtual void hydro_euler()
Computes the hydrodynamical quantities relative to the Eulerian observer from those in the fluid fram...
Tenseur omega_field
Field .
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
Tenseur nphi
Metric coefficient .
Definition etoile.h:1510
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
Coefficients storage for the multi-domain spectral method.
Definition mtbl_cf.h:186
Base_val base
Bases of the spectral expansions.
Definition mtbl_cf.h:200
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
Values and coefficients of a (real-value) function.
Definition valeur.h:287
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:712
int get_etat() const
Returns the logical state.
Definition valeur.h:726
Mtbl_cf * c_cf
Coefficients of the spectral expansion of the function.
Definition valeur.h:302
void coef() const
Computes the coeffcients of *this.
Cmp sqrt(const Cmp &)
Square root.
Definition cmp_math.C:220
Tbl max(const Cmp &)
Maximum values of a Cmp in each domain.
Definition cmp_math.C:435
Cmp abs(const Cmp &)
Absolute value.
Definition cmp_math.C:410
Lorene prototypes.
Definition app_hor.h:64