LORENE
star_bin_xcts.C
1  /*
2  * Methods for the class Star_bin_xcts
3  * (see file star.h for documentation)
4  */
5 
6 /*
7  * Copyright (c) 2010 Michal Bejger
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 char star_bin_xcts_C[] = "$Header: /cvsroot/Lorene/C++/Source/Star/star_bin_xcts.C,v 1.9 2014/10/13 08:53:39 j_novak Exp $" ;
27 
28 /*
29  * $Id: star_bin_xcts.C,v 1.9 2014/10/13 08:53:39 j_novak Exp $
30  * $Log: star_bin_xcts.C,v $
31  * Revision 1.9 2014/10/13 08:53:39 j_novak
32  * Lorene classes and functions now belong to the namespace Lorene.
33  *
34  * Revision 1.8 2010/12/09 10:45:26 m_bejger
35  * Added psi4 in the constructor, removed fait_decouple
36  *
37  * Revision 1.7 2010/10/28 13:50:27 m_bejger
38  * Added mass-shedding estimation to Star_bin_xcts::operator>>
39  *
40  * Revision 1.6 2010/10/26 20:18:34 m_bejger
41  * Correction to stdin output
42  *
43  * Revision 1.5 2010/10/26 19:57:02 m_bejger
44  * Various cleanups
45  *
46  * Revision 1.4 2010/06/17 15:06:25 m_bejger
47  * N, Psi output corrected
48  *
49  * Revision 1.3 2010/06/15 08:18:19 m_bejger
50  * Method set_chi_comp() added
51  *
52  * Revision 1.2 2010/06/04 19:59:56 m_bejger
53  * Corrected definitions of lapse and Psi4
54  *
55  * Revision 1.1 2010/05/04 07:51:05 m_bejger
56  * Initial version
57  *
58  * $Header: /cvsroot/Lorene/C++/Source/Star/star_bin_xcts.C,v 1.9 2014/10/13 08:53:39 j_novak Exp $
59  *
60  */
61 
62 // Headers C
63 #include "math.h"
64 
65 // Headers Lorene
66 #include "etoile.h"
67 #include "star.h"
68 #include "eos.h"
69 #include "unites.h"
70 
71 // Local prototype
72 namespace Lorene {
73 Cmp raccord_c1(const Cmp& uu, int l1) ;
74 
75  //--------------//
76  // Constructors //
77  //--------------//
78 
79 // Standard constructor
80 // --------------------
82  int nzet_i,
83  const Eos& eos_i,
84  bool irrot)
85  : Star(mpi, nzet_i, eos_i),
86  irrotational(irrot),
87  psi0(mpi),
88  d_psi(mpi, COV, mpi.get_bvect_cart()),
89  wit_w(mpi, CON, mpi.get_bvect_cart()),
90  loggam(mpi),
91  bsn(mpi, CON, mpi.get_bvect_cart()),
92  pot_centri(mpi),
93  chi_auto(mpi),
94  chi_comp(mpi),
95  Psi_auto(mpi),
96  Psi_comp(mpi),
97  Psi(mpi),
98  chi(mpi),
99  psi4(mpi),
100  w_beta(mpi, CON, mpi.get_bvect_cart()),
101  khi(mpi),
102  dcov_Psi(mpi, COV, mpi.get_bvect_cart()),
103  dcov_chi(mpi, COV, mpi.get_bvect_cart()),
104  flat(mpi, mpi.get_bvect_cart()),
105  beta_auto(mpi, CON, mpi.get_bvect_cart()),
106  beta_comp(mpi, CON, mpi.get_bvect_cart()),
107  haij_auto(mpi, CON, mpi.get_bvect_cart()),
108  haij_comp(mpi, CON, mpi.get_bvect_cart()),
109  hacar_auto(mpi),
110  hacar_comp(mpi),
111  ssjm1_chi(mpi),
112  ssjm1_psi(mpi),
113  ssjm1_khi(mpi),
114  ssjm1_wbeta(mpi, CON, mpi.get_bvect_cart()) {
115 
116  // Pointers of derived quantities initialized to zero :
117  set_der_0x0() ;
118 
119  // Quantities defined on a spherical triad in star.C are put on
120  // a cartesian one
124  Metric temp_met (mp.flat_met_cart()) ;
125  gamma = temp_met ;
126 
127  // Initialization of all quantities :
128  psi0 = 0 ;
129  d_psi.set_etat_zero() ;
130  wit_w.set_etat_zero() ;
131  loggam = 0 ;
132  bsn.set_etat_zero() ;
133  pot_centri = 0 ;
134 
135  Psi_auto = 0 ;
136  Psi_comp = 0 ;
137 
140  chi_auto = 0 ;
141  chi_comp = 0 ;
142 
143  Psi = 1. ;
144  chi = 1. ;
145 
147  khi = 0 ;
148 
151 
154  hacar_auto = 0 ;
155  hacar_comp = 0 ;
156  ssjm1_psi = 0 ;
157  ssjm1_chi = 0 ;
158  ssjm1_khi = 0 ;
160 
161 }
162 
163 // Copy constructor
164 // ----------------
166  : Star(star),
167  irrotational(star.irrotational),
168  psi0(star.psi0),
169  d_psi(star.d_psi),
170  wit_w(star.wit_w),
171  loggam(star.loggam),
172  bsn(star.bsn),
173  pot_centri(star.pot_centri),
174  chi_auto(star.chi_auto),
175  chi_comp(star.chi_comp),
176  Psi_auto(star.Psi_auto),
177  Psi_comp(star.Psi_comp),
178  Psi(star.Psi),
179  chi(star.chi),
180  psi4(star.psi4),
181  w_beta(star.w_beta),
182  khi(star.khi),
183  dcov_Psi(star.dcov_Psi),
184  dcov_chi(star.dcov_chi),
185  flat(star.flat),
186  beta_auto(star.beta_auto),
187  beta_comp(star.beta_comp),
188  haij_auto(star.haij_auto),
189  haij_comp(star.haij_comp),
190  hacar_auto(star.hacar_auto),
191  hacar_comp(star.hacar_comp),
192  ssjm1_chi(star.ssjm1_chi),
193  ssjm1_psi(star.ssjm1_psi),
194  ssjm1_khi(star.ssjm1_khi),
195  ssjm1_wbeta(star.ssjm1_wbeta) {
196 
197  set_der_0x0() ;
198 
199 }
200 
201 // Constructor from a file
202 // -----------------------
204  const Eos& eos_i,
205  FILE* fich)
206  : Star(mpi, eos_i, fich),
207  psi0(mpi),
208  d_psi(mpi, COV, mpi.get_bvect_cart()),
209  wit_w(mpi, CON, mpi.get_bvect_cart()),
210  loggam(mpi),
211  bsn(mpi, CON, mpi.get_bvect_cart()),
212  pot_centri(mpi),
213  chi_auto(mpi, *(mpi.get_mg()), fich),
214  chi_comp(mpi),
215  Psi_auto(mpi, *(mpi.get_mg()), fich),
216  Psi_comp(mpi),
217  Psi(mpi),
218  chi(mpi),
219  psi4(mpi),
220  w_beta(mpi, mpi.get_bvect_cart(), fich),
221  khi(mpi, *(mpi.get_mg()), fich),
222  dcov_Psi(mpi, COV, mpi.get_bvect_cart()),
223  dcov_chi(mpi, COV, mpi.get_bvect_cart()),
224  flat(mpi, mpi.get_bvect_cart()),
225  beta_auto(mpi, mpi.get_bvect_cart(), fich),
226  beta_comp(mpi, CON, mpi.get_bvect_cart()),
227  haij_auto(mpi, CON, mpi.get_bvect_cart()),
228  haij_comp(mpi, CON, mpi.get_bvect_cart()),
229  hacar_auto(mpi),
230  hacar_comp(mpi),
231  ssjm1_chi(mpi, *(mpi.get_mg()), fich),
232  ssjm1_psi(mpi, *(mpi.get_mg()), fich),
233  ssjm1_khi(mpi, *(mpi.get_mg()), fich),
234  ssjm1_wbeta(mpi, mpi.get_bvect_cart(), fich) {
235 
236  // Star parameters
237  // -----------------
238 
239  // irrotational is read in the file:
240  bool status = fread(&irrotational, sizeof(bool), 1, fich) ;
241  if(!status)
242  cout << "Star_bin_xcts::Constructor from a file: Problem with reading ! " << endl ;
243 
244  // Read of the saved fields:
245  // ------------------------
246 
247  if (irrotational) {
248  Scalar psi0_file(mp, *(mp.get_mg()), fich) ;
249  psi0 = psi0_file ;
250 
251  Scalar gam_euler_file(mp, *(mp.get_mg()), fich) ;
252  gam_euler = gam_euler_file ;
253  }
254 
255  // Quantities defined on a spherical triad in star.C are put on
256  // a cartesian one
257 
261  Metric temp_met (mp.flat_met_cart()) ;
262  gamma = temp_met ;
263 
264  // All other fields are initialized to initial values :
265  // ----------------------------------------------------
266 
267  d_psi.set_etat_zero() ;
268  wit_w.set_etat_zero() ;
269  loggam = 0 ;
270  bsn.set_etat_zero() ;
271  pot_centri = 0 ;
272  Psi_comp = 0 ;
275  chi_comp = 0 ;
277  khi = 0 ;
281  hacar_auto = 0 ;
282  hacar_comp = 0 ;
283 
284  // Pointers of derived quantities initialized to zero
285  // --------------------------------------------------
286  set_der_0x0() ;
287 
288 }
289 
290  //------------//
291  // Destructor //
292  //------------//
293 
295 
296  del_deriv() ;
297 
298 }
299 
300  //----------------------------------//
301  // Management of derived quantities //
302  //----------------------------------//
303 
305 
306  Star::del_deriv() ;
307 
308  if (p_xa_barycenter != 0x0) delete p_xa_barycenter ;
309 
310  set_der_0x0() ;
311 }
312 
314 
316 
317  p_xa_barycenter = 0x0 ;
318 
319 }
320 
322 
324 
325  del_deriv() ;
326 
327 }
328 
329 
330  //--------------//
331  // Assignment //
332  //--------------//
333 
334 // Assignment to another Star_bin_xcts
335 // --------------------------------
337 
338  // Assignment of quantities common to the derived classes of Star
339  Star::operator=(star) ;
340 
341  // Assignement of proper quantities of class Star_bin_xcts
342  irrotational = star.irrotational ;
343  psi0 = star.psi0 ;
344  d_psi = star.d_psi ;
345  wit_w = star.wit_w ;
346  loggam = star.loggam ;
347  bsn = star.bsn ;
348  pot_centri = star.pot_centri ;
349  Psi_auto = star.Psi_auto ;
350  Psi_comp = star.Psi_comp ;
351  chi_auto = star.chi_auto ;
352  chi_comp = star.chi_comp ;
353  Psi = star.Psi ;
354  chi = star.chi ;
355  psi4 = star.psi4 ;
356  w_beta = star.w_beta ;
357  khi = star.khi ;
358  dcov_Psi = star.dcov_Psi ;
359  dcov_chi = star.dcov_chi ;
360  flat = star.flat ;
361  beta_auto = star.beta_auto ;
362  beta_comp = star.beta_comp ;
363  haij_auto = star.haij_auto ;
364  haij_comp = star.haij_comp ;
365  hacar_auto = star.hacar_auto ;
366  hacar_comp = star.hacar_comp ;
367  ssjm1_psi = star.ssjm1_psi ;
368  ssjm1_chi = star.ssjm1_chi ;
369  ssjm1_khi = star.ssjm1_khi ;
370  ssjm1_wbeta = star.ssjm1_wbeta ;
371 
372  del_deriv() ; // Deletes all derived quantities
373 
374 }
375 
377 
378  del_deriv() ; // sets to 0x0 all the derived quantities
379  return pot_centri ;
380 
381 }
382 
384 
385  del_deriv() ; // sets to 0x0 all the derived quantities
386  return Psi_comp ;
387 
388 }
389 
391 
392  del_deriv() ; // sets to 0x0 all the derived quantities
393  return Psi_auto ;
394 
395 }
396 
398 
399  del_deriv() ; // sets to 0x0 all the derived quantities
400  return chi_comp ;
401 
402 }
403 
405 
406  del_deriv() ; // sets to 0x0 all the derived quantities
407  return chi_auto ;
408 
409 }
410 
412 
413  del_deriv() ; // sets to 0x0 all the derived quantities
414  return beta_auto ;
415 
416 }
417 
419 
420  del_deriv() ; // sets to 0x0 all the derived quantities
421  return beta ;
422 
423 }
424 
425 
426  //--------------//
427  // Outputs //
428  //--------------//
429 
430 // Save in a file
431 // --------------
432 void Star_bin_xcts::sauve(FILE* fich) const {
433 
434  Star::sauve(fich) ;
435 
436  chi_auto.sauve(fich) ;
437  Psi_auto.sauve(fich) ;
438 
439  w_beta.sauve(fich) ;
440  khi.sauve(fich) ;
441 
442  beta_auto.sauve(fich) ;
443 
444  ssjm1_chi.sauve(fich) ;
445  ssjm1_psi.sauve(fich) ;
446  ssjm1_khi.sauve(fich) ;
447  ssjm1_wbeta.sauve(fich) ;
448 
449  fwrite(&irrotational, sizeof(bool), 1, fich) ;
450 
451  if (irrotational) {
452 
453  psi0.sauve(fich) ;
454  gam_euler.sauve(fich) ; // required to construct d_psi from psi0
455 
456  }
457 
458 }
459 
460 // Printing
461 // --------
462 
463 ostream& Star_bin_xcts::operator>>(ostream& ost) const {
464 
465  using namespace Unites ;
466 
467 // Star::operator>>(ost) ;
468 
469  ost << endl ;
470 
471  ost << "Number of domains occupied by the star : " << nzet << endl ;
472 
473  ost << "Equation of state : " << endl ;
474  ost << eos << endl ;
475 
476  ost << endl << "Central enthalpy : " << ent.val_grid_point(0,0,0,0) << " c^2" << endl ;
477  ost << "Central proper baryon density : " << nbar.val_grid_point(0,0,0,0)
478  << " x 0.1 fm^-3" << endl ;
479  ost << "Central proper energy density : " << ener.val_grid_point(0,0,0,0)
480  << " rho_nuc c^2" << endl ;
481  ost << "Central pressure : " << press.val_grid_point(0,0,0,0)
482  << " rho_nuc c^2" << endl ;
483 
484  ost << endl ;
485 
486  Scalar psi4_local = pow(Psi_auto + Psi_comp + 1., 4.) ;
487  psi4_local.std_spectral_base() ;
488 
489  ost << "Central lapse N : " << nn.val_grid_point(0,0,0,0) << endl ;
490  ost << "Central value of Psi^4 : " << psi4_local.val_grid_point(0,0,0,0) << endl ;
491 
492  ost << endl
493  << "Coordinate equatorial radius (phi=0) a1 = "
494  << ray_eq()/km << " km" << endl ;
495  ost << "Coordinate equatorial radius (phi=pi/2) a2 = "
496  << ray_eq_pis2()/km << " km" << endl ;
497  ost << "Coordinate equatorial radius (phi=pi): "
498  << ray_eq_pi()/km << " km" << endl ;
499  ost << "Coordinate polar radius a3 = "
500  << ray_pole()/km << " km" << endl ;
501  ost << "Axis ratio a2/a1 = " << ray_eq_pis2() / ray_eq()
502  << " a3/a1 = " << ray_pole() / ray_eq() << endl ;
503 
504  double dent_eq = ent.dsdr().val_point(ray_eq(),M_PI/2.,0.) ;
505  double dent_pole = ent.dsdr().val_point(ray_pole(),0.,0.) ;
506  double mass_shedd_chi = fabs( dent_eq / dent_pole ) ;
507 
508  ost << "Mass-shedding estimator = " << mass_shedd_chi << endl ;
509 
510  ost << endl << "Baryon mass : " << mass_b() / msol << " M_sol" << endl ;
511  ost << "Gravitational mass : " << mass_g() / msol << " M_sol" << endl ;
512 
513 
514 
515  ost << endl ;
516  ost << "Star in a binary system" << endl ;
517  ost << "-----------------------" << endl ;
518 
519  if (irrotational) {
520  ost << "irrotational configuration" << endl ;
521  }
522  else {
523  ost << "corotating configuration" << endl ;
524  }
525 
526  ost << "Absolute abscidia of the stellar center: " <<
527  mp.get_ori_x() / km << " km" << endl ;
528 
529  ost << "Absolute abscidia of the barycenter of the baryon density : " <<
530  xa_barycenter() / km << " km" << endl ;
531 
532  double r_0 = 0.5 * ( ray_eq() + ray_eq_pi() ) ;
533  double d_ns = fabs( mp.get_ori_x() ) + ray_eq_pi() - r_0 ;
534  double d_tilde = 2 * d_ns / r_0 ;
535 
536  ost << "d_tilde : " << d_tilde << endl ;
537 
538  ost << "Central value of gam_euler : "
539  << gam_euler.val_grid_point(0, 0, 0, 0) << endl ;
540 
541  ost << "Central u_euler (U^r, U^t, U^p) [c] : "
542  << u_euler(1).val_grid_point(0, 0, 0, 0) << " "
543  << u_euler(2).val_grid_point(0, 0, 0, 0) << " "
544  << u_euler(3).val_grid_point(0, 0, 0, 0) << endl ;
545 
546  if (irrotational) {
547  ost << "Central d_psi (r, t, p) [c] : "
548  << d_psi(1).val_grid_point(0, 0, 0, 0) << " "
549  << d_psi(2).val_grid_point(0, 0, 0, 0) << " "
550  << d_psi(3).val_grid_point(0, 0, 0, 0) << endl ;
551 
552  ost << "Central vel. / co-orb. (W^r, W^t, W^p) [c] : "
553  << wit_w(1).val_grid_point(0, 0, 0, 0) << " "
554  << wit_w(2).val_grid_point(0, 0, 0, 0) << " "
555  << wit_w(3).val_grid_point(0, 0, 0, 0) << endl ;
556 
557  ost << "Max vel. / co-orb. (W^r, W^t, W^p) [c] : "
558  << max(max(wit_w(1))) << " "
559  << max(max(wit_w(2))) << " "
560  << max(max(wit_w(3))) << endl ;
561 
562  ost << "Min vel. / co-orb. (W^r, W^t, W^p) [c] : "
563  << min(min(wit_w(1))) << " "
564  << min(min(wit_w(2))) << " "
565  << min(min(wit_w(3))) << endl ;
566 
567  double r_surf = mp.val_r(0,1.,M_PI/4,M_PI/4) ;
568 
569  ost << "Velocity at (r_surf,pi/4,pi/4) / co-orb. [c] : "
570  << wit_w(1).val_point(r_surf,M_PI/4,M_PI/4) << " "
571  << wit_w(2).val_point(r_surf,M_PI/4,M_PI/4) << " "
572  << wit_w(3).val_point(r_surf,M_PI/4,M_PI/4) << endl ;
573 
574  ost << "Central value of loggam : "
575  << loggam.val_grid_point(0, 0, 0, 0) << endl ;
576  }
577 
578  ost << "Central value of Psi auto, comp : "
579  << Psi_auto.val_grid_point(0, 0, 0, 0) << " "
580  << Psi_comp.val_grid_point(0, 0, 0, 0) << endl ;
581 
582  ost << "Central value of beta (N^r, N^t, N^p) [c] : "
583  << beta(1).val_grid_point(0, 0, 0, 0) << " "
584  << beta(2).val_grid_point(0, 0, 0, 0) << " "
585  << beta(3).val_grid_point(0, 0, 0, 0) << endl ;
586 
587  ost << " ... beta_auto part of it [c] : "
588  << beta_auto(1).val_grid_point(0, 0, 0, 0) << " "
589  << beta_auto(2).val_grid_point(0, 0, 0, 0) << " "
590  << beta_auto(3).val_grid_point(0, 0, 0, 0) << endl ;
591 
592  ost << endl << "Central value of (B^r, B^t, B^p)/N [c] : "
593  << bsn(1).val_grid_point(0, 0, 0, 0) << " "
594  << bsn(2).val_grid_point(0, 0, 0, 0) << " "
595  << bsn(3).val_grid_point(0, 0, 0, 0) << endl ;
596 
597 
598  ost << endl << "Central \\hat{A}^{ij} [c/km] : " << endl ;
599  ost << " \\hat{A}^{xx} auto, comp : "
600  << haij_auto(1, 1).val_grid_point(0, 0, 0, 0) * km << " "
601  << haij_comp(1, 1).val_grid_point(0, 0, 0, 0) * km << endl ;
602  ost << " A^{xy} auto, comp : "
603  << haij_auto(1, 2).val_grid_point(0, 0, 0, 0) * km << " "
604  << haij_comp(1, 2).val_grid_point(0, 0, 0, 0) * km << endl ;
605  ost << " A^{xz} auto, comp : "
606  << haij_auto(1, 3).val_grid_point(0, 0, 0, 0) * km << " "
607  << haij_comp(1, 3).val_grid_point(0, 0, 0, 0) * km << endl ;
608  ost << " A^{yy} auto, comp : "
609  << haij_auto(2, 2).val_grid_point(0, 0, 0, 0) * km << " "
610  << haij_comp(2, 2).val_grid_point(0, 0, 0, 0) * km << endl ;
611  ost << " A^{yz} auto, comp : "
612  << haij_auto(2, 3).val_grid_point(0, 0, 0, 0) * km << " "
613  << haij_comp(2, 3).val_grid_point(0, 0, 0, 0) * km << endl ;
614  ost << " A^{zz} auto, comp : "
615  << haij_auto(3, 3).val_grid_point(0, 0, 0, 0) * km << " "
616  << haij_comp(3, 3).val_grid_point(0, 0, 0, 0) * km << endl ;
617 
618  ost << endl << "Central \\hat{A}_{ij}\\hat{A}^{ij} [c^2/km^2] : "
619  << endl ;
620  ost << " \\hat{A}_{ij}\\hat{A}^{ij} auto, comp : "
621  << hacar_auto.val_grid_point(0, 0, 0, 0) * km*km << " "
622  << hacar_comp.val_grid_point(0, 0, 0, 0) * km*km << endl ;
623 
624  return ost ;
625 }
626 
627  //-------------------------//
628  // Computational routines //
629  //-------------------------//
630 
632 
633  if (!irrotational) {
635  return ;
636  }
637 
638  // Specific relativistic enthalpy ---> hhh
639  //----------------------------------
640 
641  Scalar hhh = exp(ent) ; // = 1 at the Newtonian limit
642 
643  // Computation of W^i = - h Gamma_n B^i/N
644  //----------------------------------------------
645 
646  Vector www = hhh * gam_euler * bsn * psi4 ;
647 
648  // Constant value of W^i at the center of the star
649  //-------------------------------------------------
650 
651  Vector v_orb(mp, COV, mp.get_bvect_cart()) ;
652 
653  for (int i=1; i<=3; i++)
654  v_orb.set(i) = www(i).val_grid_point(0, 0, 0, 0) ;
655 
656  // Gradient of psi
657  //----------------
658 
659  Vector d_psi0 = psi0.derive_cov(flat) ;
660 
661  d_psi0.change_triad( mp.get_bvect_cart() ) ;
662  d_psi0.std_spectral_base() ;
663 
664  d_psi = d_psi0 + v_orb ;
665 
666  for (int i=1; i<=3; i++) {
667  if (d_psi(i).get_etat() == ETATZERO)
668  d_psi.set(i).annule_hard() ;
669  }
670 
671  // C^1 continuation of d_psi outside the star
672  // (to ensure a smooth enthalpy field accross the stellar surface)
673  // ----------------------------------------------------------------
674 
675  int nzm1 = mp.get_mg()->get_nzone() - 1 ;
676  d_psi.annule(nzet, nzm1) ;
677 
678  for (int i=1; i<=3; i++) {
679  Cmp d_psi_i (d_psi(i)) ;
680  d_psi_i.va.set_base( d_psi0(i).get_spectral_va().base ) ;
681  d_psi_i = raccord_c1(d_psi_i, nzet) ;
682  d_psi.set(i) = d_psi_i ;
683  }
684 
686 
687 }
688 
690  double relax_ent,
691  double relax_met,
692  int mer,
693  int fmer_met) {
694 
695  double relax_ent_jm1 = 1. - relax_ent ;
696  double relax_met_jm1 = 1. - relax_met ;
697 
698  ent = relax_ent * ent + relax_ent_jm1 * star_jm1.ent ;
699 
700  if ( (mer != 0) && (mer % fmer_met == 0)) {
701 
702  Psi_auto = relax_met * Psi_auto
703  + relax_met_jm1 * star_jm1.Psi_auto ;
704 
705  chi_auto = relax_met * chi_auto
706  + relax_met_jm1 * star_jm1.chi_auto ;
707 
708  beta_auto = relax_met * beta_auto
709  + relax_met_jm1 * star_jm1.beta_auto ;
710 
711  }
712 
713  del_deriv() ;
714 
716 
717 }
718 }
Lorene::Star::stress_euler
Sym_tensor stress_euler
Spatial part of the stress-energy tensor with respect to the Eulerian observer.
Definition: star.h:212
Lorene::Star::nbar
Scalar nbar
Baryon density in the fluid frame.
Definition: star.h:192
Lorene::Star_bin_xcts::pot_centri
Scalar pot_centri
Centrifugal potential.
Definition: star.h:1129
Lorene::Star_bin_xcts::dcov_chi
Vector dcov_chi
Covariant derivative of the function .
Definition: star.h:1174
Lorene::Star_bin_xcts::loggam
Scalar loggam
Logarithm of the Lorentz factor between the fluid and the co-orbiting observer.
Definition: star.h:1120
Lorene::Scalar::sauve
virtual void sauve(FILE *) const
Save in a file.
Definition: scalar.C:686
Lorene::Star_bin_xcts::xa_barycenter
virtual double xa_barycenter() const
Absolute coordinate X of the barycenter of the baryon density,.
Definition: star_bin_global_xcts.C:99
Lorene::Map::get_ori_x
double get_ori_x() const
Returns the x coordinate of the origin.
Definition: map.h:768
Lorene::Star_bin_xcts::set_chi_auto
Scalar & set_chi_auto()
Read/write the conformal factor .
Definition: star_bin_xcts.C:404
Lorene::Scalar::val_point
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
Lorene::Star_bin_xcts::operator=
void operator=(const Star_bin_xcts &)
Assignment to another Star_bin_xcts.
Definition: star_bin_xcts.C:336
Lorene::Star::ray_eq_pis2
double ray_eq_pis2() const
Coordinate radius at , [r_unit].
Definition: star_global.C:138
Lorene::Star::ray_eq
double ray_eq() const
Coordinate radius at , [r_unit].
Definition: star_global.C:108
Lorene::Metric
Metric for tensor calculation.
Definition: metric.h:90
Lorene::Star::del_deriv
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: star.C:298
Lorene::Star_bin_xcts::set_der_0x0
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: star_bin_xcts.C:313
Lorene::Tensor::set_etat_nondef
virtual void set_etat_nondef()
Sets the logical state of all components to ETATNONDEF (undefined state).
Definition: tensor.C:489
Lorene::Star::ener
Scalar ener
Total energy density in the fluid frame.
Definition: star.h:193
Lorene::Star::gam_euler
Scalar gam_euler
Lorentz factor between the fluid and Eulerian observers.
Definition: star.h:204
Lorene::Scalar::dsdr
const Scalar & dsdr() const
Returns of *this .
Definition: scalar_deriv.C:113
Lorene::Star_bin_xcts::beta_auto
Vector beta_auto
Part of the shift vector generated principally by the star (Spherical components with respect to the ...
Definition: star.h:1182
Lorene::Star_bin_xcts::set_beta
Vector & set_beta()
Read/write of .
Definition: star_bin_xcts.C:418
Lorene::exp
Cmp exp(const Cmp &)
Exponential.
Definition: cmp_math.C:270
Lorene::Star_bin_xcts::chi_comp
Scalar chi_comp
Scalar field generated principally by the companion star.
Definition: star.h:1139
Lorene
Lorene prototypes.
Definition: app_hor.h:64
Lorene::Star_bin_xcts::chi
Scalar chi
Total function .
Definition: star.h:1155
Lorene::Star::press
Scalar press
Fluid pressure.
Definition: star.h:194
Lorene::Star::nzet
int nzet
Number of domains of *mp occupied by the star.
Definition: star.h:183
Lorene::Tensor::set_etat_zero
virtual void set_etat_zero()
Sets the logical state of all components to ETATZERO (zero state).
Definition: tensor.C:497
Lorene::Star_bin_xcts::irrotational
bool irrotational
true for an irrotational star, false for a corotating one
Definition: star.h:1099
Lorene::Cmp::va
Valeur va
The numerical value of the Cmp
Definition: cmp.h:464
Lorene::Map::get_mg
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:765
Lorene::Tensor::annule
virtual void annule(int l_min, int l_max)
Sets the Tensor to zero in several domains.
Definition: tensor.C:671
Lorene::Scalar::derive_cov
const Vector & derive_cov(const Metric &gam) const
Returns the gradient (1-form = covariant vector) of *this
Definition: scalar_deriv.C:390
Lorene::Eos
Equation of state base class.
Definition: eos.h:190
Lorene::Map::val_r
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.
Lorene::max
Tbl max(const Cmp &)
Maximum values of a Cmp in each domain.
Definition: cmp_math.C:435
Lorene::Star_bin_xcts::hacar_comp
Scalar hacar_comp
Part of the scalar generated by beta_auto and beta_comp, i.e.
Definition: star.h:1211
Lorene::Star_bin_xcts::sauve
virtual void sauve(FILE *) const
Save in a file.
Definition: star_bin_xcts.C:432
Lorene::Star_bin_xcts::haij_comp
Sym_tensor haij_comp
Part of the extrinsic curvature tensor generated by beta_comp.
Definition: star.h:1199
Lorene::Scalar
Tensor field of valence 0 (or component of a tensorial field).
Definition: scalar.h:387
Lorene::Star::nn
Scalar nn
Lapse function N .
Definition: star.h:225
Lorene::pow
Cmp pow(const Cmp &, int)
Power .
Definition: cmp_math.C:348
Lorene::Star_bin_xcts::mass_g
virtual double mass_g() const
Gravitational mass.
Definition: star_bin_global_xcts.C:74
Lorene::Star_bin_xcts::wit_w
Vector wit_w
Spatial projection of the fluid 3-velocity with respect to the co-orbiting observer.
Definition: star.h:1115
Lorene::Star_bin_xcts::psi0
Scalar psi0
Scalar potential of the non-translational part of fluid 4-velocity (in the irrotational case)
Definition: star.h:1104
Lorene::Vector::change_triad
virtual void change_triad(const Base_vect &)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
Definition: vector_change_triad.C:75
Lorene::Star_bin_xcts::psi4
Scalar psi4
Conformal factor .
Definition: star.h:1158
Lorene::Star::u_euler
Vector u_euler
Fluid 3-velocity with respect to the Eulerian observer.
Definition: star.h:207
Unites
Standard units of space, time and mass.
Lorene::Star::ray_eq_pi
double ray_eq_pi() const
Coordinate radius at , [r_unit].
Definition: star_global.C:186
Lorene::Star_bin_xcts::w_beta
Vector w_beta
Solution for the vector part of the vector Poisson equation for .
Definition: star.h:1163
Lorene::Star_bin_xcts::Psi_comp
Scalar Psi_comp
Scalar field generated principally by the companion star.
Definition: star.h:1149
Lorene::Star::del_hydro_euler
virtual void del_hydro_euler()
Sets to ETATNONDEF (undefined state) the hydrodynamical quantities relative to the Eulerian observer.
Definition: star.C:330
Lorene::Star_bin_xcts::set_Psi_comp
Scalar & set_Psi_comp()
Read/write the conformal factor .
Definition: star_bin_xcts.C:383
Lorene::Star_bin_xcts::set_Psi_auto
Scalar & set_Psi_auto()
Read/write the conformal factor .
Definition: star_bin_xcts.C:390
Lorene::Star::equation_of_state
void equation_of_state()
Computes the proper baryon and energy density, as well as pressure from the enthalpy.
Definition: star.C:462
Lorene::Star::mp
Map & mp
Mapping associated with the star.
Definition: star.h:180
Lorene::Vector::set
Scalar & set(int)
Read/write access to a component.
Definition: vector.C:296
Lorene::Star_bin_xcts::chi_auto
Scalar chi_auto
Scalar field generated principally by the star.
Definition: star.h:1134
Lorene::Star_bin_xcts::del_hydro_euler
virtual void del_hydro_euler()
Sets to ETATNONDEF (undefined state) the hydrodynamical quantities relative to the Eulerian observer.
Definition: star_bin_xcts.C:321
Lorene::Star::gamma
Metric gamma
3-metric
Definition: star.h:235
Lorene::Star_bin_xcts::mass_b
virtual double mass_b() const
Baryon mass.
Definition: star_bin_global_xcts.C:50
Lorene::Star_bin_xcts::hacar_auto
Scalar hacar_auto
Part of the scalar generated by beta_auto, i.e.
Definition: star.h:1205
Lorene::Star_bin_xcts::set_chi_comp
Scalar & set_chi_comp()
Read/write the conformal factor .
Definition: star_bin_xcts.C:397
Lorene::Tensor::change_triad
virtual void change_triad(const Base_vect &new_triad)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
Definition: tensor_change_triad.C:77
Lorene::Star::sauve
virtual void sauve(FILE *) const
Save in a file.
Definition: star.C:397
Lorene::Star::ent
Scalar ent
Log-enthalpy.
Definition: star.h:190
Lorene::Cmp
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition: cmp.h:446
Lorene::min
Tbl min(const Cmp &)
Minimum values of a Cmp in each domain.
Definition: cmp_math.C:458
Lorene::Map::flat_met_cart
const Metric_flat & flat_met_cart() const
Returns the flat metric associated with the Cartesian coordinates and with components expressed in th...
Definition: map.C:331
Lorene::Star_bin_xcts::set_pot_centri
Scalar & set_pot_centri()
Read/write the centrifugal potential.
Definition: star_bin_xcts.C:376
Lorene::Star_bin_xcts::fait_d_psi
void fait_d_psi()
Computes the gradient of the total velocity potential .
Definition: star_bin_xcts.C:631
Lorene::Mg3d::get_nzone
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:448
Lorene::Star_bin_xcts::bsn
Vector bsn
3-vector shift, divided by N, of the rotating coordinates, .
Definition: star.h:1126
Lorene::Star
Base class for stars.
Definition: star.h:175
Lorene::Valeur::set_base
void set_base(const Base_val &)
Sets the bases for spectral expansions (member base )
Definition: valeur.C:810
Lorene::Map::get_bvect_cart
const Base_vect_cart & get_bvect_cart() const
Returns the Cartesian basis associated with the coordinates (x,y,z) of the mapping,...
Definition: map.h:791
Lorene::Star::set_der_0x0
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: star.C:316
Lorene::Vector
Tensor field of valence 1.
Definition: vector.h:188
Lorene::Star_bin_xcts::flat
Metric_flat flat
Flat metric defined on the mapping (Spherical components with respect to the mapping of the star) .
Definition: star.h:1177
Lorene::Star_bin_xcts::~Star_bin_xcts
virtual ~Star_bin_xcts()
Destructor.
Definition: star_bin_xcts.C:294
Lorene::Star_bin_xcts::ssjm1_wbeta
Vector ssjm1_wbeta
Effective source at the previous step for wbeta of the vector Poisson equation for wbeta (needed for ...
Definition: star.h:1234
Lorene::Star_bin_xcts::del_deriv
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: star_bin_xcts.C:304
Lorene::Star_bin_xcts::khi
Scalar khi
Solution for the scalar part of the vector Poisson equation for .
Definition: star.h:1168
Lorene::Star::operator=
void operator=(const Star &)
Assignment to another Star.
Definition: star.C:351
Lorene::Star_bin_xcts::set_beta_auto
Vector & set_beta_auto()
Read/write of .
Definition: star_bin_xcts.C:411
Lorene::Star_bin_xcts::Psi_auto
Scalar Psi_auto
Scalar field generated principally by the star.
Definition: star.h:1144
Lorene::Star_bin_xcts::Psi
Scalar Psi
Total conformal factor .
Definition: star.h:1152
Lorene::Scalar::std_spectral_base
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
Lorene::Star_bin_xcts
Class for stars in binary system in eXtended Conformal Thin Sandwich formulation.
Definition: star.h:1091
Lorene::Scalar::annule_hard
void annule_hard()
Sets the Scalar to zero in a hard way.
Definition: scalar.C:380
Lorene::Star_bin_xcts::Star_bin_xcts
Star_bin_xcts(Map &mp_i, int nzet_i, const Eos &eos_i, bool irrot)
Standard constructor.
Definition: star_bin_xcts.C:81
Lorene::Star_bin_xcts::operator>>
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition: star_bin_xcts.C:463
Lorene::Star_bin_xcts::ssjm1_khi
Scalar ssjm1_khi
Effective source at the previous step for the resolution of the Poisson equation for khi.
Definition: star.h:1227
Lorene::Star_bin_xcts::dcov_Psi
Vector dcov_Psi
Covariant derivative of the conformal factor .
Definition: star.h:1171
Lorene::Star::eos
const Eos & eos
Equation of state of the stellar matter.
Definition: star.h:185
Lorene::Star_bin_xcts::ssjm1_psi
Scalar ssjm1_psi
Effective source at the previous step for the resolution of the Poisson equation for \Psi_auto .
Definition: star.h:1221
Lorene::Star_bin_xcts::relaxation
void relaxation(const Star_bin_xcts &star_prev, double relax_ent, double relax_met, int mer, int fmer_met)
Performs a relaxation on ent, Psi_auto, chi_auto and beta_auto.
Definition: star_bin_xcts.C:689
Lorene::Scalar::val_grid_point
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
Lorene::Star_bin_xcts::haij_auto
Sym_tensor haij_auto
Part of the extrinsic curvature tensor generated by beta_auto.
Definition: star.h:1193
Lorene::Star_bin_xcts::d_psi
Vector d_psi
Gradient of (in the irrotational case) (Spherical components with respect to the mapping of the star...
Definition: star.h:1109
Lorene::Map
Base class for coordinate mappings.
Definition: map.h:670
Lorene::Star_bin_xcts::ssjm1_chi
Scalar ssjm1_chi
Effective source at the previous step for the resolution of the Poisson equation for \chi_auto .
Definition: star.h:1216
Lorene::Star::ray_pole
double ray_pole() const
Coordinate radius at [r_unit].
Definition: star_global.C:278
Lorene::Tensor::sauve
virtual void sauve(FILE *) const
Save in a binary file.
Definition: tensor.C:906
Lorene::Star::beta
Vector beta
Shift vector.
Definition: star.h:228
Lorene::Star_bin_xcts::p_xa_barycenter
double * p_xa_barycenter
Absolute coordinate X of the barycenter of the baryon density.
Definition: star.h:1240
Lorene::Star_bin_xcts::beta_comp
Vector beta_comp
Part of the shift vector generated principally by the star (Spherical components with respect to the ...
Definition: star.h:1187
Lorene::Vector::std_spectral_base
virtual void std_spectral_base()
Sets the standard spectal bases of decomposition for each component.
Definition: vector.C:316