LORENE
binary.h
1/*
2 * Definition of Lorene class Binary
3 *
4 */
5
6/*
7 * Copyright (c) 2004 Limousin Francois
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 version 2
13 * as published by the Free Software Foundation.
14 *
15 * LORENE is distributed in the hope that it will be useful,
16 * but WITHOUT ANY WARRANTY; without even the implied warranty of
17 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18 * GNU General Public License for more details.
19 *
20 * You should have received a copy of the GNU General Public License
21 * along with LORENE; if not, write to the Free Software
22 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
23 *
24 */
25
26#ifndef __BINARY_H_
27#define __BINARY_H_
28
29/*
30 * $Id: binary.h,v 1.10 2014/10/13 08:52:32 j_novak Exp $
31 * $Log: binary.h,v $
32 * Revision 1.10 2014/10/13 08:52:32 j_novak
33 * Lorene classes and functions now belong to the namespace Lorene.
34 *
35 * Revision 1.9 2006/04/11 14:26:12 f_limousin
36 * New version of the code : improvement of the computation of some
37 * critical sources, estimation of the dirac gauge, helical symmetry...
38 *
39 * Revision 1.8 2005/11/08 20:17:55 f_limousin
40 * Add function dirac_gauge() used to impose Dirac gauge during an iteration.
41 *
42 * Revision 1.7 2005/09/15 15:56:28 e_gourgoulhon
43 * Made the documentation compliant with Doxygen.
44 *
45 * Revision 1.6 2005/09/13 19:38:32 f_limousin
46 * Reintroduction of the resolution of the equations in cartesian coordinates.
47 *
48 * Revision 1.5 2004/07/21 11:45:20 f_limousin
49 * Add function mass_adm_vol() to compute the ADM mass of the system
50 * with a volume integral instead of a surface one.
51 *
52 * Revision 1.4 2004/05/25 14:50:06 f_limousin
53 * Add the virial theorem for conformally flat configurations.
54 *
55 * Revision 1.3 2004/01/20 15:25:39 f_limousin
56 * Nex class binary
57 *
58 *
59 * $Header: /cvsroot/Lorene/C++/Include/binary.h,v 1.10 2014/10/13 08:52:32 j_novak Exp $
60 *
61 */
62
63// Lorene headers
64#include "star.h"
65
66namespace Lorene {
73class Binary {
74
75 // Data :
76 // -----
77 protected:
78
81
84
89 Star_bin* et[2] ;
90
94 double omega ;
95
98 double x_axe ;
99
100
101 // Derived data :
102 // ------------
103 protected:
105 mutable double* p_mass_adm ;
106
108 mutable double* p_mass_kom ;
109
111 mutable Tbl* p_angu_mom ;
112
114 mutable double* p_total_ener ;
115
117 mutable double* p_virial ;
118
120 mutable double* p_ham_constr ;
121
123 mutable Tbl* p_mom_constr ;
124
125
126
127 // Constructors - Destructor
128 // -------------------------
129 public:
145 Binary(Map& mp1, int nzet1, const Eos& eos1, int irrot1,
146 Map& mp2, int nzet2, const Eos& eos2, int irrot2,
147 int conf_flat) ;
148
149
150 Binary(const Binary& ) ;
151
161 Binary(Map& mp1, const Eos& eos1, Map& mp2, const Eos& eos2,
162 FILE* fich) ;
163
164 ~Binary() ;
165
166
167
168 // Memory management
169 // -----------------
170 protected:
171
173 void del_deriv() const ;
174
176 void set_der_0x0() const ;
177
178
179 // Mutators / assignment
180 // ---------------------
181 public:
183 void operator=(const Binary&) ;
184
186 Star_bin& set(int i)
187 { assert( (i==1) || (i==2) );
188 del_deriv() ;
189 return *et[i-1] ;} ;
190
192 double& set_omega() {return omega; } ;
193
195 double& set_x_axe() {return x_axe; } ;
196
197
198 // Accessors
199 // ---------
200 public:
202 const Star_bin& operator()(int i) const
203 { assert( (i==1) || (i==2) );
204 return *et[i-1] ;} ;
205
207 double get_omega() const {return omega; } ;
208
210 double get_x_axe() const {return x_axe; } ;
211
215 double separation() const ;
216
217
218 // Outputs
219 // -------
220 public:
221 void sauve(FILE *) const ;
222
224 friend ostream& operator<<(ostream& , const Binary& ) ;
225
227 void display_poly(ostream& ) const ;
228
232 void write_global(ostream& ) const ;
233
234 private:
236 ostream& operator>>(ostream& ) const ;
237
238
239 // Computational routines
240 // ----------------------
241 public:
242
244 double mass_adm() const ;
245
247 double mass_adm_vol() const ;
248
250 double mass_kom() const ;
251
253 double mass_kom_vol() const ;
254
262 const Tbl& angu_mom() const ;
263
272 double total_ener() const ;
273
276 double virial() const ;
277
280 double ham_constr() const ;
281
284 const Tbl& mom_constr() const ;
285
312 void orbit(double fact_omeg_min, double fact_omeg_max, double& xgg1,
313 double& xgg2) ;
314
318 void analytical_omega() ;
319
324 void analytical_shift() ;
325
334 void fait_decouple () ;
335
337 void dirac_gauge() ;
338
340 void helical() ;
341
342
343};
344ostream& operator<<(ostream& , const Binary& ) ;
345
346}
347#endif
Binary systems.
Definition binary.h:73
double & set_omega()
Sets the orbital angular velocity [f_unit].
Definition binary.h:192
const Star_bin & operator()(int i) const
Returns a reference to the star no. i.
Definition binary.h:202
double get_x_axe() const
Returns the absolute coordinate X of the rotation axis [r_unit].
Definition binary.h:210
double mass_adm_vol() const
Total ADM mass (computed by a volume integral)
double ham_constr() const
Estimates the relative error on the Hamiltonian constraint.
void sauve(FILE *) const
Save in a file.
Definition binary.C:215
void display_poly(ostream &) const
Display in polytropic units.
Definition binary.C:260
double mass_kom_vol() const
Total Komar mass (computed by a volume integral)
void write_global(ostream &) const
Write global quantities in a formatted file.
Definition binary.C:478
const Tbl & mom_constr() const
Estimates the relative error on the momentum constraint.
double & set_x_axe()
Sets the absolute coordinate X of the rotation axis [r_unit].
Definition binary.h:195
void del_deriv() const
Deletes all the derived quantities.
Definition binary.C:161
Star_bin * et[2]
Array of the two stars (to perform loops on the stars): et[0] contains the address of star1 and et[1]...
Definition binary.h:89
void analytical_shift()
Sets some analytical template for the shift vector (via the members w_shift and khi_shift of the two ...
double omega
Angular velocity with respect to an asymptotically inertial observer.
Definition binary.h:94
double * p_total_ener
Total energy of the system.
Definition binary.h:114
void analytical_omega()
Sets the orbital angular velocity to some 2-PN analytical value (Keplerian value in the Newtonian cas...
Tbl * p_angu_mom
Total angular momentum of the system.
Definition binary.h:111
double total_ener() const
Total energy (excluding the rest mass energy).
void dirac_gauge()
Function used to impose Dirac gauge during an iteration.
Star_bin star2
Second star of the system.
Definition binary.h:83
const Tbl & angu_mom() const
Total angular momentum.
void fait_decouple()
Calculates decouple which is used to obtain qq_auto by the formula : qq_auto = decouple * qq.
Definition binary.C:312
Star_bin star1
First star of the system.
Definition binary.h:80
double * p_mass_adm
Total ADM mass of the system.
Definition binary.h:105
double * p_virial
Virial theorem error.
Definition binary.h:117
double mass_adm() const
Total ADM mass.
void orbit(double fact_omeg_min, double fact_omeg_max, double &xgg1, double &xgg2)
Computes the orbital angular velocity omega and the position of the rotation axis x_axe.
~Binary()
Destructor.
Definition binary.C:151
friend ostream & operator<<(ostream &, const Binary &)
Display.
Definition binary.C:227
void operator=(const Binary &)
Assignment to another Binary.
Definition binary.C:197
void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition binary.C:177
double virial() const
Estimates the relative error on the virial theorem.
Star_bin & set(int i)
Read/write of the star no. i.
Definition binary.h:186
double x_axe
Absolute X coordinate of the rotation axis.
Definition binary.h:98
Tbl * p_mom_constr
Relative error on the momentum constraint.
Definition binary.h:123
double * p_mass_kom
Total Komar mass of the system.
Definition binary.h:108
double mass_kom() const
Total Komar mass.
double separation() const
Returns the coordinate separation of the two stellar centers [r_unit].
Definition binary.C:604
void helical()
Function testing the helical symmetry.
double * p_ham_constr
Relative error on the Hamiltonian constraint.
Definition binary.h:120
ostream & operator>>(ostream &) const
Operator >> (function called by the operator <<).
Definition binary.C:233
double get_omega() const
Returns the orbital angular velocity [f_unit].
Definition binary.h:207
Equation of state base class.
Definition eos.h:190
Base class for coordinate mappings.
Definition map.h:670
Class for stars in binary system.
Definition star.h:483
Basic array class.
Definition tbl.h:161
Lorene prototypes.
Definition app_hor.h:64