LORENE
binary.C
1/*
2 * Methods of class Binary
3 *
4 */
5
6/*
7 * Copyright (c) 2004 Francois Limousin
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 */
26char Binary_C[] = "$Header: /cvsroot/Lorene/C++/Source/Binary/binary.C,v 1.16 2014/10/13 08:52:44 j_novak Exp $" ;
27
28/*
29 * $Id: binary.C,v 1.16 2014/10/13 08:52:44 j_novak Exp $
30 * $Log: binary.C,v $
31 * Revision 1.16 2014/10/13 08:52:44 j_novak
32 * Lorene classes and functions now belong to the namespace Lorene.
33 *
34 * Revision 1.15 2014/10/06 15:12:59 j_novak
35 * Modified #include directives to use c++ syntax.
36 *
37 * Revision 1.14 2005/09/15 14:39:14 e_gourgoulhon
38 * Added printing of angular momentum in display_poly.
39 *
40 * Revision 1.13 2005/09/13 19:38:31 f_limousin
41 * Reintroduction of the resolution of the equations in cartesian coordinates.
42 *
43 * Revision 1.12 2005/02/24 17:31:27 f_limousin
44 * Update of the function decouple().
45 *
46 * Revision 1.11 2005/02/18 13:14:06 j_novak
47 * Changing of malloc/free to new/delete + suppression of some unused variables
48 * (trying to avoid compilation warnings).
49 *
50 * Revision 1.10 2004/07/21 11:47:10 f_limousin
51 * Add / Delete comments.
52 *
53 * Revision 1.9 2004/05/25 14:25:12 f_limousin
54 * Add the virial theorem for conformally flat configurations.
55 *
56 * Revision 1.8 2004/03/25 10:29:01 j_novak
57 * All LORENE's units are now defined in the namespace Unites (in file unites.h).
58 *
59 * Revision 1.7 2004/03/23 10:00:47 f_limousin
60 * Minor changes.
61 *
62 * Revision 1.6 2004/02/27 09:59:33 f_limousin
63 * Modification in the routine decouple().
64 *
65 * Revision 1.5 2004/02/21 17:05:12 e_gourgoulhon
66 * Method Scalar::point renamed Scalar::val_grid_point.
67 * Method Scalar::set_point renamed Scalar::set_grid_point.
68 *
69 * Revision 1.4 2004/01/20 15:21:12 f_limousin
70 * First version
71 *
72 *
73 * $Header: /cvsroot/Lorene/C++/Source/Binary/binary.C,v 1.16 2014/10/13 08:52:44 j_novak Exp $
74 *
75 */
76
77// Headers C
78#include <cmath>
79
80// Headers Lorene
81#include "binary.h"
82#include "eos.h"
83#include "utilitaires.h"
84#include "graphique.h"
85#include "param.h"
86#include "unites.h"
87
88 //--------------//
89 // Constructors //
90 //--------------//
91
92// Standard constructor
93// --------------------
94
95namespace Lorene {
96Binary::Binary(Map& mp1, int nzet1, const Eos& eos1, int irrot1,
97 Map& mp2, int nzet2, const Eos& eos2, int irrot2,
98 int conf_flat)
99 : star1(mp1, nzet1, eos1, irrot1, conf_flat),
100 star2(mp2, nzet2, eos2, irrot2, conf_flat)
101{
102
103 et[0] = &star1 ;
104 et[1] = &star2 ;
105
106 omega = 0 ;
107 x_axe = 0 ;
108
109 // Pointers of derived quantities initialized to zero :
110 set_der_0x0() ;
111}
112
113// Copy constructor
114// ----------------
116 : star1(bibi.star1),
117 star2(bibi.star2),
118 omega(bibi.omega),
119 x_axe(bibi.x_axe)
120{
121 et[0] = &star1 ;
122 et[1] = &star2 ;
123
124 // Pointers of derived quantities initialized to zero :
125 set_der_0x0() ;
126}
127
128// Constructor from a file
129// -----------------------
130Binary::Binary(Map& mp1, const Eos& eos1, Map& mp2, const Eos& eos2,
131 FILE* fich)
132 : star1(mp1, eos1, fich),
133 star2(mp2, eos2, fich)
134{
135 et[0] = &star1 ;
136 et[1] = &star2 ;
137
138 // omega and x_axe are read in the file:
139 fread_be(&omega, sizeof(double), 1, fich) ;
140 fread_be(&x_axe, sizeof(double), 1, fich) ;
141
142 // Pointers of derived quantities initialized to zero :
143 set_der_0x0() ;
144
145}
146
147 //------------//
148 // Destructor //
149 //------------//
150
152
153 del_deriv() ;
154
155}
156
157 //----------------------------------//
158 // Management of derived quantities //
159 //----------------------------------//
160
161void Binary::del_deriv() const {
162
163 if (p_mass_adm != 0x0) delete p_mass_adm ;
164 if (p_mass_kom != 0x0) delete p_mass_kom ;
165 if (p_angu_mom != 0x0) delete p_angu_mom ;
166 if (p_total_ener != 0x0) delete p_total_ener ;
167 if (p_virial != 0x0) delete p_virial ;
168 if (p_ham_constr != 0x0) delete p_ham_constr ;
169 if (p_mom_constr != 0x0) delete p_mom_constr ;
170
171 set_der_0x0() ;
172}
173
174
175
176
178
179 p_mass_adm = 0x0 ;
180 p_mass_kom = 0x0 ;
181 p_angu_mom = 0x0 ;
182 p_total_ener = 0x0 ;
183 p_virial = 0x0 ;
184 p_ham_constr = 0x0 ;
185 p_mom_constr = 0x0 ;
186
187}
188
189
190 //--------------//
191 // Assignment //
192 //--------------//
193
194// Assignment to another Binary
195// --------------------------------
196
197void Binary::operator=(const Binary& bibi) {
198
199 star1 = bibi.star1 ;
200 star2 = bibi.star2 ;
201
202 omega = bibi.omega ;
203 x_axe = bibi.x_axe ;
204
205 del_deriv() ; // Deletes all derived quantities
206
207}
208
209 //--------------//
210 // Outputs //
211 //--------------//
212
213// Save in a file
214// --------------
215void Binary::sauve(FILE* fich) const {
216
217 star1.sauve(fich) ;
218 star2.sauve(fich) ;
219
220 fwrite_be(&omega, sizeof(double), 1, fich) ;
221 fwrite_be(&x_axe, sizeof(double), 1, fich) ;
222
223}
224
225// Printing
226// --------
227ostream& operator<<(ostream& ost, const Binary& bibi) {
228 bibi >> ost ;
229 return ost ;
230}
231
232
233ostream& Binary::operator>>(ostream& ost) const {
234
235 using namespace Unites ;
236
237 ost << endl ;
238 ost << "Binary neutron stars" << endl ;
239 ost << "=============" << endl ;
240 ost << endl <<
241 "Orbital angular velocity : " << omega * f_unit << " rad/s" << endl ;
242 ost << endl <<
243 "Coordinate separation between the two stellar centers : "
244 << separation() / km << " km" << endl ;
245 ost <<
246 "Absolute coordinate X of the rotation axis : " << x_axe / km
247 << " km" << endl ;
248 ost << endl << "Star 1 : " << endl ;
249 ost << "====== " << endl ;
250 ost << star1 << endl ;
251 ost << "Star 2 : " << endl ;
252 ost << "====== " << endl ;
253 ost << star2 << endl ;
254 return ost ;
255}
256
257// Display in polytropic units
258// ---------------------------
259
260void Binary::display_poly(ostream& ost) const {
261
262 using namespace Unites ;
263
264 const Eos* p_eos1 = &( star1.get_eos() ) ;
265 const Eos_poly* p_eos_poly = dynamic_cast<const Eos_poly*>( p_eos1 ) ;
266
267 if (p_eos_poly != 0x0) {
268
269 assert( star1.get_eos() == star2.get_eos() ) ;
270
271 double kappa = p_eos_poly->get_kap() ;
272 double gamma = p_eos_poly->get_gam() ; ;
273 double kap_ns2 = pow( kappa, 0.5 /(gamma-1) ) ;
274
275 // Polytropic unit of length in terms of r_unit :
276 double r_poly = kap_ns2 / sqrt(ggrav) ;
277
278 // Polytropic unit of time in terms of t_unit :
279 double t_poly = r_poly ;
280
281 // Polytropic unit of mass in terms of m_unit :
282 double m_poly = r_poly / ggrav ;
283
284 // Polytropic unit of angular momentum in terms of j_unit :
285 double j_poly = r_poly * r_poly / ggrav ;
286
287 ost.precision(10) ;
288 ost << endl << "Quantities in polytropic units : " << endl ;
289 ost << "==============================" << endl ;
290 ost << " ( r_poly = " << r_poly / km << " km )" << endl ;
291 ost << " d_e_max : " << separation() / r_poly << endl ;
292 ost << " d_G : "
293 << ( star2.xa_barycenter() - star1.xa_barycenter() ) / r_poly
294 << endl ;
295 ost << " Omega : " << omega * t_poly << endl ;
296 ost << " J : " << angu_mom()(2) / j_poly << endl ;
297 ost << " M_ADM : " << mass_adm() / m_poly << endl ;
298 ost << " M_Komar : " << mass_kom() / m_poly << endl ;
299 ost << " M_bar(star 1) : " << star1.mass_b() / m_poly << endl ;
300 ost << " M_bar(star 2) : " << star2.mass_b() / m_poly << endl ;
301 ost << " R_0(star 1) : " <<
302 0.5 * ( star1.ray_eq() + star1.ray_eq_pi() ) / r_poly << endl ;
303 ost << " R_0(star 2) : " <<
304 0.5 * ( star2.ray_eq() + star2.ray_eq_pi() ) / r_poly << endl ;
305
306 }
307
308
309}
310
311
313
314 int nz_un = star1.get_mp().get_mg()->get_nzone() ;
315 int nz_deux = star2.get_mp().get_mg()->get_nzone() ;
316
317 // We determine R_limite :
318 double distance = star2.get_mp().get_ori_x() - star1.get_mp().get_ori_x() ;
319 cout << "distance = " << distance << endl ;
320 double lim_un = distance/2. ;
321 double lim_deux = distance/2. ;
322 double int_un = distance/6. ;
323 double int_deux = distance/6. ;
324
325 // The functions used.
326 Scalar fonction_f_un (star1.get_mp()) ;
327 fonction_f_un = 0.5*pow(
328 cos((star1.get_mp().r-int_un)*M_PI/2./(lim_un-int_un)), 2.)+0.5 ;
329 fonction_f_un.std_spectral_base();
330
331 Scalar fonction_g_un (star1.get_mp()) ;
332 fonction_g_un = 0.5*pow
333 (sin((star1.get_mp().r-int_un)*M_PI/2./(lim_un-int_un)), 2.) ;
334 fonction_g_un.std_spectral_base();
335
336 Scalar fonction_f_deux (star2.get_mp()) ;
337 fonction_f_deux = 0.5*pow(
338 cos((star2.get_mp().r-int_deux)*M_PI/2./(lim_deux-int_deux)), 2.)+0.5 ;
339 fonction_f_deux.std_spectral_base();
340
341 Scalar fonction_g_deux (star2.get_mp()) ;
342 fonction_g_deux = 0.5*pow(
343 sin((star2.get_mp().r-int_deux)*M_PI/2./(lim_un-int_deux)), 2.) ;
344 fonction_g_deux.std_spectral_base();
345
346 // The functions total :
347 Scalar decouple_un (star1.get_mp()) ;
348 decouple_un.allocate_all() ;
349 Scalar decouple_deux (star2.get_mp()) ;
350 decouple_deux.allocate_all() ;
351
352 Mtbl xabs_un (star1.get_mp().xa) ;
353 Mtbl yabs_un (star1.get_mp().ya) ;
354 Mtbl zabs_un (star1.get_mp().za) ;
355
356 Mtbl xabs_deux (star2.get_mp().xa) ;
357 Mtbl yabs_deux (star2.get_mp().ya) ;
358 Mtbl zabs_deux (star2.get_mp().za) ;
359
360 double xabs, yabs, zabs, air_un, air_deux, theta, phi ;
361
362 for (int l=0 ; l<nz_un ; l++) {
363 int nr = star1.get_mp().get_mg()->get_nr (l) ;
364
365 if (l==nz_un-1)
366 nr -- ;
367
368 int np = star1.get_mp().get_mg()->get_np (l) ;
369 int nt = star1.get_mp().get_mg()->get_nt (l) ;
370
371 for (int k=0 ; k<np ; k++)
372 for (int j=0 ; j<nt ; j++)
373 for (int i=0 ; i<nr ; i++) {
374
375 xabs = xabs_un (l, k, j, i) ;
376 yabs = yabs_un (l, k, j, i) ;
377 zabs = zabs_un (l, k, j, i) ;
378
379 // Coordinates of the point
381 (xabs, yabs, zabs, air_un, theta, phi) ;
383 (xabs, yabs, zabs, air_deux, theta, phi) ;
384
385 if (air_un <= lim_un)
386 if (air_un < int_un)
387 decouple_un.set_grid_point(l, k, j, i) = 1 ;
388 else
389 // Close to star 1 :
390 decouple_un.set_grid_point(l, k, j, i) =
391 fonction_f_un.val_grid_point(l, k, j, i) ;
392 else
393 if (air_deux <= lim_deux)
394 if (air_deux < int_deux)
395 decouple_un.set_grid_point(l, k, j, i) = 0 ;
396 else
397 // Close to star 2 :
398 decouple_un.set_grid_point(l, k, j, i) =
399 fonction_g_deux.val_point (air_deux, theta, phi) ;
400
401 else
402 // Far from each star :
403 decouple_un.set_grid_point(l, k, j, i) = 0.5 ;
404 }
405
406 // Case infinity :
407 if (l==nz_un-1)
408 for (int k=0 ; k<np ; k++)
409 for (int j=0 ; j<nt ; j++)
410 decouple_un.set_grid_point(nz_un-1, k, j, nr)=0.5 ;
411 }
412
413 for (int l=0 ; l<nz_deux ; l++) {
414 int nr = star2.get_mp().get_mg()->get_nr (l) ;
415
416 if (l==nz_deux-1)
417 nr -- ;
418
419 int np = star2.get_mp().get_mg()->get_np (l) ;
420 int nt = star2.get_mp().get_mg()->get_nt (l) ;
421
422 for (int k=0 ; k<np ; k++)
423 for (int j=0 ; j<nt ; j++)
424 for (int i=0 ; i<nr ; i++) {
425
426 xabs = xabs_deux (l, k, j, i) ;
427 yabs = yabs_deux (l, k, j, i) ;
428 zabs = zabs_deux (l, k, j, i) ;
429
430 // coordinates of the point :
432 (xabs, yabs, zabs, air_un, theta, phi) ;
434 (xabs, yabs, zabs, air_deux, theta, phi) ;
435
436 if (air_deux <= lim_deux)
437 if (air_deux < int_deux)
438 decouple_deux.set_grid_point(l, k, j, i) = 1 ;
439 else
440 // close to star two :
441 decouple_deux.set_grid_point(l, k, j, i) =
442 fonction_f_deux.val_grid_point(l, k, j, i) ;
443 else
444 if (air_un <= lim_un)
445 if (air_un < int_un)
446 decouple_deux.set_grid_point(l, k, j, i) = 0 ;
447 else
448 // close to star one :
449 decouple_deux.set_grid_point(l, k, j, i) =
450 fonction_g_un.val_point (air_un, theta, phi) ;
451
452 else
453 // Far from each star :
454 decouple_deux.set_grid_point(l, k, j, i) = 0.5 ;
455 }
456
457 // Case infinity :
458 if (l==nz_deux-1)
459 for (int k=0 ; k<np ; k++)
460 for (int j=0 ; j<nt ; j++)
461 decouple_deux.set_grid_point(nz_un-1, k, j, nr)=0.5 ;
462 }
463
464 decouple_un.std_spectral_base() ;
465 decouple_deux.std_spectral_base() ;
466
467 int nr = star1.get_mp().get_mg()->get_nr (1) ;
468 int nt = star1.get_mp().get_mg()->get_nt (1) ;
469 int np = star1.get_mp().get_mg()->get_np (1) ;
470 cout << "decouple_un" << endl << norme(decouple_un/(nr*nt*np)) << endl ;
471 cout << "decouple_deux" << endl << norme(decouple_deux/(nr*nt*np))
472 << endl ;
473
474 star1.decouple = decouple_un ;
475 star2.decouple = decouple_deux ;
476}
477
478void Binary::write_global(ostream& ost) const {
479
480 using namespace Unites ;
481
482 const Map& mp1 = star1.get_mp() ;
483 const Mg3d* mg1 = mp1.get_mg() ;
484 int nz1 = mg1->get_nzone() ;
485
486 ost.precision(5) ;
487 ost << "# Grid 1 : " << nz1 << "x"
488 << mg1->get_nr(0) << "x" << mg1->get_nt(0) << "x" << mg1->get_np(0)
489 << " R_out(l) [km] : " ;
490 for (int l=0; l<nz1; l++) {
491 ost << " " << mp1.val_r(l, 1., M_PI/2, 0) / km ;
492 }
493 ost << endl ;
494
495 ost << "# VE(M) " << endl ;
496
497
498 ost.setf(ios::scientific) ;
499 ost.width(14) ;
500 ost << virial() << endl ;
501
502 ost << "# d [km] "
503 << " d_G [km] "
504 << " d/(a1 +a1') "
505 << " f [Hz] "
506 << " M_ADM [M_sol] "
507 << " M_ADM_vol [M_sol] "
508 << " M_Komar [M_sol] "
509 << " M_Komar_vol [M_sol] "
510 << " J [G M_sol^2/c] " << endl ;
511
512 ost.precision(14) ;
513 ost.width(20) ;
514 ost << separation() / km ; ost.width(22) ;
515 ost << ( star2.xa_barycenter() - star1.xa_barycenter() ) / km ; ost.width(22) ;
516 ost << separation() / (star1.ray_eq() + star2.ray_eq()) ; ost.width(22) ;
517 ost << omega / (2*M_PI)* f_unit ; ost.width(22) ;
518 ost << mass_adm() / msol ; ost.width(22) ;
519 ost << mass_adm_vol() / msol ; ost.width(22) ;
520 ost << mass_kom() / msol ; ost.width(22) ;
521 ost << mass_kom_vol() / msol ; ost.width(22) ;
522 ost << angu_mom()(2)/ ( qpig / (4* M_PI) * msol*msol) << endl ;
523
524 ost << "# H_c(1)[c^2] "
525 << " e_c(1)[rho_nuc] "
526 << " M_B(1) [M_sol] "
527 << " r_eq(1) [km] "
528 << " a2/a1(1) "
529 << " a3/a1(1) " << endl ;
530
531 ost.width(20) ;
532 ost << star1.get_ent().val_grid_point(0,0,0,0) ; ost.width(22) ;
533 ost << star1.get_ener().val_grid_point(0,0,0,0) ; ost.width(22) ;
534 ost << star1.mass_b() / msol ; ost.width(22) ;
535 ost << star1.ray_eq() / km ; ost.width(22) ;
536 ost << star1.ray_eq_pis2() / star1.ray_eq() ; ost.width(22) ;
537 ost << star1.ray_pole() / star1.ray_eq() << endl ;
538
539 ost << "# H_c(2)[c^2] "
540 << " e_c(2)[rho_nuc] "
541 << " M_B(2) [M_sol] "
542 << " r_eq(2) [km] "
543 << " a2/a1(2) "
544 << " a3/a1(2) " << endl ;
545
546 ost.width(20) ;
547 ost << star2.get_ent().val_grid_point(0,0,0,0) ; ost.width(22) ;
548 ost << star2.get_ener().val_grid_point(0,0,0,0) ; ost.width(22) ;
549 ost << star2.mass_b() / msol ; ost.width(22) ;
550 ost << star2.ray_eq() / km ; ost.width(22) ;
551 ost << star2.ray_eq_pis2() / star1.ray_eq() ; ost.width(22) ;
552 ost << star2.ray_pole() / star1.ray_eq() << endl ;
553
554 // Quantities in polytropic units if the EOS is a polytropic one
555 // -------------------------------------------------------------
556 const Eos* p_eos1 = &( star1.get_eos() ) ;
557 const Eos_poly* p_eos_poly = dynamic_cast<const Eos_poly*>( p_eos1 ) ;
558
559 if ((p_eos_poly != 0x0) && ( star1.get_eos() == star2.get_eos() )) {
560
561 double kappa = p_eos_poly->get_kap() ;
562 double gamma = p_eos_poly->get_gam() ; ;
563 double kap_ns2 = pow( kappa, 0.5 /(gamma-1.) ) ;
564
565 // Polytropic unit of length in terms of r_unit :
566 double r_poly = kap_ns2 / sqrt(ggrav) ;
567
568 // Polytropic unit of time in terms of t_unit :
569 double t_poly = r_poly ;
570
571 // Polytropic unit of mass in terms of m_unit :
572 double m_poly = r_poly / ggrav ;
573
574 // Polytropic unit of angular momentum in terms of j_unit :
575 double j_poly = r_poly * r_poly / ggrav ;
576
577 ost << "# d [poly] "
578 << " d_G [poly] "
579 << " Omega [poly] "
580 << " M_ADM [poly] "
581 << " J [poly] "
582 << " M_B(1) [poly] "
583 << " M_B(2) [poly] " << endl ;
584
585 ost.width(20) ;
586 ost << separation() / r_poly ; ost.width(22) ;
587 ost << ( star2.xa_barycenter() - star1.xa_barycenter() ) / r_poly ; ost.width(22) ;
588 ost << omega * t_poly ; ost.width(22) ;
589 ost << mass_adm() / m_poly ; ost.width(22) ;
590 ost << angu_mom()(2) / j_poly ; ost.width(22) ;
591 ost << star1.mass_b() / m_poly ; ost.width(22) ;
592 ost << star2.mass_b() / m_poly << endl ;
593
594 }
595
596}
597
598
599
600 //-------------------------------//
601 // Miscellaneous //
602 //-------------------------------//
603
604double Binary::separation() const {
605
606 double dx = star1.mp.get_ori_x() - star2.mp.get_ori_x() ;
607 double dy = star1.mp.get_ori_y() - star2.mp.get_ori_y() ;
608 double dz = star1.mp.get_ori_z() - star2.mp.get_ori_z() ;
609
610 return sqrt( dx*dx + dy*dy + dz*dz ) ;
611
612}
613}
Binary systems.
Definition binary.h:73
double mass_adm_vol() const
Total ADM mass (computed by a volume integral)
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
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
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
Tbl * p_angu_mom
Total angular momentum of the system.
Definition binary.h:111
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.
~Binary()
Destructor.
Definition binary.C:151
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.
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
Binary(Map &mp1, int nzet1, const Eos &eos1, int irrot1, Map &mp2, int nzet2, const Eos &eos2, int irrot2, int conf_flat)
Standard constructor.
Definition binary.C:96
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
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
Polytropic equation of state (relativistic case).
Definition eos.h:757
double get_gam() const
Returns the adiabatic index (cf. Eq. (3))
Definition eos_poly.C:256
double get_kap() const
Returns the pressure coefficient (cf.
Definition eos_poly.C:260
Equation of state base class.
Definition eos.h:190
Base class for coordinate mappings.
Definition map.h:670
double get_ori_z() const
Returns the z coordinate of the origin.
Definition map.h:772
Coord ya
Absolute y coordinate.
Definition map.h:731
Coord r
r coordinate centered on the grid
Definition map.h:718
void convert_absolute(double xx, double yy, double zz, double &rr, double &theta, double &pphi) const
Determines the coordinates corresponding to given absolute Cartesian coordinates (X,...
Definition map.C:302
double get_ori_y() const
Returns the y coordinate of the origin.
Definition map.h:770
Coord za
Absolute z coordinate.
Definition map.h:732
double get_ori_x() const
Returns the x coordinate of the origin.
Definition map.h:768
Coord xa
Absolute x coordinate.
Definition map.h:730
virtual double val_r(int l, double xi, double theta, double pphi) const =0
Returns the value of the radial coordinate r for a given in a given domain.
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition map.h:765
Multi-domain grid.
Definition grilles.h:273
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
Definition grilles.h:462
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
Definition grilles.h:457
int get_nzone() const
Returns the number of domains.
Definition grilles.h:448
int get_nr(int l) const
Returns the number of points in the radial direction ( ) in domain no. l.
Definition grilles.h:452
Multi-domain array.
Definition mtbl.h:118
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
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
virtual void allocate_all()
Sets the logical state to ETATQCQ (ordinary state) and performs the memory allocation of all the elem...
Definition scalar.C:365
double & set_grid_point(int l, int k, int j, int i)
Setting the value of the field at a given grid point.
Definition scalar.h:684
double val_point(double r, double theta, double phi) const
Computes the value of the field at an arbitrary point , by means of the spectral expansion.
Definition scalar.C:890
Scalar decouple
Function used to construct the part generated by the star from the total .
Definition star.h:676
virtual double xa_barycenter() const
Absolute coordinate X of the barycenter of the baryon density,.
virtual double mass_b() const
Baryon mass.
virtual void sauve(FILE *) const
Save in a file.
Definition star_bin.C:489
const Map & get_mp() const
Returns the mapping.
Definition star.h:355
const Eos & get_eos() const
Returns the equation of state.
Definition star.h:361
const Scalar & get_ener() const
Returns the proper total energy density.
Definition star.h:370
double ray_eq() const
Coordinate radius at , [r_unit].
double ray_eq_pis2() const
Coordinate radius at , [r_unit].
const Scalar & get_ent() const
Returns the enthalpy field.
Definition star.h:364
double ray_eq_pi() const
Coordinate radius at , [r_unit].
Map & mp
Mapping associated with the star.
Definition star.h:180
double ray_pole() const
Coordinate radius at [r_unit].
Cmp sqrt(const Cmp &)
Square root.
Definition cmp_math.C:220
Cmp sin(const Cmp &)
Sine.
Definition cmp_math.C:69
Tbl norme(const Cmp &)
Sums of the absolute values of all the values of the Cmp in each domain.
Definition cmp_math.C:481
Cmp pow(const Cmp &, int)
Power .
Definition cmp_math.C:348
Cmp cos(const Cmp &)
Cosine.
Definition cmp_math.C:94
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.