LORENE
star_rot_hydro.C
1/*
2 * Method Star_rot::hydro_euler
3 *
4 * (see file star_rot.h for documentation)
5 *
6 */
7
8/*
9 * Copyright (c) 2010 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 star_rot_hydro_C[] = "$Header: /cvsroot/Lorene/C++/Source/Star/star_rot_hydro.C,v 1.3 2014/10/13 08:53:39 j_novak Exp $" ;
31
32/*
33 * $Id: star_rot_hydro.C,v 1.3 2014/10/13 08:53:39 j_novak Exp $
34 * $Log: star_rot_hydro.C,v $
35 * Revision 1.3 2014/10/13 08:53:39 j_novak
36 * Lorene classes and functions now belong to the namespace Lorene.
37 *
38 * Revision 1.2 2014/10/06 15:13:17 j_novak
39 * Modified #include directives to use c++ syntax.
40 *
41 * Revision 1.1 2010/01/25 18:15:52 e_gourgoulhon
42 * First version.
43 *
44 *
45 * $Header: /cvsroot/Lorene/C++/Source/Star/star_rot_hydro.C,v 1.3 2014/10/13 08:53:39 j_novak Exp $
46 *
47 */
48
49// Headers C
50#include <cstdlib>
51
52// Headers Lorene
53#include "star_rot.h"
54#include "utilitaires.h"
55
56namespace Lorene {
58
59 int nz = mp.get_mg()->get_nzone() ;
60 int nzm1 = nz - 1 ;
61
62 // Computation of u_euler
63 // ----------------------
64
65 const Coord& x = mp.x ;
66 const Coord& y = mp.y ;
67
69
70 u_euler.set(1) = - omega * y ; // Cartesian components of solid rotation
71 u_euler.set(2) = omega * x ;
72 u_euler.set(3) = 0 ;
74
75 u_euler.set_triad( mp.get_bvect_cart() ) ; // Triad = Cartesian triad
76
77 u_euler = ( u_euler + beta ) / nn ;
78
79 u_euler.std_spectral_base() ; // sets the standard bases for spectral expansions
80
81 if ( (u_euler(1).get_etat() == ETATZERO) &&
82 (u_euler(2).get_etat() == ETATZERO) &&
83 (u_euler(3).get_etat() == ETATZERO) ) {
84
86 }
87
88
89 // Computation of uuu (norme of u_euler)
90 // ------------------
91
92 // The scalar product is performed on the spherical components:
93
94 Vector us = u_euler ;
95 us.change_triad( mp.get_bvect_spher() ) ;
96
97 Scalar uuu2 = a_car * ( us(1)*us(1) + us(2)*us(2) ) + b_car * us(3)*us(3) ;
98
99 uuu = sqrt( uuu2 ) ;
100
101 if (uuu.get_etat() == ETATQCQ) {
102 // Same basis as (Omega -N^phi) r sin(theta)
103 (uuu.set_spectral_va()).set_base( us(3).get_spectral_va().get_base() ) ;
104 }
105
106 // Lorentz factor
107 // --------------
108
109 Scalar u2 = unsurc2 * uuu2 ;
110
111 Scalar gam2 = double(1) / (double(1) - u2) ;
112
113 gam_euler = sqrt(gam2) ;
114
115 gam_euler.std_spectral_base() ; // sets the standard spectral bases for
116 // a scalar field
117
118 // Energy density E with respect to the Eulerian observer
119 //--------------------------------------------------------
120
121 ener_euler = gam2 * ( ener + press ) - press ;
122
123 // Trace of the stress tensor with respect to the Eulerian observer
124 //------------------------------------
125
126 s_euler = 3 * press + ( ener_euler + press ) * u2 ;
127
128 // The derived quantities are obsolete
129 // -----------------------------------
130
131 del_deriv() ;
132
133}
134}
Active physical coordinates and mapping derivatives.
Definition coord.h:90
Time evolution with partial storage (*** under development ***).
Definition evolution.h:371
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
Tensor field of valence 0 (or component of a tensorial field).
Definition scalar.h:387
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
Valeur & set_spectral_va()
Returns va (read/write version)
Definition scalar.h:604
int get_etat() const
Returns the logical state ETATNONDEF (undefined), ETATZERO (null) or ETATQCQ (ordinary).
Definition scalar.h:554
double unsurc2
: unsurc2=1 for a relativistic star, 0 for a Newtonian one.
Definition star_rot.h:99
Scalar b_car
Square of the metric factor B.
Definition star_rot.h:110
virtual void hydro_euler()
Computes the hydrodynamical quantities relative to the Eulerian observer from those in the fluid fram...
double omega
Rotation angular velocity ([f_unit] )
Definition star_rot.h:101
Scalar uuu
Norm of u_euler.
Definition star_rot.h:121
Scalar a_car
Square of the metric factor A.
Definition star_rot.h:104
virtual void del_deriv() const
Deletes all the derived quantities.
Definition star_rot.C:297
Scalar ener
Total energy density in the fluid frame.
Definition star.h:193
Scalar nn
Lapse function N .
Definition star.h:225
Scalar ener_euler
Total energy density in the Eulerian frame.
Definition star.h:198
Scalar gam_euler
Lorentz factor between the fluid and Eulerian observers.
Definition star.h:204
Scalar s_euler
Trace of the stress scalar in the Eulerian frame.
Definition star.h:201
Scalar press
Fluid pressure.
Definition star.h:194
Vector u_euler
Fluid 3-velocity with respect to the Eulerian observer.
Definition star.h:207
Map & mp
Mapping associated with the star.
Definition star.h:180
Vector beta
Shift vector.
Definition star.h:228
Tensor field of valence 1.
Definition vector.h:188
virtual void std_spectral_base()
Sets the standard spectal bases of decomposition for each component.
Definition vector.C:316
Scalar & set(int)
Read/write access to a component.
Definition vector.C:296
Cmp sqrt(const Cmp &)
Square root.
Definition cmp_math.C:220
void set_triad(const Base_vect &new_triad)
Assigns a new vectorial basis (triad) of decomposition.
Definition tensor.C:519
void annule_domain(int l)
Sets the Tensor to zero in a given domain.
Definition tensor.C:666
virtual void set_etat_qcq()
Sets the logical state of all components to ETATQCQ (ordinary state).
Definition tensor.C:481
virtual void set_etat_zero()
Sets the logical state of all components to ETATZERO (zero state).
Definition tensor.C:497
Lorene prototypes.
Definition app_hor.h:64