LORENE
star_rot_dirac_diff.C
1/*
2 * Methods of class Star_rot_Dirac_diff
3 *
4 * (see file star_rot_dirac_diff.h for documentation).
5 *
6 */
7
8/*
9 * Copyright (c) 2005 Motoyuki Saijo
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 version 2
15 * as published by the Free Software Foundation.
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
28char star_rot_dirac_diff_C[] = "$Header: /cvsroot/Lorene/C++/Source/Star/star_rot_dirac_diff.C,v 1.4 2014/10/13 08:53:39 j_novak Exp $" ;
29
30/*
31 * $Header: /cvsroot/Lorene/C++/Source/Star/star_rot_dirac_diff.C,v 1.4 2014/10/13 08:53:39 j_novak Exp $
32 *
33 */
34
35
36// C headers
37#include <cmath>
38#include <cassert>
39
40// Lorene headers
41#include "star_rot_dirac_diff.h"
42#include "unites.h"
43#include "utilitaires.h"
44
45
46 //--------------//
47 // Constructors //
48 //--------------//
49
50// Standard constructor
51//-------------------------
52namespace Lorene {
54 double (*frot_i)(double, const Tbl&),
55 double (*primfrot_i)(double, const Tbl&),
56 const Tbl& par_frot_i)
58 frot(frot_i),
59 primfrot(primfrot_i),
60 par_frot(par_frot_i),
61 omega_field(mpi),
62 prim_field(mpi)
63 {
64
65 // Initialization to a static state :
66 omega_field = 0 ;
67 prim_field = 0 ;
68 omega_min = 0 ;
69 omega_max = 0 ;
70
71}
72
73
74// Copy constructor
75//-----------------
77 : Star_rot_Dirac(star),
78 frot(star.frot),
79 primfrot(star.primfrot),
80 par_frot(star.par_frot),
81 omega_field(star.omega_field),
82 omega_min(star.omega_min),
83 omega_max(star.omega_max),
84 prim_field(star.prim_field)
85 {}
86
87
88//Constructor from a file //## to be more general...
89//------------------------
91 double (*frot_i)(double, const Tbl&),
92 double (*primfrot_i)(double, const Tbl&) )
94 frot(frot_i),
95 primfrot(primfrot_i),
96 par_frot(fich),
97 omega_field(mpi),
98 prim_field(mpi)
99 {
100
103 fait_prim_field() ;
104
105 // omega_min and omega_max are read in the file:
106 fread_be(&omega_min, sizeof(double), 1, fich) ;
107 fread_be(&omega_max, sizeof(double), 1, fich) ;
108
109}
110
111
112 //------------//
113 // Destructor //
114 //------------//
115
117
118
119 //---------------//
120 // Assignment //
121 //---------------//
122
123// Assignment to another Star_rot_Dirac_diff
124// ------------------------------------
125
127
128 // Assignment of quantities common to all the derived classes of
129 // Star_rot_Dirac
131
132 // Assignment of proper quantities of class Star_rot_Dirac
133 frot = star.frot ;
134 primfrot = star.primfrot ;
135 par_frot = star.par_frot ;
136 omega_field = star.omega_field ;
137 prim_field = star.prim_field ;
138 omega_min = star.omega_min ;
139 omega_max = star.omega_max ;
140
141}
142
143
144 //-----------//
145 // Outputs //
146 //-----------//
147
148// Save in a file
149// --------------
150
152
154
156
158
159 fwrite_be(&omega_min, sizeof(double), 1, fich) ;
160 fwrite_be(&omega_max, sizeof(double), 1, fich) ;
161
162 // What else to save? //## to be more general ...
163
164}
165
166
167// Printing
168// ---------
169
171
172 using namespace Unites ;
173
175
176 ost << "Differentially rotating star in Dirac gauge" << '\n';
177
178 // Only differentially rotating star for the moment....
179 ost << '\n';
180 ost << "Differentially rotating star" << '\n';
181 ost << "-----------------------" << '\n';
182
183 ost << '\n'<< "Parameters of F(Omega) : " << '\n';
184 ost << par_frot << '\n';
185
186 ost << "Min, Max of Omega/(2pi) : " << omega_min / (2*M_PI) * f_unit
187 << " Hz, " << omega_max / (2*M_PI) * f_unit << " Hz" << '\n';
188 int lsurf = nzet - 1;
189 int nt = mp.get_mg()->get_nt(lsurf) ;
190 int nr = mp.get_mg()->get_nr(lsurf) ;
191 ost << "Central, equatorial value of Omega: "
192 << omega_field.val_grid_point(0, 0, 0, 0) * f_unit
193 << " rad/s, "
194 << omega_field.val_grid_point(nzet-1, 0, nt-1, nr-1) * f_unit
195 << " rad/s" << '\n';
196
197 ost << "Central, equatorial value of Omega/(2 Pi): "
198 << omega_field.val_grid_point(0, 0, 0, 0) * f_unit / (2*M_PI)
199 << " Hz, "
200 << omega_field.val_grid_point(nzet-1, 0, nt-1, nr-1) *
201 f_unit / (2*M_PI)
202 << " Hz" << '\n';
203
204 ost << "Error on the virial identity GRV2 : " << '\n';
205 ost << "GRV2 = " << grv2() << '\n';
206 ost << "Error on the virial identity GRV3 : " << '\n';
207 ost << "GRV3 = " << grv3() << '\n';
208
209 ost << "Angular momentum J : "
210 << angu_mom()/( qpig / (4*M_PI) *msol*msol) << " G M_sol^2 / c"
211 << '\n';
212 ost << "c J / (G M^2) : "
213 << angu_mom()/( qpig / (4*M_PI) * pow(mass_g(), 2.) ) << '\n';
214
215 if (omega != 0.) {
216 double mom_iner = angu_mom() / omega ;
217 double mom_iner_38si = mom_iner * rho_unit * (pow(r_unit, double(5.))
218 / double(1.e38) ) ;
219 ost << "Moment of inertia: " << mom_iner_38si << " 10^38 kg m^2"
220 << '\n';
221 }
222
223 ost << "Ratio T/W : " << tsw() << '\n';
224
225// ost << "Omega_c / e_max^{1/2} : " << get_omega_c() * len_rho << '\n';
226
227 ost << "Circumferential equatorial radius R_circ : "
228 << r_circ()/km << " km" << '\n';
229 ost << "Circumferential polar radius Rp_circ : "
230 << rp_circ()/km << " km" << endl ;
231 ost << "Coordinate equatorial radius r_eq : " << ray_eq()/km << " km"
232 << endl ;
233 ost << "Flattening r_pole/r_eq : " << aplat() << endl ;
234 ost << "Ellipticity sqrt(1-(Rp_circ/R_circ)^2) : " << ellipt() << endl ;
235 double compact = qpig/(4.*M_PI) * mass_g() / r_circ() ;
236 ost << "Compaction parameter M_g / R_circ : " << compact << '\n';
237
238
239 // More to come here.....
240
241 return ost ;
242
243}
244
245
246 //-----------------------//
247 // Miscellaneous //
248 //-----------------------//
249
251
252 return frot(omeg, par_frot) ;
253
254}
255
257
258 return primfrot(omeg, par_frot) ;
259
260}
261
263
264 return omega_field.val_grid_point(0, 0, 0, 0) ;
265
266}
267}
Equation of state base class.
Definition eos.h:190
Time evolution with partial storage (*** under development ***).
Definition evolution.h:371
Base class for coordinate mappings.
Definition map.h:670
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 sauve(FILE *) const
Save in a file.
Definition scalar.C:686
double val_grid_point(int l, int k, int j, int i) const
Returns the value of the field at a specified grid point.
Definition scalar.h:637
Class for relativistic differentially rotating stars in Dirac gauge and maximal slicing.
virtual ~Star_rot_Dirac_diff()
Destructor.
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
double funct_omega(double omeg) const
Evaluates , where F is the function defining the rotation profile.
virtual void sauve(FILE *) const
Save in a file.
double prim_funct_omega(double omeg) const
Evaluates the primitive of , where F is the function defining the rotation profile.
Star_rot_Dirac_diff(Map &mp_i, int nzet_i, const Eos &eos_i, double(*frot_i)(double, const Tbl &), double(*primfrot_i)(double, const Tbl &), const Tbl &par_frot_i)
Standard constructor.
void fait_prim_field()
Computes the member prim_field from omega_field .
double(* frot)(double, const Tbl &)
Function defining the rotation profile.
double omega_min
Minimum value of .
virtual double get_omega_c() const
Returns the central value of the rotation angular velocity ([f_unit] )
double(* primfrot)(double, const Tbl &)
Primitive of the function , which vanishes at the stellar center.
virtual double tsw() const
Ratio T/W.
Tbl par_frot
Parameters of the function .
void operator=(const Star_rot_Dirac_diff &)
Assignment to another Star_rot_Dirac_diff.
double omega_max
Maximum value of .
Class for relativistic rotating stars in Dirac gauge and maximal slicing.
virtual double mass_g() const
Gravitational mass.
virtual double ellipt() const
Ellipticity e.
virtual double grv3() const
Error on the virial identity GRV3.
virtual double angu_mom() const
Angular momentum.
virtual double r_circ() const
Circumferential equatorial radius.
double omega
Rotation angular velocity ([f_unit] )
virtual double rp_circ() const
Circumferential polar radius.
virtual double aplat() const
Flattening r_pole/r_eq.
void operator=(const Star_rot_Dirac &)
Assignment to another Star_rot_Dirac.
virtual double grv2() const
Error on the virial identity GRV2.
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition star.C:417
double ray_eq() const
Coordinate radius at , [r_unit].
virtual void sauve(FILE *) const
Save in a file.
Definition star.C:397
Map & mp
Mapping associated with the star.
Definition star.h:180
int nzet
Number of domains of *mp occupied by the star.
Definition star.h:183
Basic array class.
Definition tbl.h:161
void sauve(FILE *) const
Save in a file.
Definition tbl.C:326
Cmp pow(const Cmp &, int)
Power .
Definition cmp_math.C:348
int fread_be(int *aa, int size, int nb, FILE *fich)
Reads integer(s) from a binary file according to the big endian convention.
Definition fread_be.C:69
int fwrite_be(const int *aa, int size, int nb, FILE *fich)
Writes integer(s) into a binary file according to the big endian convention.
Definition fwrite_be.C:70
Lorene prototypes.
Definition app_hor.h:64
Standard units of space, time and mass.