LORENE
star_bin.C
1  /*
2  * Methods for the class Star_bin
3  *
4  * (see file star.h for documentation)
5  */
6 
7 /*
8  * Copyright (c) 2004 Francois Limousin
9  *
10  * This file is part of LORENE.
11  *
12  * LORENE is free software; you can redistribute it and/or modify
13  * it under the terms of the GNU General Public License as published by
14  * the Free Software Foundation; either version 2 of the License, or
15  * (at your option) any later version.
16  *
17  * LORENE is distributed in the hope that it will be useful,
18  * but WITHOUT ANY WARRANTY; without even the implied warranty of
19  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
20  * GNU General Public License for more details.
21  *
22  * You should have received a copy of the GNU General Public License
23  * along with LORENE; if not, write to the Free Software
24  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
25  *
26  */
27 
28 
29 char star_bin_C[] = "$Header: /cvsroot/Lorene/C++/Source/Star/star_bin.C,v 1.19 2014/10/13 08:53:37 j_novak Exp $" ;
30 
31 /*
32  * $Id: star_bin.C,v 1.19 2014/10/13 08:53:37 j_novak Exp $
33  * $Log: star_bin.C,v $
34  * Revision 1.19 2014/10/13 08:53:37 j_novak
35  * Lorene classes and functions now belong to the namespace Lorene.
36  *
37  * Revision 1.18 2006/04/11 14:24:44 f_limousin
38  * New version of the code : improvement of the computation of some
39  * critical sources, estimation of the dirac gauge, helical symmetry...
40  *
41  * Revision 1.17 2005/09/13 19:38:31 f_limousin
42  * Reintroduction of the resolution of the equations in cartesian coordinates.
43  *
44  * Revision 1.16 2005/04/08 12:36:44 f_limousin
45  * Just to avoid warnings...
46  *
47  * Revision 1.15 2005/02/24 16:03:01 f_limousin
48  * Change the name of some variables (for instance dcov_logn --> dlogn).
49  *
50  * Revision 1.14 2005/02/17 17:29:28 f_limousin
51  * Change the name of some quantities to be consistent with other classes
52  * (for instance nnn is changed to nn, shift to beta, beta to lnq...)
53  *
54  * Revision 1.13 2005/02/11 18:11:57 f_limousin
55  * Introduction of a new member Map_af.
56  *
57  * Revision 1.12 2004/11/11 16:29:49 j_novak
58  * set_der_0x0 is no longer virtual (to be coherent with Tensor/Scalar classes).
59  *
60  * Revision 1.11 2004/07/21 11:49:03 f_limousin
61  * Remove function sprod.
62  *
63  * Revision 1.10 2004/06/22 12:48:52 f_limousin
64  * Change qq, qq_auto and qq_comp to beta, beta_auto and beta_comp.
65  *
66  * Revision 1.9 2004/04/08 16:32:28 f_limousin
67  * The new variable is ln(Q) instead of Q=psi^2*N. It improves the
68  * convergence of the code.
69  *
70  * Revision 1.8 2004/03/25 10:29:26 j_novak
71  * All LORENE's units are now defined in the namespace Unites (in file unites.h).
72  *
73  * Revision 1.7 2004/03/23 09:54:54 f_limousin
74  * Add comments
75  *
76  * Revision 1.6 2004/02/27 09:50:57 f_limousin
77  * Scalars ssjm1_logn, ssjm1_qq ... have been added for all metric
78  * quantities for the resolution of Poisson equations.
79  * Members bsn and d_psi are now constructed on a cartesian triad.
80  *
81  * Revision 1.5 2004/02/21 17:05:13 e_gourgoulhon
82  * Method Scalar::point renamed Scalar::val_grid_point.
83  * Method Scalar::set_point renamed Scalar::set_grid_point.
84  *
85  * Revision 1.4 2004/01/22 10:07:18 f_limousin
86  * Add methods set_logn_comp() and set_shift_auto().
87  *
88  * Revision 1.3 2004/01/20 15:17:34 f_limousin
89  * First version
90  *
91  *
92  * $Header: /cvsroot/Lorene/C++/Source/Star/star_bin.C,v 1.19 2014/10/13 08:53:37 j_novak Exp $
93  *
94  */
95 
96 // Headers C
97 #include "math.h"
98 
99 // Headers Lorene
100 #include "etoile.h"
101 #include "star.h"
102 #include "eos.h"
103 #include "unites.h"
104 
105 // Local prototype
106 namespace Lorene {
107 Cmp raccord_c1(const Cmp& uu, int l1) ;
108 
109  //--------------//
110  // Constructors //
111  //--------------//
112 
113 // Standard constructor
114 // --------------------
115 Star_bin::Star_bin(Map& mpi, int nzet_i, const Eos& eos_i,
116  bool irrot, bool conf_flat0)
117  : Star(mpi, nzet_i, eos_i),
118  irrotational(irrot),
119  psi0(mpi),
120  d_psi(mpi, COV, mpi.get_bvect_cart()),
121  wit_w(mpi, CON, mpi.get_bvect_cart()),
122  loggam(mpi),
123  bsn(mpi, CON, mpi.get_bvect_cart()),
124  pot_centri(mpi),
125  logn_auto(mpi),
126  logn_comp(mpi),
127  dcov_logn(mpi, COV, mpi.get_bvect_cart()),
128  dcon_logn(mpi, CON, mpi.get_bvect_cart()),
129  lnq_auto(mpi),
130  lnq_comp(mpi),
131  psi4(mpi),
132  dcov_phi(mpi, COV, mpi.get_bvect_cart()),
133  dcon_phi(mpi, CON, mpi.get_bvect_cart()),
134  flat(mpi, mpi.get_bvect_cart()),
135  gtilde(flat),
136  beta_auto(mpi, CON, mpi.get_bvect_cart()),
137  beta_comp(mpi, CON, mpi.get_bvect_cart()),
138  hij(mpi, CON, mpi.get_bvect_cart()),
139  hij_auto(mpi, CON, mpi.get_bvect_cart()),
140  hij_comp(mpi, CON, mpi.get_bvect_cart()),
141  tkij_auto(mpi, CON, mpi.get_bvect_cart()),
142  tkij_comp(mpi, CON, mpi.get_bvect_cart()),
143  kcar_auto(mpi),
144  kcar_comp(mpi),
145  ssjm1_logn(mpi),
146  ssjm1_lnq(mpi),
147  ssjm1_khi(mpi),
148  ssjm1_wbeta(mpi, CON, mpi.get_bvect_cart()),
149  ssjm1_h11(mpi),
150  ssjm1_h21(mpi),
151  ssjm1_h31(mpi),
152  ssjm1_h22(mpi),
153  ssjm1_h32(mpi),
154  ssjm1_h33(mpi),
155  decouple(mpi),
156  conf_flat(conf_flat0){
157 
158  // Pointers of derived quantities initialized to zero :
159  set_der_0x0() ;
160 
161  // Quantities defined on a spherical triad in star.C are put on
162  // a cartesian one
166  Metric temp_met (mp.flat_met_cart()) ;
167  gamma = temp_met ;
168 
169  // All quantities are initialized to zero :
170  psi0 = 0 ;
171  d_psi.set_etat_zero() ;
172  wit_w.set_etat_zero() ;
173  loggam = 0 ;
174  bsn.set_etat_zero() ;
175  pot_centri = 0 ;
176 
177  logn_auto = 0 ;
178  logn_comp = 0 ;
183  lnq_auto = 0 ;
184  lnq_comp = 0 ;
185  psi4 = 1 ;
188  hij.set_etat_zero() ;
191 
194  kcar_auto = 0 ;
195  kcar_comp = 0 ;
196  ssjm1_logn = 0 ;
197  ssjm1_lnq = 0 ;
198  ssjm1_khi = 0 ;
199  ssjm1_wbeta.set_etat_zero() ;
200  ssjm1_h11 = 0 ;
201  ssjm1_h21 = 0 ;
202  ssjm1_h31 = 0 ;
203  ssjm1_h22 = 0 ;
204  ssjm1_h32 = 0 ;
205  ssjm1_h33 = 0 ;
206 }
207 
208 // Copy constructor
209 // ----------------
211  : Star(star),
212  irrotational(star.irrotational),
213  psi0(star.psi0),
214  d_psi(star.d_psi),
215  wit_w(star.wit_w),
216  loggam(star.loggam),
217  bsn(star.bsn),
218  pot_centri(star.pot_centri),
219  logn_auto(star.logn_auto),
220  logn_comp(star.logn_comp),
221  dcov_logn(star.dcov_logn),
222  dcon_logn(star.dcon_logn),
223  lnq_auto(star.lnq_auto),
224  lnq_comp(star.lnq_comp),
225  psi4(star.psi4),
226  dcov_phi(star.dcov_phi),
227  dcon_phi(star.dcon_phi),
228  flat(star.flat),
229  gtilde(star.gtilde),
230  beta_auto(star.beta_auto),
231  beta_comp(star.beta_comp),
232  hij(star.hij),
233  hij_auto(star.hij_auto),
234  hij_comp(star.hij_comp),
235  tkij_auto(star.tkij_auto),
236  tkij_comp(star.tkij_comp),
237  kcar_auto(star.kcar_auto),
238  kcar_comp(star.kcar_comp),
239  ssjm1_logn(star.ssjm1_logn),
240  ssjm1_lnq(star.ssjm1_lnq),
241  ssjm1_khi(star.ssjm1_khi),
242  ssjm1_wbeta(star.ssjm1_wbeta),
243  ssjm1_h11(star.ssjm1_h11),
244  ssjm1_h21(star.ssjm1_h21),
245  ssjm1_h31(star.ssjm1_h31),
246  ssjm1_h22(star.ssjm1_h22),
247  ssjm1_h32(star.ssjm1_h32),
248  ssjm1_h33(star.ssjm1_h33),
249  decouple(star.decouple),
250  conf_flat(star.conf_flat)
251 {
252  set_der_0x0() ;
253 
254 }
255 
256 // Constructor from a file
257 // -----------------------
258 Star_bin::Star_bin(Map& mpi, const Eos& eos_i, FILE* fich)
259  : Star(mpi, eos_i, fich),
260  psi0(mpi),
261  d_psi(mpi, COV, mpi.get_bvect_cart()),
262  wit_w(mpi, CON, mpi.get_bvect_cart()),
263  loggam(mpi),
264  bsn(mpi, CON, mpi.get_bvect_cart()),
265  pot_centri(mpi),
266  logn_auto(mpi, *(mpi.get_mg()), fich),
267  logn_comp(mpi),
268  dcov_logn(mpi, COV, mpi.get_bvect_cart()),
269  dcon_logn(mpi, CON, mpi.get_bvect_cart()),
270  lnq_auto(mpi, *(mpi.get_mg()), fich),
271  lnq_comp(mpi),
272  psi4(mpi),
273  dcov_phi(mpi, COV, mpi.get_bvect_cart()),
274  dcon_phi(mpi, CON, mpi.get_bvect_cart()),
275  flat(mpi, mpi.get_bvect_cart()),
276  gtilde(flat),
277  beta_auto(mpi, mpi.get_bvect_cart(), fich),
278  beta_comp(mpi, CON, mpi.get_bvect_cart()),
279  hij(mpi, CON, mpi.get_bvect_cart()),
280  hij_auto(mpi, mpi.get_bvect_cart(), fich),
281  hij_comp(mpi, CON, mpi.get_bvect_cart()),
282  tkij_auto(mpi, CON, mpi.get_bvect_cart()),
283  tkij_comp(mpi, CON, mpi.get_bvect_cart()),
284  kcar_auto(mpi),
285  kcar_comp(mpi),
286  ssjm1_logn(mpi, *(mpi.get_mg()), fich),
287  ssjm1_lnq(mpi, *(mpi.get_mg()), fich),
288  ssjm1_khi(mpi, *(mpi.get_mg()), fich),
289  ssjm1_wbeta(mpi, mpi.get_bvect_cart(), fich),
290  ssjm1_h11(mpi, *(mpi.get_mg()), fich),
291  ssjm1_h21(mpi, *(mpi.get_mg()), fich),
292  ssjm1_h31(mpi, *(mpi.get_mg()), fich),
293  ssjm1_h22(mpi, *(mpi.get_mg()), fich),
294  ssjm1_h32(mpi, *(mpi.get_mg()), fich),
295  ssjm1_h33(mpi, *(mpi.get_mg()), fich),
296  decouple(mpi){
297 
298  // Etoile parameters
299  // -----------------
300 
301  // irrotational and conf_flat is read in the file:
302  fread(&irrotational, sizeof(bool), 1, fich) ;
303  fread(&conf_flat, sizeof(bool), 1, fich) ;
304 
305 
306  // Read of the saved fields:
307  // ------------------------
308 
309  if (irrotational) {
310  Scalar gam_euler_file(mp, *(mp.get_mg()), fich) ;
311  gam_euler = gam_euler_file ;
312 
313  Scalar psi0_file(mp, *(mp.get_mg()), fich) ;
314  psi0 = psi0_file ;
315  }
316 
317  // Quantities defined on a spherical triad in star.C are put on
318  // a cartesian one
322  Metric temp_met (mp.flat_met_cart()) ;
323  gamma = temp_met ;
324 
325  // All other fields are initialized to zero :
326  // ----------------------------------------
327 
328  d_psi.set_etat_zero() ;
329  wit_w.set_etat_zero() ;
330  loggam = 0 ;
331  bsn.set_etat_zero() ;
332  pot_centri = 0 ;
333  logn_comp = 0 ;
337  lnq_comp = 0 ;
338  psi4 = 1 ;
341  hij.set_etat_zero() ;
343 
346  kcar_auto = 0 ;
347  kcar_comp = 0 ;
348 
349  // Pointers of derived quantities initialized to zero
350  // --------------------------------------------------
351  set_der_0x0() ;
352 
353 }
354 
355  //------------//
356  // Destructor //
357  //------------//
358 
360 
361  del_deriv() ;
362 
363 }
364 
365  //----------------------------------//
366  // Management of derived quantities //
367  //----------------------------------//
368 
369 void Star_bin::del_deriv() const {
370 
371  Star::del_deriv() ;
372 
373  if (p_xa_barycenter != 0x0) delete p_xa_barycenter ;
374 
375  set_der_0x0() ;
376 }
377 
378 
379 
380 
381 void Star_bin::set_der_0x0() const {
382 
384 
385  p_xa_barycenter = 0x0 ;
386 
387 }
388 
390 
392 
393  del_deriv() ;
394 
395 }
396 
397 
398  //--------------//
399  // Assignment //
400  //--------------//
401 
402 // Assignment to another Star_bin
403 // --------------------------------
404 void Star_bin::operator=(const Star_bin& star) {
405 
406  // Assignment of quantities common to the derived classes of Star
407  Star::operator=(star) ;
408 
409  // Assignement of proper quantities of class Star_bin
410  irrotational = star.irrotational ;
411  psi0 = star.psi0 ;
412  d_psi = star.d_psi ;
413  wit_w = star.wit_w ;
414  loggam = star.loggam ;
415  bsn = star.bsn ;
416  pot_centri = star.pot_centri ;
417  logn_auto = star.logn_auto ;
418  logn_comp = star.logn_comp ;
419  dcov_logn = star.dcov_logn ;
420  dcon_logn = star.dcon_logn ;
421  lnq_auto = star.lnq_auto ;
422  lnq_comp = star.lnq_comp ;
423  psi4 = star.psi4 ;
424  dcov_phi = star.dcov_phi ;
425  dcon_phi = star.dcon_phi ;
426  flat = star.flat ;
427  gtilde = star.gtilde ;
428  beta_auto = star.beta_auto ;
429  beta_comp = star.beta_comp ;
430  hij = star.hij ;
431  hij_auto = star.hij_auto ;
432  hij_comp = star.hij_comp ;
433  tkij_auto = star.tkij_auto ;
434  tkij_comp = star.tkij_comp ;
435  kcar_auto = star.kcar_auto ;
436  kcar_comp = star.kcar_comp ;
437  ssjm1_logn = star.ssjm1_logn ;
438  ssjm1_lnq = star.ssjm1_lnq ;
439  ssjm1_khi = star.ssjm1_khi ;
440  ssjm1_wbeta = star.ssjm1_wbeta ;
441  ssjm1_h11 = star.ssjm1_h11 ;
442  ssjm1_h21 = star.ssjm1_h21 ;
443  ssjm1_h31 = star.ssjm1_h31 ;
444  ssjm1_h22 = star.ssjm1_h22 ;
445  ssjm1_h32 = star.ssjm1_h32 ;
446  ssjm1_h33 = star.ssjm1_h33 ;
447  decouple = star.decouple ;
448  conf_flat = star.conf_flat ;
449 
450  del_deriv() ; // Deletes all derived quantities
451 
452 }
453 
455 
456  del_deriv() ; // sets to 0x0 all the derived quantities
457  return pot_centri ;
458 
459 }
460 
462 
463  del_deriv() ; // sets to 0x0 all the derived quantities
464  return logn_comp ;
465 
466 }
467 
469 
470  del_deriv() ; // sets to 0x0 all the derived quantities
471  return beta_auto ;
472 
473 }
474 
476 
477  del_deriv() ; // sets to 0x0 all the derived quantities
478  return beta ;
479 
480 }
481 
482 
483  //--------------//
484  // Outputs //
485  //--------------//
486 
487 // Save in a file
488 // --------------
489 void Star_bin::sauve(FILE* fich) const {
490 
491  Star::sauve(fich) ;
492 
493  logn_auto.sauve(fich) ;
494  lnq_auto.sauve(fich) ;
495  beta_auto.sauve(fich) ;
496  hij_auto.sauve(fich) ;
497 
498  ssjm1_logn.sauve(fich) ;
499  ssjm1_lnq.sauve(fich) ;
500  ssjm1_khi.sauve(fich) ;
501  ssjm1_wbeta.sauve(fich) ;
502  ssjm1_h11.sauve(fich) ;
503  ssjm1_h21.sauve(fich) ;
504  ssjm1_h31.sauve(fich) ;
505  ssjm1_h22.sauve(fich) ;
506  ssjm1_h32.sauve(fich) ;
507  ssjm1_h33.sauve(fich) ;
508 
509  fwrite(&irrotational, sizeof(bool), 1, fich) ;
510  fwrite(&conf_flat, sizeof(bool), 1, fich) ;
511 
512  if (irrotational) {
513  gam_euler.sauve(fich) ; // required to construct d_psi from psi0
514  psi0.sauve(fich) ;
515  }
516 
517 }
518 
519 // Printing
520 // --------
521 
522 ostream& Star_bin::operator>>(ostream& ost) const {
523 
524  using namespace Unites ;
525 
526  Star::operator>>(ost) ;
527 
528  ost << endl ;
529  ost << "Star in a binary system" << endl ;
530  ost << "-----------------------" << endl ;
531 
532  if (irrotational) {
533  ost << "irrotational configuration" << endl ;
534  }
535  else {
536  ost << "corotating configuration" << endl ;
537  }
538 
539  ost << "Absolute abscidia of the stellar center: " <<
540  mp.get_ori_x() / km << " km" << endl ;
541 
542  ost << "Absolute abscidia of the barycenter of the baryon density : " <<
543  xa_barycenter() / km << " km" << endl ;
544 
545  double r_0 = 0.5 * ( ray_eq() + ray_eq_pi() ) ;
546  double d_ns = fabs( mp.get_ori_x() ) + ray_eq_pi() - r_0 ;
547  double d_tilde = 2 * d_ns / r_0 ;
548 
549  ost << "d_tilde : " << d_tilde << endl ;
550 
551  ost << "Central value of gam_euler : "
552  << gam_euler.val_grid_point(0, 0, 0, 0) << endl ;
553 
554  ost << "Central u_euler (U^r, U^t, U^p) [c] : "
555  << u_euler(1).val_grid_point(0, 0, 0, 0) << " "
556  << u_euler(2).val_grid_point(0, 0, 0, 0) << " "
557  << u_euler(3).val_grid_point(0, 0, 0, 0) << endl ;
558 
559  if (irrotational) {
560  ost << "Central d_psi (r, t, p) [c] : "
561  << d_psi(1).val_grid_point(0, 0, 0, 0) << " "
562  << d_psi(2).val_grid_point(0, 0, 0, 0) << " "
563  << d_psi(3).val_grid_point(0, 0, 0, 0) << endl ;
564 
565  ost << "Central vel. / co-orb. (W^r, W^t, W^p) [c] : "
566  << wit_w(1).val_grid_point(0, 0, 0, 0) << " "
567  << wit_w(2).val_grid_point(0, 0, 0, 0) << " "
568  << wit_w(3).val_grid_point(0, 0, 0, 0) << endl ;
569 
570  ost << "Max vel. / co-orb. (W^r, W^t, W^p) [c] : "
571  << max(max(wit_w(1))) << " "
572  << max(max(wit_w(2))) << " "
573  << max(max(wit_w(3))) << endl ;
574 
575  ost << "Min vel. / co-orb. (W^r, W^t, W^p) [c] : "
576  << min(min(wit_w(1))) << " "
577  << min(min(wit_w(2))) << " "
578  << min(min(wit_w(3))) << endl ;
579 
580  double r_surf = mp.val_r(0,1.,M_PI/4,M_PI/4) ;
581 
582  ost << "Velocity at (r_surf,pi/4,pi/4) / co-orb. [c] : "
583  << wit_w(1).val_point(r_surf,M_PI/4,M_PI/4) << " "
584  << wit_w(2).val_point(r_surf,M_PI/4,M_PI/4) << " "
585  << wit_w(3).val_point(r_surf,M_PI/4,M_PI/4) << endl ;
586 
587  ost << "Central value of loggam : "
588  << loggam.val_grid_point(0, 0, 0, 0) << endl ;
589  }
590 
591 
592  ost << "Central value of log(N) auto, comp : "
593  << logn_auto.val_grid_point(0, 0, 0, 0) << " "
594  << logn_comp.val_grid_point(0, 0, 0, 0) << endl ;
595 
596  ost << "Central value of beta (N^r, N^t, N^p) [c] : "
597  << beta(1).val_grid_point(0, 0, 0, 0) << " "
598  << beta(2).val_grid_point(0, 0, 0, 0) << " "
599  << beta(3).val_grid_point(0, 0, 0, 0) << endl ;
600 
601  ost << " ... beta_auto part of it [c] : "
602  << beta_auto(1).val_grid_point(0, 0, 0, 0) << " "
603  << beta_auto(2).val_grid_point(0, 0, 0, 0) << " "
604  << beta_auto(3).val_grid_point(0, 0, 0, 0) << endl ;
605 
606  ost << endl << "Central value of (B^r, B^t, B^p)/N [c] : "
607  << bsn(1).val_grid_point(0, 0, 0, 0) << " "
608  << bsn(2).val_grid_point(0, 0, 0, 0) << " "
609  << bsn(3).val_grid_point(0, 0, 0, 0) << endl ;
610 
611 
612  ost << endl << "Central A^{ij} [c/km] : " << endl ;
613  ost << " A^{xx} auto, comp : "
614  << tkij_auto(1, 1).val_grid_point(0, 0, 0, 0) * km << " "
615  << tkij_comp(1, 1).val_grid_point(0, 0, 0, 0) * km << endl ;
616  ost << " A^{xy} auto, comp : "
617  << tkij_auto(1, 2).val_grid_point(0, 0, 0, 0) * km << " "
618  << tkij_comp(1, 2).val_grid_point(0, 0, 0, 0) * km << endl ;
619  ost << " A^{xz} auto, comp : "
620  << tkij_auto(1, 3).val_grid_point(0, 0, 0, 0) * km << " "
621  << tkij_comp(1, 3).val_grid_point(0, 0, 0, 0) * km << endl ;
622  ost << " A^{yy} auto, comp : "
623  << tkij_auto(2, 2).val_grid_point(0, 0, 0, 0) * km << " "
624  << tkij_comp(2, 2).val_grid_point(0, 0, 0, 0) * km << endl ;
625  ost << " A^{yz} auto, comp : "
626  << tkij_auto(2, 3).val_grid_point(0, 0, 0, 0) * km << " "
627  << tkij_comp(2, 3).val_grid_point(0, 0, 0, 0) * km << endl ;
628  ost << " A^{zz} auto, comp : "
629  << tkij_auto(3, 3).val_grid_point(0, 0, 0, 0) * km << " "
630  << tkij_comp(3, 3).val_grid_point(0, 0, 0, 0) * km << endl ;
631 
632  ost << endl << "Central A_{ij} A^{ij} [c^2/km^2] : " << endl ;
633  ost << " A_{ij} A^{ij} auto, comp : "
634  << kcar_auto.val_grid_point(0, 0, 0, 0) * km*km << " "
635  << kcar_comp.val_grid_point(0, 0, 0, 0) * km*km << endl ;
636 
637 
638  return ost ;
639 }
640 
641  //-------------------------//
642  // Computational routines //
643  //-------------------------//
644 
646 
647  if (!irrotational) {
649  return ;
650  }
651 
652  // Specific relativistic enthalpy ---> hhh
653  //----------------------------------
654 
655  Scalar hhh = exp(ent) ; // = 1 at the Newtonian limit
656 
657  // Computation of W^i = - h Gamma_n B^i/N
658  //----------------------------------------------
659 
660  Vector www = hhh * gam_euler * bsn * psi4 ;
661 
662  // Constant value of W^i at the center of the star
663  //-------------------------------------------------
664 
665  Vector v_orb(mp, COV, mp.get_bvect_cart()) ;
666 
667  for (int i=1; i<=3; i++) {
668  v_orb.set(i) = www(i).val_grid_point(0, 0, 0, 0) ;
669  }
670 
671  // Gradient of psi
672  //----------------
673 
674  Vector d_psi0 = psi0.derive_cov(flat) ;
675 
676  d_psi0.change_triad( mp.get_bvect_cart() ) ;
677  d_psi0.std_spectral_base() ;
678 
679  d_psi = d_psi0 + v_orb ;
680  for (int i=1; i<=3; i++) {
681  if (d_psi(i).get_etat() == ETATZERO)
682  d_psi.set(i).annule_hard() ;
683  }
684 
685  // C^1 continuation of d_psi outside the star
686  // (to ensure a smooth enthalpy field accross the stellar surface)
687  // ----------------------------------------------------------------
688 
689  int nzm1 = mp.get_mg()->get_nzone() - 1 ;
690  d_psi.annule(nzet, nzm1) ;
691  for (int i=1; i<=3; i++) {
692  Cmp d_psi_i (d_psi(i)) ;
693  d_psi_i.va.set_base( d_psi0(i).get_spectral_va().base ) ;
694  d_psi_i = raccord_c1(d_psi_i, nzet) ;
695  d_psi.set(i) = d_psi_i ;
696  }
697 
698 }
699 
700 
701 void Star_bin::relaxation(const Star_bin& star_jm1, double relax_ent,
702  double relax_met, int mer, int fmer_met) {
703 
704  double relax_ent_jm1 = 1. - relax_ent ;
705  double relax_met_jm1 = 1. - relax_met ;
706 
707  ent = relax_ent * ent + relax_ent_jm1 * star_jm1.ent ;
708 
709  if ( (mer != 0) && (mer % fmer_met == 0)) {
710 
711  logn_auto = relax_met * logn_auto + relax_met_jm1 * star_jm1.logn_auto ;
712  lnq_auto = relax_met * lnq_auto + relax_met_jm1 * star_jm1.lnq_auto ;
713 
714  beta_auto = relax_met * beta_auto
715  + relax_met_jm1 * star_jm1.beta_auto ;
716 
717  hij_auto = relax_met * hij_auto + relax_met_jm1 * star_jm1.hij_auto ;
718 
719  }
720 
721  del_deriv() ;
722 
723  equation_of_state() ;
724 
725 }
726 
727 void Star_bin::test_K_Hi() const {
728 
729  int nr = mp.get_mg()->get_nr(0) ;
730  int nt = mp.get_mg()->get_nt(0) ;
731  int np = mp.get_mg()->get_np(0) ;
732 
733  cout << "La jauge de Dirac est elle bien satisfaite ??" << endl ;
734  cout << "Vector Hi" << endl ;
735  for (int i=1; i<=3; i++)
736  cout << " Comp. " << i << " : " << norme((gtilde.con()
737  .divergence(flat))(i)/(nr*nt*np)) << endl ;
738 
739 
740  cout << "Pour comparaison valeur de D_i(g^1i)" << endl ;
741  for (int i=1; i<=3; i++)
742  cout << " i = " << i << " : " << norme((gtilde.con().derive_cov(flat))
743  (1, i, i)/(nr*nt*np)) << endl ;
744 
745 
746 
747 }
748 }
Sym_tensor hij_comp
Deviation of the inverse conformal metric from the inverse flat metric generated principally by the ...
Definition: star.h:594
Vector dcov_phi
Covariant derivative of the logarithm of the conformal factor.
Definition: star.h:555
bool irrotational
true for an irrotational star, false for a corotating one
Definition: star.h:491
Metric for tensor calculation.
Definition: metric.h:90
virtual void sauve(FILE *) const
Save in a binary file.
Definition: tensor_sym.C:372
virtual const Sym_tensor & con() const
Read-only access to the contravariant representation.
Definition: metric.C:290
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: star.C:316
void operator=(const Star_bin &)
Assignment to another Star_bin.
Definition: star_bin.C:404
Scalar psi0
Scalar potential of the non-translational part of fluid 4-velocity (in the irrotational case) ...
Definition: star.h:496
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition: cmp.h:446
Cmp exp(const Cmp &)
Exponential.
Definition: cmp_math.C:270
Vector wit_w
Spatial projection of the fluid 3-velocity with respect to the co-orbiting observer.
Definition: star.h:507
Map & mp
Mapping associated with the star.
Definition: star.h:180
Scalar ssjm1_lnq
Effective source at the previous step for the resolution of the Poisson equation for lnq_auto...
Definition: star.h:628
Sym_tensor hij_auto
Deviation of the inverse conformal metric from the inverse flat metric generated principally by the ...
Definition: star.h:588
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
Definition: grilles.h:462
Metric gamma
3-metric
Definition: star.h:235
Vector d_psi
Gradient of (in the irrotational case) (Spherical components with respect to the mapping of the star...
Definition: star.h:501
Scalar ssjm1_h31
Effective source at the previous step for the resolution of the Poisson equation for h20_auto...
Definition: star.h:651
Lorene prototypes.
Definition: app_hor.h:64
Scalar ssjm1_h11
Effective source at the previous step for the resolution of the Poisson equation for h00_auto...
Definition: star.h:641
Standard units of space, time and mass.
Equation of state base class.
Definition: eos.h:190
Sym_tensor stress_euler
Spatial part of the stress-energy tensor with respect to the Eulerian observer.
Definition: star.h:212
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:765
Tensor field of valence 0 (or component of a tensorial field).
Definition: scalar.h:387
virtual void set_etat_nondef()
Sets the logical state of all components to ETATNONDEF (undefined state).
Definition: tensor.C:489
Scalar & set_logn_comp()
Read/write of the logarithm of the lapse generated principally by the companion.
Definition: star_bin.C:461
Base class for coordinate mappings.
Definition: map.h:670
double get_ori_x() const
Returns the x coordinate of the origin.
Definition: map.h:768
Scalar ssjm1_h33
Effective source at the previous step for the resolution of the Poisson equation for h22_auto...
Definition: star.h:666
Scalar logn_comp
Part of the lapse logarithm (gravitational potential at the Newtonian limit) generated principally by...
Definition: star.h:532
Scalar loggam
Logarithm of the Lorentz factor between the fluid and the co-orbiting observer.
Definition: star.h:512
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: star_bin.C:381
Scalar kcar_auto
Part of the scalar generated by beta_auto, i.e.
Definition: star.h:612
virtual void sauve(FILE *) const
Save in a binary file.
Definition: tensor.C:906
virtual void change_triad(const Base_vect &)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
Scalar ent
Log-enthalpy.
Definition: star.h:190
Base class for stars.
Definition: star.h:175
double * p_xa_barycenter
Absolute coordinate X of the barycenter of the baryon density.
Definition: star.h:687
Tbl min(const Cmp &)
Minimum values of a Cmp in each domain.
Definition: cmp_math.C:458
Vector beta
Shift vector.
Definition: star.h:228
Tensor field of valence 1.
Definition: vector.h:188
Scalar logn_auto
Part of the lapse logarithm (gravitational potential at the Newtonian limit) generated principally by...
Definition: star.h:527
const Tensor_sym & derive_cov(const Metric &gam) const
Returns the covariant derivative of this with respect to some metric .
void set_base(const Base_val &)
Sets the bases for spectral expansions (member base )
Definition: valeur.C:810
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
Scalar ssjm1_h21
Effective source at the previous step for the resolution of the Poisson equation for h10_auto...
Definition: star.h:646
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 std_spectral_base()
Sets the standard spectal bases of decomposition for each component.
Definition: vector.C:316
void annule_hard()
Sets the Scalar to zero in a hard way.
Definition: scalar.C:380
void operator=(const Star &)
Assignment to another Star.
Definition: star.C:351
int nzet
Number of domains of *mp occupied by the star.
Definition: star.h:183
void test_K_Hi() const
Test if the gauge conditions we impose are well satisfied.
Definition: star_bin.C:727
Scalar gam_euler
Lorentz factor between the fluid and Eulerian observers.
Definition: star.h:204
virtual void sauve(FILE *) const
Save in a file.
Definition: scalar.C:686
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.
Scalar decouple
Function used to construct the part generated by the star from the total .
Definition: star.h:676
const Vector & divergence(const Metric &) const
Returns the divergence of this with respect to a Metric .
Definition: sym_tensor.C:349
Star_bin(Map &mp_i, int nzet_i, const Eos &eos_i, bool irrot, bool conf_flat)
Standard constructor.
Definition: star_bin.C:115
Tbl norme(const Cmp &)
Sums of the absolute values of all the values of the Cmp in each domain.
Definition: cmp_math.C:481
Sym_tensor hij
Total deviation of the inverse conformal metric from the inverse flat metric.
Definition: star.h:581
Scalar pot_centri
Centrifugal potential.
Definition: star.h:521
virtual void del_hydro_euler()
Sets to ETATNONDEF (undefined state) the hydrodynamical quantities relative to the Eulerian observer...
Definition: star_bin.C:389
virtual ~Star_bin()
Destructor.
Definition: star_bin.C:359
Metric_flat flat
Flat metric defined on the mapping (Spherical components with respect to the mapping of the star) ...
Definition: star.h:562
Scalar lnq_auto
Scalar field generated principally by the star.
Definition: star.h:543
Scalar lnq_comp
Scalar field generated principally by the companion star.
Definition: star.h:548
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:448
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: star_bin.C:369
Tbl max(const Cmp &)
Maximum values of a Cmp in each domain.
Definition: cmp_math.C:435
Vector u_euler
Fluid 3-velocity with respect to the Eulerian observer.
Definition: star.h:207
Sym_tensor tkij_comp
Part of the extrinsic curvature tensor generated by beta_comp.
Definition: star.h:606
double ray_eq_pi() const
Coordinate radius at , [r_unit].
Definition: star_global.C:186
virtual void sauve(FILE *) const
Save in a file.
Definition: star_bin.C:489
void fait_d_psi()
Computes the gradient of the total velocity potential .
Definition: star_bin.C:645
Vector dcov_logn
Covariant derivative of the total logarithm of the lapse.
Definition: star.h:535
Class for stars in binary system.
Definition: star.h:483
Metric gtilde
Conformal metric .
Definition: star.h:565
Vector dcon_logn
Contravariant derivative of the total logarithm of the lapse.
Definition: star.h:538
int get_nr(int l) const
Returns the number of points in the radial direction ( ) in domain no. l.
Definition: grilles.h:452
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: star.C:298
virtual void del_hydro_euler()
Sets to ETATNONDEF (undefined state) the hydrodynamical quantities relative to the Eulerian observer...
Definition: star.C:330
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:791
Scalar ssjm1_h22
Effective source at the previous step for the resolution of the Poisson equation for h11_auto...
Definition: star.h:656
Vector beta_comp
Part of the shift vector generated principally by the star (Spherical components with respect to the ...
Definition: star.h:575
Sym_tensor tkij_auto
Part of the extrinsic curvature tensor generated by beta_auto.
Definition: star.h:600
bool conf_flat
true if the 3-metric is conformally flat, false for a more general metric.
Definition: star.h:681
double ray_eq() const
Coordinate radius at , [r_unit].
Definition: star_global.C:108
Scalar ssjm1_h32
Effective source at the previous step for the resolution of the Poisson equation for h21_auto...
Definition: star.h:661
Vector beta_auto
Part of the shift vector generated principally by the star (Spherical components with respect to the ...
Definition: star.h:570
Scalar kcar_comp
Part of the scalar generated by beta_auto and beta_comp, i.e.
Definition: star.h:618
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition: star_bin.C:522
virtual void change_triad(const Base_vect &new_triad)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition: star.C:417
Scalar & set_pot_centri()
Read/write the centrifugal potential.
Definition: star_bin.C:454
virtual void sauve(FILE *) const
Save in a file.
Definition: star.C:397
virtual void set_etat_zero()
Sets the logical state of all components to ETATZERO (zero state).
Definition: tensor.C:497
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
Definition: grilles.h:457
void relaxation(const Star_bin &star_prev, double relax_ent, double relax_met, int mer, int fmer_met)
Performs a relaxation on ent, logn_auto, lnq_auto, beta_auto and hij_auto.
Definition: star_bin.C:701
Scalar ssjm1_khi
Effective source at the previous step for the resolution of the Poisson equation for khi...
Definition: star.h:634
Scalar & set(int)
Read/write access to a component.
Definition: vector.C:296
Vector & set_beta()
Read/write of .
Definition: star_bin.C:475
void equation_of_state()
Computes the proper baryon and energy density, as well as pressure from the enthalpy.
Definition: star.C:462
const Vector & derive_cov(const Metric &gam) const
Returns the gradient (1-form = covariant vector) of *this.
Definition: scalar_deriv.C:390
virtual void annule(int l_min, int l_max)
Sets the Tensor to zero in several domains.
Definition: tensor.C:671
Vector bsn
3-vector shift, divided by N, of the rotating coordinates, .
Definition: star.h:518
Valeur va
The numerical value of the Cmp.
Definition: cmp.h:464
Scalar psi4
Conformal factor .
Definition: star.h:552
Vector dcon_phi
Contravariant derivative of the logarithm of the conformal factor.
Definition: star.h:557
Vector & set_beta_auto()
Read/write of .
Definition: star_bin.C:468
virtual double xa_barycenter() const
Absolute coordinate X of the barycenter of the baryon density,.
Scalar ssjm1_logn
Effective source at the previous step for the resolution of the Poisson equation for logn_auto...
Definition: star.h:623