LORENE
binaire.h
1/*
2 * Definition of Lorene class Binaire
3 *
4 */
5
6/*
7 * Copyright (c) 2000-2003 Eric Gourgoulhon
8 *
9 * This file is part of LORENE.
10 *
11 * LORENE is free software; you can redistribute it and/or modify
12 * it under the terms of the GNU General Public License as published by
13 * the Free Software Foundation; either version 2 of the License, or
14 * (at your option) any later version.
15 *
16 * LORENE is distributed in the hope that it will be useful,
17 * but WITHOUT ANY WARRANTY; without even the implied warranty of
18 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19 * GNU General Public License for more details.
20 *
21 * You should have received a copy of the GNU General Public License
22 * along with LORENE; if not, write to the Free Software
23 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
24 *
25 */
26
27
28#ifndef __BINAIRE_H_
29#define __BINAIRE_H_
30
31/*
32 * $Id: binaire.h,v 1.6 2014/10/13 08:52:32 j_novak Exp $
33 * $Log: binaire.h,v $
34 * Revision 1.6 2014/10/13 08:52:32 j_novak
35 * Lorene classes and functions now belong to the namespace Lorene.
36 *
37 * Revision 1.5 2009/06/18 18:42:13 k_taniguchi
38 * Defined a slightly modified code to determine
39 * the orbital angular velocity.
40 *
41 * Revision 1.4 2003/09/15 15:09:47 e_gourgoulhon
42 * Added the member function write_global.
43 *
44 * Revision 1.3 2002/06/17 14:05:16 j_novak
45 * friend functions are now also declared outside the class definition
46 *
47 * Revision 1.2 2001/12/20 13:00:31 k_taniguchi
48 * Addition of the Komar mass, the virial error by Gourgoulhon and Bonazzola, and the virial error by Friedman, Uryu, and Shibata.
49 *
50 * Revision 1.1.1.1 2001/11/20 15:19:27 e_gourgoulhon
51 * LORENE
52 *
53 * Revision 2.11 2001/02/28 14:27:32 keisuke
54 * *** empty log message ***
55 *
56 * Revision 2.10 2000/07/07 14:10:04 eric
57 * AJout de display_poly.
58 *
59 * Revision 2.9 2000/03/17 15:26:53 eric
60 * Ajout de la fonction analytical_omega.
61 *
62 * Revision 2.8 2000/03/15 16:43:12 eric
63 * Ajout de la fonction analytical_shift().
64 *
65 * Revision 2.7 2000/03/13 14:24:50 eric
66 * Ajout des membres p_ham_constr et p_mom_constr ainsi que des
67 * fonctions de calcul associees (verification des equations de contrainte).
68 *
69 * Revision 2.6 2000/02/18 14:52:05 eric
70 * Ajout des membres p_virial et p_total_ener et des fonctions de calcul
71 * associees.
72 *
73 * Revision 2.5 2000/02/12 17:36:39 eric
74 * Ajout de la fonction separation().
75 *
76 * Revision 2.4 2000/02/12 17:09:02 eric
77 * Ajout de la fonction orbit.
78 *
79 * Revision 2.3 2000/02/04 17:15:05 eric
80 * Ajout du membre ref_triad.
81 *
82 * Revision 2.2 2000/02/02 10:13:08 eric
83 * Ajout des fonctions de lecture/ecriture omega, x_axe.
84 *
85 * Revision 2.1 2000/01/31 17:01:46 eric
86 * *** empty log message ***
87 *
88 * Revision 2.0 2000/01/31 15:57:39 eric
89 * *** empty log message ***
90 *
91 *
92 * $Header: /cvsroot/Lorene/C++/Include/binaire.h,v 1.6 2014/10/13 08:52:32 j_novak Exp $
93 *
94 */
95
96#include "etoile.h"
97
98namespace Lorene {
104class Binaire {
105
106 // Data :
107 // -----
108 private:
113
116
119
125
129 double omega ;
130
133 double x_axe ;
134
135 // Derived data :
136 // ------------
137 private:
139 mutable double* p_mass_adm ;
140
142 mutable double* p_mass_kom ;
143
145 mutable Tbl* p_angu_mom ;
146
148 mutable double* p_total_ener ;
149
151 mutable double* p_virial ;
152
154 mutable double* p_virial_gb ;
155
157 mutable double* p_virial_fus ;
158
160 mutable double* p_ham_constr ;
161
163 mutable Tbl* p_mom_constr ;
164
165
166 // Constructors - Destructor
167 // -------------------------
168 public:
185 Binaire(Map& mp1, int nzet1, const Eos& eos1, int irrot1,
186 Map& mp2, int nzet2, const Eos& eos2, int irrot2,
187 int relat) ;
188
189
190 Binaire(const Binaire& ) ;
191
201 Binaire(Map& mp1, const Eos& eos1, Map& mp2, const Eos& eos2,
202 FILE* fich) ;
203
204 ~Binaire() ;
205
206
207 // Memory management
208 // -----------------
209 private:
211 void del_deriv() const ;
212
214 void set_der_0x0() const ;
215
216 // Mutators / assignment
217 // ---------------------
218 public:
220 void operator=(const Binaire&) ;
221
224 { assert( (i==1) || (i==2) );
225 del_deriv() ;
226 return *et[i-1] ;} ;
227
229 double& set_omega() {return omega; } ;
230
232 double& set_x_axe() {return x_axe; } ;
233
234
235 // Accessors
236 // ---------
237 public:
239 const Etoile_bin& operator()(int i) const
240 { assert( (i==1) || (i==2) );
241 return *et[i-1] ;} ;
242
244 double get_omega() const {return omega; } ;
245
247 double get_x_axe() const {return x_axe; } ;
248
252 double separation() const ;
253
254 // Outputs
255 // -------
256 public:
257 void sauve(FILE *) const ;
258
260 friend ostream& operator<<(ostream& , const Binaire& ) ;
261
263 void display_poly(ostream& ) const ;
264
268 void write_global(ostream& ) const ;
269
270 private:
272 ostream& operator>>(ostream& ) const ;
273
274
275 // Computational routines
276 // ----------------------
277 public:
279 double mass_adm() const ;
280
282 double mass_kom() const ;
283
291 const Tbl& angu_mom() const ;
292
301 double total_ener() const ;
302
307 double virial() const ;
308
314 double virial_gb() const ;
315
323 double virial_fus() const ;
324
332 double ham_constr() const ;
333
340 const Tbl& mom_constr() const ;
341
368 void orbit(double fact_omeg_min, double fact_omeg_max, double& xgg1,
369 double& xgg2) ;
370
401 void orbit_eqmass(double fact_omeg_min, double fact_omeg_max,
402 double mass1, double mass2,
403 double& xgg1, double& xgg2) ;
404
408 void analytical_omega() ;
409
414 void analytical_shift() ;
415
416};
417ostream& operator<<(ostream& , const Binaire& ) ;
418
419}
420#endif
Cartesian vectorial bases (triads).
Definition base_vect.h:201
Binary systems.
Definition binaire.h:104
void write_global(ostream &) const
Write global quantities in a formatted file.
Definition binaire.C:334
double * p_ham_constr
Relative error on the Hamiltonian constraint.
Definition binaire.h:160
Tbl * p_angu_mom
Total angular momentum of the system.
Definition binaire.h:145
void analytical_shift()
Sets some analytical template for the shift vector (via the members {\tt w_shift} and {\tt khi_shift}...
double * p_virial_gb
Virial theorem error by E.Gourgoulhon and S.Bonazzola.
Definition binaire.h:154
double mass_adm() const
Total ADM mass.
void set_der_0x0() const
Sets to {\tt 0x0} all the pointers on derived quantities.
Definition binaire.C:195
void analytical_omega()
Sets the orbital angular velocity to some 2-PN analytical value (Keplerian value in the Newtonian cas...
double virial_gb() const
Estimates the relative error on the virial theorem calculated by E.Gourgoulhon and S....
Etoile_bin * et[2]
Array of the two stars (to perform loops on the stars): {\tt et[0]} contains the address of {\tt star...
Definition binaire.h:124
const Tbl & angu_mom() const
Total angular momentum.
double * p_mass_adm
Total ADM mass of the system.
Definition binaire.h:139
double separation() const
Returns the coordinate separation of the two stellar centers [{\tt r_unit}].
Definition binaire.C:457
void del_deriv() const
Destructor.
Definition binaire.C:177
friend ostream & operator<<(ostream &, const Binaire &)
Save in a file.
Definition binaire.C:251
double * p_total_ener
Total energy of the system.
Definition binaire.h:148
double * p_virial
Virial theorem error.
Definition binaire.h:151
Tbl * p_mom_constr
Relative error on the momentum constraint.
Definition binaire.h:163
double & set_omega()
Sets the orbital angular velocity [{\tt f_unit}].
Definition binaire.h:229
void orbit_eqmass(double fact_omeg_min, double fact_omeg_max, double mass1, double mass2, double &xgg1, double &xgg2)
Computes the orbital angular velocity {\tt omega} and the position of the rotation axis {\tt x_axe}.
double get_omega() const
Returns the orbital angular velocity [{\tt f_unit}].
Definition binaire.h:244
const Etoile_bin & operator()(int i) const
Returns a reference to the star no. i.
Definition binaire.h:239
Etoile_bin star2
Second star of the system.
Definition binaire.h:118
double mass_kom() const
Total Komar mass.
double * p_mass_kom
Total Komar mass of the system.
Definition binaire.h:142
const Base_vect_cart ref_triad
Cartesian triad of the absolute reference frame.
Definition binaire.h:112
double virial() const
Estimates the relative error on the virial theorem (for a relativistic one, it returns $|1 - M_{\rm K...
void display_poly(ostream &) const
Display in polytropic units.
Definition binaire.C:284
const Tbl & mom_constr() const
Estimates the relative error on the momentum constraint equation by comparing ${\overline\nabla}_j K^...
double & set_x_axe()
Sets the absolute coordinate X of the rotation axis [{\tt r_unit}].
Definition binaire.h:232
double omega
Angular velocity with respect to an asymptotically inertial observer.
Definition binaire.h:129
double * p_virial_fus
Virial theorem error by J.L.Friedman, K.Uryu, and M.Shibata.
Definition binaire.h:157
void operator=(const Binaire &)
Assignment to another {\tt Binaire}.
Definition binaire.C:217
double get_x_axe() const
Returns the absolute coordinate X of the rotation axis [{\tt r_unit}].
Definition binaire.h:247
Etoile_bin star1
First star of the system.
Definition binaire.h:115
void orbit(double fact_omeg_min, double fact_omeg_max, double &xgg1, double &xgg2)
Computes the orbital angular velocity {\tt omega} and the position of the rotation axis {\tt x_axe}.
double total_ener() const
Total energy (excluding the rest mass energy).
ostream & operator>>(ostream &) const
Operator >> (function called by the operator <<).
Definition binaire.C:257
double virial_fus() const
Estimates the relative error on the virial theorem calculated by J.L.Friedman, K.Uryu,...
Etoile_bin & set(int i)
Read/write of the star no. i.
Definition binaire.h:223
double x_axe
Absolute X coordinate of the rotation axis.
Definition binaire.h:133
double ham_constr() const
Estimates the relative error on the Hamiltonian constraint equation by comparing $\underline\Delta\ln...
Equation of state base class.
Definition eos.h:190
Class for stars in binary system.
Definition etoile.h:814
Base class for coordinate mappings.
Definition map.h:670
Basic array class.
Definition tbl.h:161
Lorene prototypes.
Definition app_hor.h:64