LORENE
et_rot_bifluid.C
1 /*
2  * Methods for two fluids rotating relativistic stars.
3  *
4  * See the file et_rot_bifluid.h for documentation
5  *
6  */
7 
8 /*
9  * Copyright (c) 2001 Jerome Novak
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 
30 
31 
32 /*
33  * $Id: et_rot_bifluid.C,v 1.17 2017/10/06 12:36:34 a_sourie Exp $
34  * $Log: et_rot_bifluid.C,v $
35  * Revision 1.17 2017/10/06 12:36:34 a_sourie
36  * Cleaning of tabulated 2-fluid EoS class + superfluid rotating star model.
37  *
38  * Revision 1.16 2016/12/05 16:17:53 j_novak
39  * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
40  *
41  * Revision 1.15 2015/06/11 13:50:19 j_novak
42  * Minor corrections
43  *
44  * Revision 1.14 2015/06/10 14:39:17 a_sourie
45  * New class Eos_bf_tabul for tabulated 2-fluid EoSs and associated functions for the computation of rotating stars with such EoSs.
46  *
47  * Revision 1.13 2014/10/13 08:52:57 j_novak
48  * Lorene classes and functions now belong to the namespace Lorene.
49  *
50  * Revision 1.12 2004/03/25 10:29:04 j_novak
51  * All LORENE's units are now defined in the namespace Unites (in file unites.h).
52  *
53  * Revision 1.11 2003/12/04 14:28:26 r_prix
54  * allow for the case of "slow-rot-style" EOS inversion, in which we need to adapt
55  * the inner domain to n_outer=0 instead of mu_outer=0 ...
56  * (this should only be used for comparison to analytic slow-rot solution!)
57  *
58  * Revision 1.10 2003/11/20 14:01:26 r_prix
59  * changed member names to better conform to Lorene coding standards:
60  * J_euler -> j_euler, EpS_euler -> enerps_euler, Delta_car -> delta_car
61  *
62  * Revision 1.9 2003/11/18 18:38:11 r_prix
63  * use of new member EpS_euler: matter sources in equilibrium() and global quantities
64  * no longer distinguish Newtonian/relativistic, as all terms should have the right limit...
65  *
66  * Revision 1.8 2003/11/17 13:49:43 r_prix
67  * - moved superluminal check into hydro_euler()
68  * - removed some warnings
69  *
70  * Revision 1.7 2003/11/13 12:07:57 r_prix
71  * *) changed xxx2 -> Delta_car
72  * *) added (non 2-fluid specific!) members sphph_euler J_euler
73  * *) more or less rewritten hydro_euler() to see if I understand it ;)
74  * - somewhat simplified and more adapted to the notation used in our notes/paper.
75  * - Main difference: u_euler is no longer used!!, the "output" instead
76  * consists of ener_euler, s_euler, sphph_euler and J_euler, which are
77  * the general 3+1 components for Tmunu.
78  *
79  * Revision 1.6 2003/09/17 08:27:50 j_novak
80  * New methods: mass_b1() and mass_b2().
81  *
82  * Revision 1.5 2002/10/18 08:42:58 j_novak
83  * Take into account the sign for uuu and uuu2
84  *
85  * Revision 1.4 2002/01/16 15:03:28 j_novak
86  * *** empty log message ***
87  *
88  * Revision 1.3 2002/01/11 14:09:34 j_novak
89  * Added newtonian version for 2-fluid stars
90  *
91  * Revision 1.2 2001/12/04 21:27:53 e_gourgoulhon
92  *
93  * All writing/reading to a binary file are now performed according to
94  * the big endian convention, whatever the system is big endian or
95  * small endian, thanks to the functions fwrite_be and fread_be
96  *
97  * Revision 1.1.1.1 2001/11/20 15:19:28 e_gourgoulhon
98  * LORENE
99  *
100  * Revision 1.3 2001/08/28 16:04:22 novak
101  * Use of new definition of relative velocity and new declarations for EOS
102  *
103  * Revision 1.2 2001/08/27 09:58:43 novak
104  * *** empty log message ***
105  *
106  * Revision 1.1 2001/06/22 15:39:17 novak
107  * Initial revision
108  *
109  *
110  * $Header: /cvsroot/Lorene/C++/Source/Etoile/et_rot_bifluid.C,v 1.17 2017/10/06 12:36:34 a_sourie Exp $
111  *
112  */
113 // Headers C
114 #include "math.h"
115 
116 // Headers Lorene
117 #include "et_rot_bifluid.h"
118 #include "nbr_spx.h"
119 #include "utilitaires.h"
120 #include "unites.h"
121 
122  //--------------//
123  // Constructors //
124  //--------------//
125 // Standard constructor
126 // --------------------
127 namespace Lorene {
128 Et_rot_bifluid::Et_rot_bifluid(Map& mpi, int nzet_i, bool relat, const Eos_bifluid& eos_i):
129  Etoile_rot(mpi, nzet_i, relat, *eos_i.trans2Eos()),
130  eos(eos_i),
131  ent2(mpi),
132  nbar2(mpi),
133  K_nn(mpi),
134  K_np(mpi),
135  K_pp(mpi),
136  alpha_eos(mpi),
137  sphph_euler(mpi),
138  j_euler(mpi, 1, CON, mp.get_bvect_cart()),
139  j_euler1 (mpi, 1, CON, mp.get_bvect_cart()),
140  j_euler2(mpi, 1, CON, mp.get_bvect_cart()),
141  enerps_euler(mpi),
142  uuu2(mpi),
143  gam_euler2(mpi),
144  delta_car(mpi)
145 {
146  // All the matter quantities are initialized to zero :
147  nbar2 = 0 ;
148  ent2 = 0 ;
149  K_nn = 0 ;
150  K_np = 0 ;
151  K_pp = 0 ;
152  alpha_eos = 0.;
153  sphph_euler = 0;
154  j_euler = 0;
155  j_euler1 = 0 ;
156  j_euler2 = 0;
157  enerps_euler = 0;
159 
160  // Initialization to a static state :
161  omega2 = 0 ;
162  uuu2 = 0 ;
163  gam_euler2 = 1 ;
164  delta_car = 0 ;
165 
166  // Pointers of derived quantities initialized to zero :
167  set_der_0x0() ;
168 
169 }
170 
171 // Copy constructor
172 // ----------------
173 
175  Etoile_rot(et),
176  eos(et.eos),
177  ent2(et.ent2),
178  nbar2(et.nbar2),
179  K_nn(et.K_nn),
180  K_np(et.K_np),
181  K_pp(et.K_pp),
182  alpha_eos(et.alpha_eos),
183  sphph_euler(et.sphph_euler),
184  j_euler(et.j_euler),
185  j_euler1(et.j_euler1),
186  j_euler2(et.j_euler2),
187  enerps_euler(et.enerps_euler),
188  uuu2(et.uuu2),
189  gam_euler2(et.gam_euler2),
190  delta_car(et.delta_car)
191 {
192  omega2 = et.omega2 ;
193 
194  // Pointers of derived quantities initialized to zero :
195  set_der_0x0() ;
196 }
197 
198 
199 // Constructor from a file
200 // ------------------------
201 Et_rot_bifluid::Et_rot_bifluid(Map& mpi, const Eos_bifluid& eos_i, FILE* fich):
202  Etoile_rot(mpi, *eos_i.trans2Eos(), fich),
203  eos(eos_i),
204  ent2(mpi),
205  nbar2(mpi),
206  K_nn(mpi),
207  K_np(mpi),
208  K_pp(mpi),
209  alpha_eos(mpi),
210  sphph_euler(mpi),
211  j_euler(mpi, 1, CON, mp.get_bvect_cart()),
212  j_euler1(mpi, 1, CON, mp.get_bvect_cart()),
213  j_euler2(mpi, 1, CON, mp.get_bvect_cart()),
214  enerps_euler(mpi),
215  uuu2(mpi),
216  gam_euler2(mpi),
217  delta_car(mpi)
218 {
219 
220  // Etoile parameters
221  // -----------------
222  // omega2 is read in the file:
223  fread_be(&omega2, sizeof(double), 1, fich) ;
224 
225 
226  // Read of the saved fields:
227  // ------------------------
228 
229  Tenseur ent2_file(mp, fich) ;
230  ent2 = ent2_file ;
231 
232  // All other fields are initialized to zero :
233  // ----------------------------------------
234  uuu2 = 0 ;
235  delta_car = 0 ;
236 
237  // Pointers of derived quantities initialized to zero
238  // --------------------------------------------------
239  set_der_0x0() ;
240 
241 }
242 
243  //------------//
244  // Destructor //
245  //------------//
246 
248 
249  del_deriv() ;
250 
251 }
252 
253  //----------------------------------//
254  // Management of derived quantities //
255  //----------------------------------//
256 
258 
260 
261  if (p_ray_eq2 != 0x0) delete p_ray_eq2 ;
262  if (p_ray_eq2_pis2 != 0x0) delete p_ray_eq2_pis2 ;
263  if (p_ray_eq2_pi != 0x0) delete p_ray_eq2_pi ;
264  if (p_ray_pole2 != 0x0) delete p_ray_pole2 ;
265  if (p_l_surf2 != 0x0) delete p_l_surf2 ;
266  if (p_xi_surf2 != 0x0) delete p_xi_surf2 ;
267  if (p_r_circ2 != 0x0) delete p_r_circ2 ;
268  if (p_area2 != 0x0) delete p_area2 ;
269  if (p_aplat2 != 0x0) delete p_aplat2 ;
270  if (p_mass_b1 != 0x0) delete p_mass_b1 ;
271  if (p_mass_b2 != 0x0) delete p_mass_b2 ;
272  if (p_angu_mom_1 != 0x0) delete p_angu_mom_1 ;
273  if (p_angu_mom_2 != 0x0) delete p_angu_mom_2 ;
274  if (p_coupling_mominert_1 != 0x0) delete p_coupling_mominert_1 ;
275  if (p_coupling_mominert_2 != 0x0) delete p_coupling_mominert_2 ;
276  if (p_coupling_entr != 0x0) delete p_coupling_entr ;
277  if (p_coupling_LT_1 != 0x0) delete p_coupling_LT_1 ;
278  if (p_coupling_LT_2 != 0x0) delete p_coupling_LT_2 ;
279 
280  set_der_0x0() ;
281 }
282 
283 
284 
285 
287 
289 
290  p_ray_eq2 = 0x0 ;
291  p_ray_eq2_pis2 = 0x0 ;
292  p_ray_eq2_pi = 0x0 ;
293  p_ray_pole2 = 0x0 ;
294  p_l_surf2 = 0x0 ;
295  p_xi_surf2 = 0x0 ;
296  p_r_circ2 = 0x0 ;
297  p_area2 = 0x0 ;
298  p_aplat2 = 0x0 ;
299  p_mass_b1 = 0x0;
300  p_mass_b2 = 0x0;
301  p_angu_mom_1 = 0x0;
302  p_angu_mom_2 = 0x0;
303  p_coupling_mominert_1 = 0x0;
304  p_coupling_mominert_2 = 0x0;
305  p_coupling_entr = 0x0;
306  p_coupling_LT_1 = 0x0;
307  p_coupling_LT_2 = 0x0;
308 
309 }
310 
312 
326  del_deriv() ;
327 }
328 
329 // Assignment of the enthalpy field
330 // --------------------------------
331 
332 void Et_rot_bifluid::set_enthalpies(const Cmp& ent_i, const Cmp& ent2_i) {
333 
334  ent = ent_i ;
335  ent2 = ent2_i ;
336 
337  // Update of (nbar, ener, press) :
338  equation_of_state() ;
339 
340  // The derived quantities are obsolete:
341  del_deriv() ;
342 
343 }
344 
345 
346  //--------------//
347  // Assignment //
348  //--------------//
349 
350 // Assignment to another Et_rot_bifluid
351 // --------------------------------
353 
354  // Assignment of quantities common to all the derived classes of Etoile
355  Etoile_rot::operator=(et) ;
356 
357  assert( &(et.eos) == &eos ) ; // Same EOS
358  // Assignement of proper quantities of class Et_rot_bifluid
359  omega2 = et.omega2 ;
360 
361  ent2 = et.ent2 ;
362  nbar2 = et.nbar2 ;
363  K_nn = et.K_nn ;
364  K_np = et.K_np ;
365  K_pp = et.K_pp ;
366  alpha_eos = et.alpha_eos ;
368  j_euler = et.j_euler;
369  j_euler1 = et.j_euler1;
370  j_euler2 = et.j_euler2;
372  uuu2 = et.uuu2 ;
373  gam_euler2 = et.gam_euler2 ;
374  delta_car = et.delta_car ;
375 
376  del_deriv() ; // Deletes all derived quantities
377 
378 }
379 
380  //--------------//
381  // Outputs //
382  //--------------//
383 
384 // Save in a file
385 // --------------
386 void Et_rot_bifluid::sauve(FILE* fich) const {
387 
388  Etoile_rot::sauve(fich) ;
389 
390  fwrite_be(&omega2, sizeof(double), 1, fich) ;
391 
392  ent2.sauve(fich) ;
393 
394 
395 }
396 
397 // Printing
398 // --------
399 
400 ostream& Et_rot_bifluid::operator>>(ostream& ost) const {
401 
402  using namespace Unites ;
403 
404  Etoile::operator>>(ost) ;
405 
406  ost << endl ;
407  ost << "Bifluid rotating star" << endl ;
408  ost << "-------------" << endl ;
409  ost << setprecision(16);
410  double freq = omega / (2.*M_PI) ;
411  ost << "Omega1 : " << omega * f_unit
412  << " rad/s f : " << freq * f_unit << " Hz" << endl ;
413  ost << "Rotation period 1: " << 1000. / (freq * f_unit) << " ms"
414  << endl ;
415 
416  double freq2 = omega2 / (2.*M_PI) ;
417  ost << "Omega2 : " << omega2 * f_unit
418  << " rad/s f : " << freq2 * f_unit << " Hz" << endl ;
419  ost << "Rotation period 2: " << 1000. / (freq2 * f_unit) << " ms"
420  << endl ;
421 
422 
423  ost << "Total angular momentum J : "
424  << angu_mom()/( qpig / (4* M_PI) * msol*msol) << " G M_sol^2 / c"
425  << endl ;
426  ost << "c J / (G M^2) : "
427  << angu_mom()/( qpig / (4* M_PI) * pow(mass_g(), 2.) ) << endl ;
428 
429  double mom_iner = fabs(angu_mom() / omega2) ;
430  ost << "Total moment of inertia I = J/Omega2 : "
431  << mom_iner << " Lorene units"
432  << endl;
433  double mom_iner_38si = mom_iner * rho_unit * (pow(r_unit, double(5.)) / double(1.e38) ) ;
434  ost << "Total moment of inertia I = J/Omega2 : "
435  << mom_iner_38si << " 10^38 kg m^2"
436  << endl ;
437 
438  // partial angular momenta
439  ost << "Angular momentum J of fluid 1 : "
440  << angu_mom_1()/( qpig / (4* M_PI) * msol*msol) << " G M_sol^2 / c"
441  << endl ;
442  ost << "Angular momentum J of fluid 2 : "
443  << angu_mom_2()/( qpig / (4* M_PI) * msol*msol) << " G M_sol^2 / c"
444  << endl ;
445 
446 
447  // partial moments of inertia defined as Jn/Omega and Jp/Omega
448  // Remark : such definitions only makes sense in corotation
449  ost.precision(16) ;
450  if (omega == omega2) {
451 
452  double mom_iner_1 = 0.; //< In
453  double mom_iner_2 = 0.; //< Ip
454 
455  if (omega != __infinity) {
456 
457  ost << "Partial moments of inertia (defined in corotation only) :" << endl;
458 
459  mom_iner_1 = fabs(angu_mom_1() / omega) ;
460  ost << "Moment of inertia of fluid 1 : " << mom_iner_1 << " Lorene units " << endl ;
461  double mom_iner_1_38si = mom_iner_1 * rho_unit * (pow(r_unit, double(5.)) / double(1.e38) ) ;
462  ost << "Moment of inertia of fluid 1 : " << mom_iner_1_38si << " 10^38 kg m^2"
463  << endl ;
464 
465  mom_iner_2 = fabs(angu_mom_2() / omega) ;
466  ost << "Moment of inertia of fluid 2 : " << mom_iner_2 << " Lorene units " << endl ;
467  double mom_iner_2_38si = mom_iner_2 * rho_unit * (pow(r_unit, double(5.)) / double(1.e38) ) ;
468  ost << "Moment of inertia of fluid 2 : " << mom_iner_2_38si << " 10^38 kg m^2"
469  << endl;
470  }
471  }
472 
473  ost << "***** Fluid coupling quantities *****" << endl;
474 
475  double mominert_1_tilde_si = coupling_mominert_1() * rho_unit * (pow(r_unit, double(5.)) / double(1.e38) ) ;
476  double mominert_2_tilde_si = coupling_mominert_2() * rho_unit * (pow(r_unit, double(5.)) / double(1.e38) ) ;
477 
478  ost << "tilde{I}_n : " << coupling_mominert_1() << " SI: " << mominert_1_tilde_si << " 10^38 kg m^2"
479  << endl;
480  ost << "tilde{I}_p : " << coupling_mominert_2() << " SI: " << mominert_2_tilde_si << " 10^38 kg m^2"
481  << endl;
482  ost << "tilde{varepsilon}_n : " << coupling_entr() / coupling_mominert_1() << endl;
483  ost << "tilde{varepsilon}_p : " << coupling_entr() / coupling_mominert_2() << endl;
484  ost << "tilde{omega}_n : " << coupling_LT_1() / coupling_mominert_1() << endl;
485  ost << "tilde{omega}_p : " << coupling_LT_2() / coupling_mominert_2() << endl;
486  ost << " Verif : Jn = " << coupling_mominert_1() * omega - coupling_LT_1() + coupling_entr() * (omega2 - omega)
487  << " Jp = " << coupling_mominert_2() * omega2 - coupling_LT_2() + coupling_entr() * (omega - omega2)
488  << endl ;
489  ost << " Num : Jn = " << angu_mom_1() << " Jp = " << angu_mom_2()
490  << endl;
491 
492  double nphi_c = nphi()(0, 0, 0, 0) ;
493  if ( (omega==0) && (nphi_c==0) ) {
494  ost << "Central N^phi : " << nphi_c << endl ;
495  }
496  else{
497  ost << "Central N^phi/Omega : " << nphi_c / omega << endl ;
498  }
499  if (omega2!=0)
500  ost << "Central N^phi/Omega2 : " << nphi_c / omega2 << endl ;
501 
502  ost << "Error on the virial identity GRV2 : " << endl ;
503  ost << "GRV2 = " << grv2() << endl ;
504  ost << "Error on the virial identity GRV3 : " << endl ;
505  double xgrv3 = grv3(&ost) ;
506  ost << "GRV3 = " << xgrv3 << endl ;
507 
508  ost << "Circumferential equatorial radius R_circ : "
509  << r_circ()/km << " km" << endl ;
510  ost << "Mean radius R_mean : "
511  << mean_radius()/km << " km" << endl ;
512  ost << "Coordinate equatorial radius r_eq : " << ray_eq()/km << " km"
513  << endl ;
514  ost << "Flattening r_pole/r_eq : " << aplat() << endl ;
515  ost << "Circumferential equatorial radius R_circ2 : "
516  << r_circ2()/km << " km" << endl ;
517  ost << "Mean radius R_mean2 : "
518  << mean_radius2()/km << " km" << endl ;
519  ost << "Coordinate equatorial radius r_eq2 : " << ray_eq2()/km << " km"
520  << endl ;
521  ost << "Flattening r_pole2/r_eq2 : " << aplat2() << endl ;
522 
523  int lsurf = nzet - 1;
524  int nt = mp.get_mg()->get_nt(lsurf) ;
525  int nr = mp.get_mg()->get_nr(lsurf) ;
526  ost << "Equatorial value of the velocity U: "
527  << uuu()(nzet-1, 0, nt-1, nr-1) << " c" << endl ;
528  ost << "Equatorial value of the velocity U2: "
529  << uuu2()(nzet-1, 0, nt-1, nr-1) << " c" << endl ;
530  ost << "Redshift at the equator (forward) : " << z_eqf() << endl ;
531  ost << "Redshift at the equator (backward): " << z_eqb() << endl ;
532  ost << "Redshift at the pole : " << z_pole() << endl ;
533 
534 
535  ost << "Central value of log(N) : "
536  << logn()(0, 0, 0, 0) << endl ;
537 
538  ost << "Central value of dzeta=log(AN) : "
539  << dzeta()(0, 0, 0, 0) << endl ;
540 
541  ost << "Central value of B^2 : " << b_car()(0,0,0,0) << endl ;
542 
543  Tbl diff_a_b = diffrel( a_car(), b_car() ) ;
544  ost <<
545  "Relative discrepancy in each domain between the metric coef. A^2 and B^2 : "
546  << endl ;
547  for (int l=0; l<diff_a_b.get_taille(); l++) {
548  ost << diff_a_b(l) << " " ;
549  }
550  ost << endl;
551  ost << "Quadrupole moment : " << mom_quad() << endl ;
552  double mom_quad_38si = mom_quad() * rho_unit * (pow(r_unit, double(5.)) / double(1.e38) ) ;
553  ost << "Quadrupole moment Q : " << mom_quad_38si << " 10^38 kg m^2"
554  << endl ;
555  ost << "Old quadrupole moment : " << mom_quad_old() << endl ;
556  ost << "Coefficient b : " << mom_quad_Bo() / pow(mass_g(), 2.) << endl ;
557  ost << "q = c^4 Q / (G^2 M^3) : "
558  << mom_quad() / ( ggrav * ggrav * pow(mass_g(), 3.) ) << endl ;
559  ost << "j = c J / (G M^2) : "
560  << angu_mom()/( ggrav * pow(mass_g(), 2.) ) << endl ;
561 
562 
563  ost << "Baryon mass 1 : " << mass_b1() / msol << " Msol" << endl ;
564  ost << "Baryon mass 2 : " << mass_b2() / msol << " Msol" << endl ;
565 
566  ost << endl ;
567 
568  return ost ;
569 
570 }
571 
572 
573 void Et_rot_bifluid::partial_display(ostream& ost) const {
574 
575  using namespace Unites ;
576 
577  double freq = omega / (2.*M_PI) ;
578  ost << "Omega : " << omega * f_unit
579  << " rad/s f : " << freq * f_unit << " Hz" << endl ;
580  ost << "Rotation period : " << 1000. / (freq * f_unit) << " ms"
581  << endl ;
582  ost << endl << "Central enthalpy : " << ent()(0,0,0,0) << " c^2" << endl ;
583  ost << "Central proper baryon density : " << nbar()(0,0,0,0)
584  << " x 0.1 fm^-3" << endl ;
585  double freq2 = omega2 / (2.*M_PI) ;
586  ost << "Omega2 : " << omega2 * f_unit
587  << " rad/s f : " << freq2 * f_unit << " Hz" << endl ;
588  ost << "Rotation period 2: " << 1000. / (freq2 * f_unit) << " ms"
589  << endl ;
590  ost << endl << "Central enthalpy 2: " << ent2()(0,0,0,0) << " c^2" << endl ;
591  ost << "Central proper baryon density 2: " << nbar2()(0,0,0,0)
592  << " x 0.1 fm^-3" << endl ;
593  ost << "Central proper energy density : " << ener()(0,0,0,0)
594  << " rho_nuc c^2" << endl ;
595  ost << "Central pressure : " << press()(0,0,0,0)
596  << " rho_nuc c^2" << endl ;
597 
598  ost << "Central value of log(N) : "
599  << logn()(0, 0, 0, 0) << endl ;
600  ost << "Central lapse N : " << nnn()(0,0,0,0) << endl ;
601  ost << "Central value of dzeta=log(AN) : "
602  << dzeta()(0, 0, 0, 0) << endl ;
603  ost << "Central value of A^2 : " << a_car()(0,0,0,0) << endl ;
604  ost << "Central value of B^2 : " << b_car()(0,0,0,0) << endl ;
605 
606  double nphi_c = nphi()(0, 0, 0, 0) ;
607  if ( (omega==0) && (nphi_c==0) ) {
608  ost << "Central N^phi : " << nphi_c << endl ;
609  }
610  else{
611  ost << "Central N^phi/Omega : " << nphi_c / omega << endl ;
612  }
613 
614 
615  int lsurf = nzet - 1;
616  int nt = mp.get_mg()->get_nt(lsurf) ;
617  int nr = mp.get_mg()->get_nr(lsurf) ;
618  ost << "Equatorial value of the velocity U: "
619  << uuu()(nzet-1, 0, nt-1, nr-1) << " c" << endl ;
620 
621  ost << "Equatorial value of the velocity U2: "
622  << uuu2()(nzet-1, 0, nt-1, nr-1) << " c" << endl ;
623 
624  ost << endl
625  << "Coordinate equatorial radius r_eq = "
626  << ray_eq()/km << " km" << endl ;
627  ost << "Flattening r_pole/r_eq : " << aplat() << endl ;
628  ost << endl
629  << "Coordinate equatorial radius r_eq2 = "
630  << ray_eq2()/km << " km" << endl ;
631  ost << "Flattening r_pole2/r_eq2 : " << aplat2() << endl ;
632 
633 }
634 
635 //
636 // Equation of state
637 //
638 
640 
641  Cmp ent_eos = ent() ;
642  Cmp ent2_eos = ent2() ;
643  Tenseur rel_vel(delta_car) ;
644 
645  if (nzet > 1) {
646  // Slight rescale of the enthalpy field in case of 2 domains inside the
647  // star
648 
649  if (nzet > 2) {
650  cout << "Et_rot_bifluid::equation_of_state: not ready yet for nzet > 2 !" << endl ;
651  }
652 
653  double epsilon = 1.e-12 ;
654 
655  const Mg3d* mg = mp.get_mg() ;
656  int nz = mg->get_nzone() ;
657 
658  Mtbl xi(mg) ;
659  xi.set_etat_qcq() ;
660  for (int l=0; l<nz; l++) {
661  xi.t[l]->set_etat_qcq() ;
662  for (int k=0; k<mg->get_np(l); k++) {
663  for (int j=0; j<mg->get_nt(l); j++) {
664  for (int i=0; i<mg->get_nr(l); i++) {
665  xi.set(l,k,j,i) =
666  mg->get_grille3d(l)->x[i] ;
667  }
668  }
669  }
670 
671  }
672 
673  Cmp fact_ent(mp) ;
674  fact_ent.allocate_all() ;
675 
676  fact_ent.set(0) = 1 + epsilon * xi(0) * xi(0) ;
677  fact_ent.set(1) = 1 - 0.25 * epsilon * (xi(1) - 1) * (xi(1) - 1) ;
678 
679  for (int l=nzet; l<nz; l++) {
680  fact_ent.set(l) = 1 ;
681  }
682 
683  ent_eos = fact_ent * ent_eos ;
684  ent2_eos = fact_ent * ent2_eos ;
685  ent_eos.std_base_scal() ;
686  ent2_eos.std_base_scal() ;
687 
688  } // if nzet > 1
689 
690 
691  // Call to the EOS
692  nbar.set_etat_qcq() ;
693  nbar2.set_etat_qcq() ;
694  ener.set_etat_qcq() ;
695  press.set_etat_qcq() ;
696  K_nn.set_etat_qcq() ;
697  K_np.set_etat_qcq() ;
698  K_pp.set_etat_qcq() ;
700 
701  const Eos_bf_tabul* eos_t = dynamic_cast<const Eos_bf_tabul*>(&eos) ;
702 
703  if (eos_t != 0x0) { // The EoS is tabulated
704  eos_t->calcule_interpol(ent_eos, ent2_eos, rel_vel(), nbar.set(), nbar2.set(),
705  ener.set(), press.set(), K_nn.set(), K_np.set(), K_pp.set(), alpha_eos.set(),
706  nzet) ;
707  }
708  else { // The EoS is analytic
709  eos.calcule_tout(ent_eos, ent2_eos, rel_vel(), nbar.set(), nbar2.set(),
710  ener.set(), press.set(), nzet) ;
711 
712  K_nn.set() = eos.get_Knn(nbar(), nbar2(), delta_car(), nzet);
713  K_pp.set() = eos.get_Kpp(nbar(), nbar2(), delta_car(), nzet);
714  K_np.set() = eos.get_Knp(nbar(), nbar2(), delta_car(), nzet);
715  alpha_eos.set() = eos.get_Knp(nbar(), nbar2(), delta_car(), nzet) * nbar() * nbar2() * pow(1. - unsurc2 * rel_vel(), -1.5) / 2. ;
716  }
717  // Set the bases for spectral expansion
718  nbar.set_std_base() ;
719  nbar2.set_std_base() ;
720  ener.set_std_base() ;
721  press.set_std_base() ;
722  K_pp.set_std_base() ;
723  K_nn.set_std_base() ;
724  K_np.set_std_base() ;
726  // The derived quantities are obsolete
727  del_deriv() ;
728 
729 }
730 
731 //
732 // Computation of hydro quantities
733 //
734 
736 
737  const Mg3d* mg = mp.get_mg();
738  int nz = mg->get_nzone() ;
739  int nzm1 = nz - 1 ;
740 
741  // RP: I prefer to use the 3-vector j_euler instead of u_euler
742  // for better physical "encapsulation"
743  // (i.e. --> use same form of Poisson-equations for all etoile sub-classes!)
744  u_euler.set_etat_nondef(); // make sure it's not used
745 
746  // (Norm of) Euler-velocity of the first fluid
747  //------------------------------
748  uuu.set_etat_qcq();
749 
750  uuu.set() = bbb() * (omega - nphi() ) / nnn();
751  uuu.annule(nzm1) ;
752 
753  // gosh, we have to exclude the thing being zero here... :(
754  if( uuu.get_etat() != ETATZERO )
755  {
756  (uuu.set()).std_base_scal() ;
757  (uuu.set()).mult_rsint();
758  }
759  uuu.set_std_base();
760 
761 
762  // (Norm of) Euler-velocity of the second fluid
763  //----------------------------------------
764  uuu2.set_etat_qcq();
765 
766  uuu2.set() = bbb() * (omega2 - nphi() ) / nnn();
767  uuu2.annule(nzm1) ;
768 
769  if( uuu2.get_etat() != ETATZERO )
770  {
771  (uuu2.set()).std_base_scal();
772  (uuu2.set()).mult_rsint();
773  }
774  uuu2.set_std_base();
775 
776  // Sanity check:
777  // Is one of the new velocities larger than c in the equatorial plane ?
778  //----------------------------------------
779 
780  // Index of the point at phi=0, theta=pi/2 at the surface of the star:
781  int l_b = nzet - 1 ;
782  int j_b = mg->get_nt(l_b) - 1 ;
783 
784  bool superlum = false ;
785  if (relativistic) {
786  for (int l=0; l<nzet; l++) {
787  for (int i=0; i<mg->get_nr(l); i++) {
788 
789  double u1 = uuu()(l, 0, j_b, i) ;
790  double u2 = uuu2()(l, 0, j_b, i) ;
791  if ((u1 >= 1.) || (u2>=1.)) { // superluminal velocity
792  superlum = true ;
793  cout << "U > c for l, i : " << l << " " << i
794  << " U1 = " << u1 << endl ;
795  cout << " U2 = " << u2 << endl ;
796  }
797  }
798  }
799  if ( superlum ) {
800  cout << "**** VELOCITY OF LIGHT REACHED ****" << endl ;
801  abort() ;
802  }
803  }
804 
805 
806  Tenseur uuu_car = uuu * uuu;
807  Tenseur uuu2_car = uuu2 * uuu2;
808 
809 
810  // Lorentz factors
811  // --------------
812  Tenseur gam_car = 1.0 / (1.0 - unsurc2 * uuu_car) ;
813  gam_euler = sqrt(gam_car) ;
814  gam_euler.set_std_base() ; // sets the standard spectral bases for a scalar field
815 
816 
817  Tenseur gam2_car = 1.0 / (1.0 - unsurc2 * uuu2_car) ;
818  gam_euler2 = sqrt(gam2_car) ;
820 
821  // Update of "relative velocity" squared: $\Delta^2$
822  // ---------------------------
823 
824  delta_car = (uuu - uuu2)*(uuu - uuu2) / ( (1 - unsurc2* uuu*uuu2) *(1 - unsurc2* uuu*uuu2 ) ) ;
826 
827  Tenseur Ann(ent) ;
828  Tenseur Anp(ent) ;
829  Tenseur App(ent) ;
830 
831  Ann = gam_car() * nbar() * nbar() * K_nn() ;
832  Anp = gam_euler() * gam_euler2() * nbar() * nbar2() * K_np();
833  App = gam2_car() * nbar2() * nbar2() * K_pp();
834 
835 
836  // Energy density E with respect to the Eulerian observer
837  //------------------------------------
838  // use of ener_euler is deprecated, because it's useless in Newtonian limit!
839  ener_euler.set_etat_nondef(); // make sure, it's not used
840 
841  Tenseur E_euler(mp);
842  E_euler = Ann + 2. * Anp + App - press ;
843  E_euler.set_std_base() ;
844 
845 
846  // S^phi_phi component of stress-tensor S^i_j
847  //------------------------------------
848  sphph_euler = press() + Ann() * uuu_car() + 2. * Anp() * uuu() * uuu2() + App() * uuu2_car();
850 
851 
852  // Trace of the stress tensor with respect to the Eulerian observer
853  //------------------------------------
854  s_euler = 2. * press() + sphph_euler();
855  s_euler.set_std_base() ;
856 
857  // The combination enerps_euler := (E + S_i^i) which has Newtonian limit -> rho
858  if (relativistic)
859  enerps_euler = E_euler + s_euler;
860  else
861  enerps_euler = eos.get_m1() * nbar() + eos.get_m2() * nbar2();
862 
863 
864  // the (flat-space) angular-momentum 3-vector j_euler^i
865  //-----------------------------------
866  Tenseur Jph(mp); // the normalized phi-component of J^i: Sqrt[g_phiphi]*J^phi
867  Jph = Ann*uuu + Anp*(uuu + uuu2) + App*uuu2 ;
868 
870 
871  j_euler.set(0) = 0; // r tetrad component
872  j_euler.set(1) = 0; // theta tetrad component
873  j_euler.set(2) = Jph()/ bbb(); // phi tetrad component ... = J^phi r sin(th)
876 
877  // RP: it seems that j_euler _HAS_ to have cartesian triad set on exit from here...!!
878  j_euler.change_triad( mp.get_bvect_cart() ) ; // Triad = Cartesian triad
879 
880  if( (j_euler(0).get_etat() == ETATZERO)&&(j_euler(1).get_etat() == ETATZERO)&&(j_euler(2).get_etat()==ETATZERO))
881  j_euler = 0;
882 
883  // the (flat-space) angular-momentum 3-vector j_euler^i for fluid 1
884  //-----------------------------------
885  Tenseur Jph1(mp); // the normalized phi-component of J_n^i: Sqrt[g_phiphi]*J_n^phi
886  Jph1 = Ann*uuu + Anp*uuu2 ;
887 
889 
890  j_euler1.set(0) = 0;
891  j_euler1.set(1) = 0;
892  j_euler1.set(2) = Jph1()/ bbb();
895 
897 
898  if( (j_euler1(0).get_etat() == ETATZERO)&&(j_euler1(1).get_etat() == ETATZERO)&&(j_euler1(2).get_etat()==ETATZERO))
899  j_euler1 = 0;
900 
901  // the (flat-space) angular-momentum 3-vector j_euler^i for fluid 2
902  //-----------------------------------
903  Tenseur Jph2(mp); // the normalized phi-component of J_p^i: Sqrt[g_phiphi]*J_p^phi
904  Jph2 = Anp*uuu + App*uuu2 ;
905 
907 
908  j_euler2.set(0) = 0;
909  j_euler2.set(1) = 0;
910  j_euler2.set(2) = Jph2()/ bbb();
913 
914  j_euler2.change_triad( mp.get_bvect_cart() ) ; // Triad = Cartesian triad
915 
916  if( (j_euler2(0).get_etat() == ETATZERO)&&(j_euler2(1).get_etat() == ETATZERO)&&(j_euler2(2).get_etat()==ETATZERO))
917  j_euler2 = 0;
918 
919  // The derived quantities are obsolete
920  // -----------------------------------
921  del_deriv() ;
922 
923 } // hydro_euler()
924 
925 }
926 
virtual double z_eqf() const
Forward redshift factor at equator.
Cmp get_Kpp(const Cmp &nbar1, const Cmp &nbar2, const Cmp &x2, int nzet, int l_min=0) const
Computes the derivatives of the energy/(baryonic density 2) .
Definition: eos_bifluid.C:772
double get_m2() const
Return the individual particule mass .
Definition: eos_bifluid.h:269
double * p_mass_b2
Baryon mass of fluid 2.
virtual double r_circ() const
Circumferential equatorial radius.
const Grille3d * get_grille3d(int l) const
Returns a pointer on the 3D mono-grid for domain no. l.
Definition: grilles.h:517
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition: etoile.C:514
void annule(int l)
Sets the Tenseur to zero in a given domain.
Definition: tenseur.C:916
void calcule_interpol(const Cmp &ent1, const Cmp &ent2, const Cmp &delta2, Cmp &nbar1, Cmp &nbar2, Cmp &ener, Cmp &press, Cmp &K_nn, Cmp &K_np, Cmp &K_pp, Cmp &alpha_eos, int nzet, int l_min=0) const
General computational method for Cmp &#39;s, it computes both baryon densities, energy and pressure profi...
Definition: eos_bf_tabul.C:757
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition: cmp.h:446
virtual double mass_g() const
Gravitational mass.
Tenseur j_euler2
To compute .
Tenseur j_euler1
To compute .
Tenseur K_pp
Coefficient .
void set_triad(const Base_vect &new_triad)
Assigns a new vectorial basis (triad) of decomposition.
Definition: tenseur.C:690
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
Definition: grilles.h:479
Cmp sqrt(const Cmp &)
Square root.
Definition: cmp_math.C:223
double mass_b1() const
Baryon mass of fluid 1.
void set_enthalpies(const Cmp &, const Cmp &)
Sets both enthalpy profiles.
void set_std_base()
Set the standard spectal basis of decomposition for each component.
Definition: tenseur.C:1186
Multi-domain array.
Definition: mtbl.h:118
const Base_vect_spher & get_bvect_spher() const
Returns the orthonormal vectorial basis associated with the coordinates of the mapping.
Definition: map.h:801
double * p_r_circ2
Circumferential radius of fluid no.2.
Tbl * p_xi_surf2
Description of the surface of fluid 2: 2-D Tbl containing the values of the radial coordinate on the...
Lorene prototypes.
Definition: app_hor.h:67
Standard units of space, time and mass.
Tenseur uuu2
Norm of the (fluid no.2) 3-velocity with respect to the eulerian observer.
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:783
Tenseur nnn
Total lapse function.
Definition: etoile.h:512
virtual void equation_of_state()
Computes the proper baryon and energy densities, as well as pressure and the coefficients Knn...
Tenseur nphi
Metric coefficient .
Definition: etoile.h:1513
Tenseur s_euler
Trace of the stress tensor in the Eulerian frame.
Definition: etoile.h:471
double * p_angu_mom_1
Angular momentum of fluid 1.
double ray_eq() const
Coordinate radius at , [r_unit].
Base class for coordinate mappings.
Definition: map.h:688
Class for isolated rotating stars *** DEPRECATED : use class Star_rot instead ***.
Definition: etoile.h:1499
Tenseur b_car
Square of the metric factor B.
Definition: etoile.h:1510
void operator=(const Etoile_rot &)
Assignment to another Etoile_rot.
Definition: etoile_rot.C:415
virtual void sauve(FILE *) const
Save in a file.
double unsurc2
: unsurc2=1 for a relativistic star, 0 for a Newtonian one.
Definition: etoile.h:445
double * p_ray_eq2
Coordinate radius at , .
Tenseur press
Fluid pressure.
Definition: etoile.h:464
double * x
Array of values of at the nr collocation points.
Definition: grilles.h:215
Class for two-fluid rotating relativistic stars.
virtual double mean_radius() const
Mean radius.
virtual double mom_quad_Bo() const
Part of the quadrupole moment.
Tbl diffrel(const Cmp &a, const Cmp &b)
Relative difference between two Cmp (norme version).
Definition: cmp_math.C:507
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: etoile_rot.C:344
2-fluids equation of state base class.
Definition: eos_bifluid.h:180
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: tbl.C:364
Tenseur j_euler
Total angular momentum (flat-space!) 3-vector , which is related to of the "3+1" decomposition...
Cmp get_Knn(const Cmp &nbar1, const Cmp &nbar2, const Cmp &x2, int nzet, int l_min=0) const
Computes the derivatives of the energy/(baryonic density 1) .
Definition: eos_bifluid.C:750
double * p_mass_b1
Baryon mass of fluid 1.
Tenseur gam_euler2
Lorentz factor between the fluid 2 and Eulerian observers.
Cmp & set()
Read/write for a scalar (see also operator=(const Cmp&) ).
Definition: tenseur.C:840
Tenseur u_euler
Fluid 3-velocity with respect to the Eulerian observer.
Definition: etoile.h:477
void change_triad(const Base_vect &new_triad)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
Definition: tenseur.C:684
void sauve(FILE *) const
Save in a file.
Definition: tenseur.C:1341
const Eos_bifluid & eos
Equation of state for two-fluids model.
virtual double z_pole() const
Redshift factor at North pole.
Tenseur nbar
Baryon density in the fluid frame.
Definition: etoile.h:462
double * p_coupling_mominert_1
Quantities used to describe the different couplings between the fluids.
double omega2
Rotation angular velocity for fluid 2 ([f_unit] )
virtual void sauve(FILE *) const
Save in a file.
Definition: etoile_rot.C:452
Tenseur gam_euler
Lorentz factor between the fluid and Eulerian observers.
Definition: etoile.h:474
virtual double grv3(ostream *ost=0x0) const
Error on the virial identity GRV3.
virtual double angu_mom_2() const
Angular momentum of fluid 2.
Map & mp
Mapping associated with the star.
Definition: etoile.h:432
double * p_area2
Surface area of fluid no.2.
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:465
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
virtual double z_eqb() const
Backward redshift factor at equator.
int get_etat() const
Returns the logical state.
Definition: tenseur.h:710
Tenseur enerps_euler
the combination : useful because in the Newtonian limit .
Tenseur bbb
Metric factor B.
Definition: etoile.h:1507
double * p_ray_eq2_pi
Coordinate radius at , .
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:73
Tenseur K_np
Coefficient .
Cmp pow(const Cmp &, int)
Power .
Definition: cmp_math.C:351
void operator=(const Et_rot_bifluid &)
Assignment to another Et_rot_bifluid.
virtual double mom_quad_old() const
Part of the quadrupole moment.
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: mtbl.C:302
virtual double grv2() const
Error on the virial identity GRV2.
Tenseur uuu
Norm of u_euler.
Definition: etoile.h:1521
Tbl & set(int l)
Read/write of the Tbl in a given domain.
Definition: mtbl.h:221
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
double omega
Rotation angular velocity ([f_unit] )
Definition: etoile.h:1504
virtual double mean_radius2() const
Mean radius for fluid 2.
void std_base_scal()
Sets the spectral bases of the Valeur va to the standard ones for a scalar.
Definition: cmp.C:647
int nzet
Number of domains of *mp occupied by the star.
Definition: etoile.h:435
virtual void del_deriv() const
Deletes all the derived quantities.
virtual void partial_display(ostream &) const
Printing of some informations, excluding all global quantities.
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:72
int get_nr(int l) const
Returns the number of points in the radial direction ( ) in domain no. l.
Definition: grilles.h:469
virtual double angu_mom_1() const
Angular momentum of fluid 1.
Tenseur a_car
Total conformal factor .
Definition: etoile.h:518
double * p_angu_mom_2
Angular momentum of fluid 2.
Tenseur ent2
Log-enthalpy for the second fluid.
void allocate_all()
Sets the logical state to ETATQCQ (ordinary state) and performs the memory allocation of all the elem...
Definition: cmp.C:326
bool relativistic
Indicator of relativity: true for a relativistic star, false for a Newtonian one. ...
Definition: etoile.h:440
Multi-domain grid.
Definition: grilles.h:279
Tenseur sphph_euler
The component of the stress tensor .
virtual void del_hydro_euler()
Sets to ETATNONDEF (undefined state) the hydrodynamical quantities relative to the Eulerian observer...
Tbl & set(int l)
Read/write of the value in a given domain.
Definition: cmp.h:724
virtual double mom_quad() const
Quadrupole moment.
Tenseur ener
Total energy density in the fluid frame.
Definition: etoile.h:463
const Base_vect_cart & get_bvect_cart() const
Returns the Cartesian basis associated with the coordinates (x,y,z) of the mapping, i.e.
Definition: map.h:809
virtual void calcule_tout(const Cmp &ent1, const Cmp &ent2, const Cmp &delta2, Cmp &nbar1, Cmp &nbar2, Cmp &ener, Cmp &press, int nzet, int l_min=0) const
General computational method for Cmp &#39;s, it computes both baryon densities, energy and pressure profi...
Definition: eos_bifluid.C:286
double * p_ray_eq2_pis2
Coordinate radius at , .
double * p_ray_pole2
Coordinate radius at .
Tenseur ent
Log-enthalpy (relativistic case) or specific enthalpy (Newtonian case)
Definition: etoile.h:460
virtual double r_circ2() const
Circumferential radius for fluid 2.
Tenseur nbar2
Baryon density in the fluid frame, for fluid 2.
Tenseur K_nn
Coefficient .
Class for a two-fluid (tabulated) equation of state.
Definition: eos_bifluid.h:1436
Tenseur & logn
Metric potential = logn_auto.
Definition: etoile.h:1524
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: etoile_rot.C:374
virtual double coupling_mominert_1() const
Quantities used to study the different fluid couplings: , and .
Basic array class.
Definition: tbl.h:164
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
Definition: grilles.h:474
virtual double angu_mom() const
Angular momentum.
double mass_b2() const
Baryon mass of fluid 2.
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: tenseur.C:652
virtual ~Et_rot_bifluid()
Destructor.
Cmp get_Knp(const Cmp &nbar1, const Cmp &nbar2, const Cmp &x2, int nzet, int l_min=0) const
Computes the derivatives of the energy with respect to .
Definition: eos_bifluid.C:761
virtual void del_hydro_euler()
Sets to ETATNONDEF (undefined state) the hydrodynamical quantities relative to the Eulerian observer...
Definition: etoile_rot.C:400
Tenseur alpha_eos
Coefficient relative to entrainment effects.
Tbl ** t
Array (size nzone ) of pointers on the Tbl &#39;s.
Definition: mtbl.h:132
Tenseur delta_car
The "relative velocity" (squared) of the two fluids.
Tenseur ener_euler
Total energy density in the Eulerian frame.
Definition: etoile.h:468
double ray_eq2() const
Coordinate radius for fluid 2 at , [r_unit].
Tenseur & dzeta
Metric potential = beta_auto.
Definition: etoile.h:1537
virtual double aplat2() const
Flatening r_pole/r_eq for fluid 2.
double * p_aplat2
Flatening r_pole/r_eq of fluid no.2.
virtual void hydro_euler()
Computes the hydrodynamical quantities relative to the Eulerian observer from those in the fluid fram...
double get_m1() const
Return the individual particule mass .
Definition: eos_bifluid.h:263
Tensor handling *** DEPRECATED : use class Tensor instead ***.
Definition: tenseur.h:304
void set_etat_nondef()
Sets the logical state to ETATNONDEF (undefined state).
Definition: tenseur.C:666
Et_rot_bifluid(Map &mp_i, int nzet_i, bool relat, const Eos_bifluid &eos_i)
Standard constructor.
Itbl * p_l_surf2
Description of the surface of fluid 2: 2-D Itbl containing the values of the domain index l on the su...
virtual double aplat() const
Flatening r_pole/r_eq.