LORENE
star_bhns.C
1 /*
2  * Methods of class Star_bhns
3  *
4  * (see file star_bhns.h for documentation).
5  *
6  */
7 
8 /*
9  * Copyright (c) 2005-2007 Keisuke Taniguchi
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 version 2
15  * as published by the Free Software Foundation.
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 char star_bhns_C[] = "$Header: /cvsroot/Lorene/C++/Source/Star_bhns/star_bhns.C,v 1.4 2014/10/13 08:53:40 j_novak Exp $" ;
29 
30 /*
31  * $Id: star_bhns.C,v 1.4 2014/10/13 08:53:40 j_novak Exp $
32  * $Log: star_bhns.C,v $
33  * Revision 1.4 2014/10/13 08:53:40 j_novak
34  * Lorene classes and functions now belong to the namespace Lorene.
35  *
36  * Revision 1.3 2014/10/06 15:13:16 j_novak
37  * Modified #include directives to use c++ syntax.
38  *
39  * Revision 1.2 2008/05/15 19:12:38 k_taniguchi
40  * Change of some parameters.
41  *
42  * Revision 1.1 2007/06/22 01:30:10 k_taniguchi
43  * *** empty log message ***
44  *
45  *
46  * $Header: /cvsroot/Lorene/C++/Source/Star_bhns/star_bhns.C,v 1.4 2014/10/13 08:53:40 j_novak Exp $
47  *
48  */
49 
50 // C++ headers
51 //#include <>
52 
53 // C headers
54 #include <cmath>
55 
56 // Lorene headers
57 #include "etoile.h"
58 #include "star.h"
59 #include "star_bhns.h"
60 #include "eos.h"
61 #include "unites.h"
62 
63 // Local prototype
64 namespace Lorene {
65 Cmp raccord_c1(const Cmp& uu, int l1) ;
66 
67  //---------------------//
68  // Constructor //
69  //---------------------//
70 
71 // Standard constructor
72 // --------------------
73 Star_bhns::Star_bhns(Map& mp_i, int nzet_i, const Eos& eos_i, bool irrot_i)
74  : Star(mp_i, nzet_i, eos_i),
75  mp_aff(mp_i),
76  irrotational(irrot_i),
77  psi0(mp_i),
78  d_psi(mp_i, COV, mp_i.get_bvect_cart()),
79  wit_w(mp_i, CON, mp_i.get_bvect_cart()),
80  loggam(mp_i),
81  bsn(mp_i, CON, mp_i.get_bvect_cart()),
82  gam(mp_i),
83  gam0(mp_i),
84  pot_centri(mp_i),
85  lapconf_auto(mp_i),
86  lapconf_comp(mp_i),
87  lapconf_tot(mp_i),
88  lapse_auto(mp_i),
89  lapse_tot(mp_i),
90  d_lapconf_auto(mp_i, COV, mp_i.get_bvect_cart()),
91  d_lapconf_comp(mp_i, COV, mp_i.get_bvect_cart()),
92  shift_auto(mp_i, CON, mp_i.get_bvect_cart()),
93  shift_comp(mp_i, CON, mp_i.get_bvect_cart()),
94  shift_tot(mp_i, CON, mp_i.get_bvect_cart()),
95  d_shift_auto(mp_i, 2, CON, mp_i.get_bvect_cart()),
96  d_shift_comp(mp_i, 2, CON, mp_i.get_bvect_cart()),
97  confo_auto(mp_i),
98  confo_comp(mp_i),
99  confo_tot(mp_i),
100  d_confo_auto(mp_i, COV, mp_i.get_bvect_cart()),
101  d_confo_comp(mp_i, COV, mp_i.get_bvect_cart()),
102  psi4(mp_i),
103  taij_auto(mp_i, CON, mp_i.get_bvect_cart()),
104  taij_quad_auto(mp_i),
105  flat(mp_i, mp_i.get_bvect_cart()),
106  ssjm1_lapconf(mp_i),
107  ssjm1_confo(mp_i),
108  ssjm1_khi(mp_i),
109  ssjm1_wshift(mp_i, CON, mp_i.get_bvect_cart())
110 {
111 
112  // Pointers of derived quantities initialized to zero :
113  set_der_0x0() ;
114 
115  // Quantities defined on a spherical triad in star.C are put on
116  // a cartesian one
118 
119  // All quantities are initialized to zero
120  psi0 = 0. ;
122  d_psi.set_etat_zero() ;
123  wit_w.set_etat_zero() ;
124  loggam = 0. ;
126  bsn.set_etat_zero() ;
127  gam = 0. ;
129  gam0 = 0. ;
131  pot_centri = 0. ;
133 
134  lapconf_auto = 1. ;
136  lapconf_comp = 0. ;
138  lapconf_tot = 1. ;
140  lapse_auto = 1. ;
142  lapse_tot = 1. ;
151  confo_auto = 1. ;
153  confo_comp = 0. ;
155  confo_tot = 1. ;
159  psi4 = 1. ;
161 
163  taij_quad_auto = 0. ;
165 
166  ssjm1_lapconf = 0. ;
168  ssjm1_confo = 0. ;
170  ssjm1_khi = 0. ;
173 
174 }
175 
176 // Copy constructor
177 // ----------------
179  : Star(star),
180  mp_aff(star.mp_aff),
182  psi0(star.psi0),
183  d_psi(star.d_psi),
184  wit_w(star.wit_w),
185  loggam(star.loggam),
186  bsn(star.bsn),
187  gam(star.gam),
188  gam0(star.gam0),
189  pot_centri(star.pot_centri),
192  lapconf_tot(star.lapconf_tot),
193  lapse_auto(star.lapse_auto),
194  lapse_tot(star.lapse_tot),
197  shift_auto(star.shift_auto),
198  shift_comp(star.shift_comp),
199  shift_tot(star.shift_tot),
202  confo_auto(star.confo_auto),
203  confo_comp(star.confo_comp),
204  confo_tot(star.confo_tot),
207  psi4(star.psi4),
208  taij_auto(star.taij_auto),
210  flat(star.flat),
212  ssjm1_confo(star.ssjm1_confo),
213  ssjm1_khi(star.ssjm1_khi),
215 {
216  set_der_0x0() ;
217 }
218 
219 // Constructor from a file
220 // -----------------------
221 Star_bhns::Star_bhns(Map& mp_i, const Eos& eos_i, FILE* fich)
222  : Star(mp_i, eos_i, fich),
223  mp_aff(mp_i),
224  psi0(mp_i),
225  d_psi(mp_i, COV, mp_i.get_bvect_cart()),
226  wit_w(mp_i, CON, mp_i.get_bvect_cart()),
227  loggam(mp_i),
228  bsn(mp_i, CON, mp_i.get_bvect_cart()),
229  gam(mp_i),
230  gam0(mp_i),
231  pot_centri(mp_i),
232  lapconf_auto(mp_i, *(mp_i.get_mg()), fich),
233  lapconf_comp(mp_i),
234  lapconf_tot(mp_i),
235  lapse_auto(mp_i),
236  lapse_tot(mp_i),
237  d_lapconf_auto(mp_i, COV, mp_i.get_bvect_cart()),
238  d_lapconf_comp(mp_i, COV, mp_i.get_bvect_cart()),
239  shift_auto(mp_i, mp_i.get_bvect_cart(), fich),
240  shift_comp(mp_i, CON, mp_i.get_bvect_cart()),
241  shift_tot(mp_i, CON, mp_i.get_bvect_cart()),
242  d_shift_auto(mp_i, 2, CON, mp_i.get_bvect_cart()),
243  d_shift_comp(mp_i, 2, CON, mp_i.get_bvect_cart()),
244  confo_auto(mp_i, *(mp_i.get_mg()), fich),
245  confo_comp(mp_i),
246  confo_tot(mp_i),
247  d_confo_auto(mp_i, COV, mp_i.get_bvect_cart()),
248  d_confo_comp(mp_i, COV, mp_i.get_bvect_cart()),
249  psi4(mp_i),
250  taij_auto(mp_i, CON, mp_i.get_bvect_cart()),
251  taij_quad_auto(mp_i),
252  flat(mp_i, mp_i.get_bvect_cart()),
253  ssjm1_lapconf(mp_i, *(mp_i.get_mg()), fich),
254  ssjm1_confo(mp_i, *(mp_i.get_mg()), fich),
255  ssjm1_khi(mp_i, *(mp_i.get_mg()), fich),
256  ssjm1_wshift(mp_i, mp_i.get_bvect_cart(), fich)
257 {
258 
259  // Star parameter
260  fread(&irrotational, sizeof(bool), 1, fich) ;
261 
262  // Read of the saved fields
263  // ------------------------
264 
265  if (irrotational) {
266  Scalar gam_euler_file(mp, *(mp.get_mg()), fich) ;
267  gam_euler = gam_euler_file ;
268 
269  Scalar psi0_file(mp, *(mp.get_mg()), fich) ;
270  psi0 = psi0_file ;
271  }
272 
273  // Quantities defined on a spherical triad in star.C are put on
274  // a cartesian one
276 
277  // All other quantities are initialized to zero
278  // --------------------------------------------
279 
280  d_psi.set_etat_zero() ;
281  wit_w.set_etat_zero() ;
282  loggam = 0. ;
284  bsn.set_etat_zero() ;
285  gam = 0. ;
287  gam0 = 0. ;
289  pot_centri = 0. ;
291 
292  lapconf_comp = 0. ;
294  lapconf_tot = 1. ;
296  lapse_auto = 1. ;
298  lapse_tot = 1. ;
306  confo_comp = 0. ;
308  confo_tot = 1. ;
312  psi4 = 1. ;
315  taij_quad_auto = 0. ;
317 
318  // Pointers of derived quantities initialized to zero
319  // --------------------------------------------------
320  set_der_0x0() ;
321 
322 }
323 
324 
325  //--------------------//
326  // Destructor //
327  //--------------------//
328 
330 {
331 
332  del_deriv() ;
333 
334 }
335 
336 
337  //------------------------------------------//
338  // Management of derived quantities //
339  //------------------------------------------//
340 
341 void Star_bhns::del_deriv() const {
342 
343  Star::del_deriv() ;
344 
345  if (p_mass_b_bhns != 0x0) delete p_mass_b_bhns ;
346  if (p_mass_g_bhns != 0x0) delete p_mass_g_bhns ;
347 
348  set_der_0x0() ;
349 
350 }
351 
353 
355 
356  p_mass_b_bhns = 0x0 ;
357  p_mass_g_bhns = 0x0 ;
358 
359 }
360 
361 
362  //--------------------//
363  // Assignment //
364  //--------------------//
365 
366 // Assignment to another Star_bhns
367 // -------------------------------
368 void Star_bhns::operator=(const Star_bhns& star) {
369 
370  // Assignment of quantities common to the derived classes of Star
371  Star::operator=(star) ;
372 
373  // Assignment of proper quantities of class Star_bhns
374  mp_aff = star.mp_aff ;
375  irrotational = star.irrotational ;
376  psi0 = star.psi0 ;
377  d_psi = star.d_psi ;
378  wit_w = star.wit_w ;
379  loggam = star.loggam ;
380  bsn = star.bsn ;
381  gam = star.gam ;
382  gam0 = star.gam0 ;
383  pot_centri = star.pot_centri ;
384  lapconf_auto = star.lapconf_auto ;
385  lapconf_comp = star.lapconf_comp ;
386  lapconf_tot = star.lapconf_tot ;
387  lapse_auto = star.lapse_auto ;
388  lapse_tot = star.lapse_tot ;
391  shift_auto = star.shift_auto ;
392  shift_comp = star.shift_comp ;
393  shift_tot = star.shift_tot ;
394  d_shift_auto = star.d_shift_auto ;
395  d_shift_comp = star.d_shift_comp ;
396  confo_auto = star.confo_auto ;
397  confo_comp = star.confo_comp ;
398  confo_tot = star.confo_tot ;
399  d_confo_auto = star.d_confo_auto ;
400  d_confo_comp = star.d_confo_comp ;
401  psi4 = star.psi4 ;
402  taij_auto = star.taij_auto ;
404  flat = star.flat ;
406  ssjm1_confo = star.ssjm1_confo ;
407  ssjm1_khi = star.ssjm1_khi ;
408  ssjm1_wshift = star.ssjm1_wshift ;
409 
410  del_deriv() ;
411 
412 }
413 
415 
416  del_deriv() ;
417  return pot_centri ;
418 
419 }
420 
422 
423  del_deriv() ;
424  return lapconf_auto ;
425 
426 }
427 
429 
430  del_deriv() ;
431  return lapconf_comp ;
432 
433 }
434 
436 
437  del_deriv() ;
438  return shift_auto ;
439 
440 }
441 
443 
444  del_deriv() ;
445  return shift_comp ;
446 
447 }
448 
450 
451  del_deriv() ;
452  return confo_auto ;
453 
454 }
455 
457 
458  del_deriv() ;
459  return confo_comp ;
460 
461 }
462 
463 
464  //-----------------//
465  // Outputs //
466  //-----------------//
467 
468 // Save in a file
469 // --------------
470 void Star_bhns::sauve(FILE* fich) const {
471 
472  Star::sauve(fich) ;
473 
474  lapconf_auto.sauve(fich) ;
475  shift_auto.sauve(fich) ;
476  confo_auto.sauve(fich) ;
477 
478  ssjm1_lapconf.sauve(fich) ;
479  ssjm1_confo.sauve(fich) ;
480  ssjm1_khi.sauve(fich) ;
481  ssjm1_wshift.sauve(fich) ;
482 
483  fwrite(&irrotational, sizeof(bool), 1, fich) ;
484 
485  if (irrotational) {
486  gam_euler.sauve(fich) ; // required to construct d_psi from psi0
487  psi0.sauve(fich) ;
488  }
489 
490 }
491 
492 // Printing
493 // --------
494 
495 ostream& Star_bhns::operator>>(ostream& ost) const {
496 
497  using namespace Unites ;
498 
499  // Star::operator>>(ost) ;
500 
501  ost << endl ;
502  ost << "Neutron star in a BHNS binary" << endl ;
503  ost << "-----------------------------" << endl ;
504 
505  ost << "Coordinate radius R_eq_tow : "
506  << ray_eq_pi() / km << " [km]" << endl ;
507  ost << "Coordinate radius R_eq_opp : "
508  << ray_eq() / km << " [km]" << endl ;
509  ost << "Coordinate radius R_eq_orb : "
510  << ray_eq_pis2() / km << " [km]" << endl ;
511  ost << "Coordinate radius R_pole : "
512  << ray_pole() / km << " [km]" << endl ;
513  ost << "Central enthalph H_ent : "
514  << ent.val_grid_point(0,0,0,0) << endl ;
515  ost << "Lapse function at the center of NS : "
516  << lapse_tot.val_grid_point(0,0,0,0) << endl ;
517  ost << "Conformal factor at the center of NS : "
518  << confo_tot.val_grid_point(0,0,0,0) << endl ;
519  ost << "shift(1) at the center of NS : "
520  << shift_tot(1).val_grid_point(0,0,0,0) << endl ;
521  ost << "shift(2) at the center of NS : "
522  << shift_tot(2).val_grid_point(0,0,0,0) << endl ;
523  ost << "shift(3) at the center of NS : "
524  << shift_tot(3).val_grid_point(0,0,0,0) << endl ;
525 
526  return ost ;
527 
528 }
529 
530 
531  //--------------------------------//
532  // Computational routines //
533  //--------------------------------//
534 
536 
537  if (!irrotational) {
539  return ;
540  }
541 
542  // Specific relativistic enthalpy ---> hhh
543  // ------------------------------
544 
545  Scalar hhh = exp(ent) ; // = 1 at the Newtonian limit
546 
547  // Computation of W^i = h Gamma_n B^i / N
548  // --------------------------------------
549 
550  Vector www = hhh * gam_euler * bsn * psi4 ;
551 
552  // Constant value of W^i at the center of the star
553  // -----------------------------------------------
554 
555  Vector v_orb(mp, COV, mp.get_bvect_cart()) ;
556 
557  for (int i=1; i<=3; i++) {
558  v_orb.set(i) = www(i).val_grid_point(0,0,0,0) ;
559  }
560 
561  // Gradient of psi
562  // ---------------
563 
564  Vector d_psi0(mp, COV, mp.get_bvect_cart()) ;
565  d_psi0.set_etat_qcq() ;
566  for (int i=1; i<=3; i++)
567  d_psi0.set(i) = psi0.deriv(i) ;
568 
569  d_psi0.std_spectral_base() ;
570 
571  d_psi = d_psi0 + v_orb ;
572 
573  for (int i=1; i<=3; i++) {
574  if (d_psi(i).get_etat() == ETATZERO)
575  d_psi.set(i).annule_hard() ;
576  }
577 
578  // C^1 continuation of d_psi outside the star
579  // (to ensure a smooth enthalpy field accross the stellar surface)
580  // ---------------------------------------------------------------
581 
582  int nzm1 = mp.get_mg()->get_nzone() - 1 ;
583  d_psi.annule(nzet, nzm1) ;
584  for (int i=1; i<=3; i++) {
585  Cmp d_psi_i (d_psi(i)) ;
586  d_psi_i.va.set_base( d_psi0(i).get_spectral_va().base ) ;
587  d_psi_i = raccord_c1(d_psi_i, nzet) ;
588  d_psi.set(i) = d_psi_i ;
589  }
590 
591 }
592 
593 
594 void Star_bhns::relax_bhns(const Star_bhns& star_prev, double relax_ent,
595  double relax_met, int mer, int fmer_met) {
596 
597  double relax_ent_jm1 = 1. - relax_ent ;
598  double relax_met_jm1 = 1. - relax_met ;
599 
600  ent = relax_ent * ent + relax_ent_jm1 * star_prev.ent ;
601 
602  if ( (mer != 0) && (mer % fmer_met == 0)) {
603 
604  lapconf_auto = relax_met * lapconf_auto
605  + relax_met_jm1 * star_prev.lapconf_auto ;
606 
607  shift_auto = relax_met * shift_auto
608  + relax_met_jm1 * star_prev.shift_auto ;
609 
610  confo_auto = relax_met * confo_auto
611  + relax_met_jm1 * star_prev.confo_auto ;
612 
613  }
614 
615  del_deriv() ;
616 
618 
619 }
620 }
Scalar psi4
Fourth power of the total conformal factor.
Definition: star_bhns.h:176
void fait_d_psi_bhns()
Computes the gradient of the total velocity potential .
Definition: star_bhns.C:535
Vector d_lapconf_auto
Derivative of the lapconf function generated by the star .
Definition: star_bhns.h:130
Vector shift_tot
Total shift vector.
Definition: star_bhns.h:144
Scalar lapconf_tot
Total lapconf function.
Definition: star_bhns.h:119
Class for stars in black hole-neutron star binaries.
Definition: star_bhns.h:59
virtual void set_etat_qcq()
Sets the logical state of all components to ETATQCQ (ordinary state).
Definition: tensor.C:481
virtual void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: star.C:316
double * p_mass_g_bhns
Gravitational mass.
Definition: star_bhns.h:226
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Definition: cmp.h:443
Cmp exp(const Cmp &)
Exponential.
Definition: cmp_math.C:270
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition: star_bhns.C:495
Map & mp
Mapping associated with the star.
Definition: star.h:180
Scalar ssjm1_khi
Effective source at the previous step for the resolution of the Poisson equation for the scalar by m...
Definition: star_bhns.h:210
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: star_bhns.C:341
Vector bsn
3-vector shift, divided by N , of the rotating coordinates, .
Definition: star_bhns.h:99
Lorene prototypes.
Definition: app_hor.h:64
Standard units of space, time and mass.
Equation of state base class.
Definition: eos.h:190
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:765
Vector d_lapconf_comp
Derivative of the lapconf function generated by the companion black hole.
Definition: star_bhns.h:135
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
Vector wit_w
Spatial projection of the fluid 3-velocity with respect to the co-orbiting observer.
Definition: star_bhns.h:88
Base class for coordinate mappings.
Definition: map.h:670
Scalar loggam
Logarithm of the Lorentz factor between the fluid and the co-orbiting observer.
Definition: star_bhns.h:93
double * p_mass_b_bhns
Baryon mass.
Definition: star_bhns.h:225
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
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 & set_pot_centri()
Read/write the centrifugal potential.
Definition: star_bhns.C:414
Scalar ent
Log-enthalpy.
Definition: star.h:190
Base class for stars.
Definition: star.h:175
Vector shift_auto
Shift vector generated by the star.
Definition: star_bhns.h:138
Tensor d_shift_auto
Derivative of the shift vector generated by the star .
Definition: star_bhns.h:149
Tensor field of valence 1.
Definition: vector.h:188
void set_base(const Base_val &)
Sets the bases for spectral expansions (member base )
Definition: valeur.C:810
Scalar lapse_tot
Total lapse function.
Definition: star_bhns.h:125
Scalar gam
Lorentz factor between the fluid and the co-orbiting observer.
Definition: star_bhns.h:102
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
Map_af mp_aff
Affine mapping for solving poisson&#39;s equations of metric quantities.
Definition: star_bhns.h:67
Vector d_psi
Gradient of (in the irrotational case) (Spherical components with respect to the mapping of the star...
Definition: star_bhns.h:82
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
Metric_flat flat
Flat metric defined on the mapping (Spherical components with respect to the mapping of the star )...
Definition: star_bhns.h:192
void operator=(const Star_bhns &)
Assignment to another Star_bhns.
Definition: star_bhns.C:368
Scalar confo_tot
Total conformal factor.
Definition: star_bhns.h:163
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 ~Star_bhns()
Destructor.
Definition: star_bhns.C:329
double ray_eq_pis2() const
Coordinate radius at , [r_unit].
Definition: star_global.C:138
Vector d_confo_comp
Derivative of the conformal factor generated by the companion black hole.
Definition: star_bhns.h:173
Vector & set_shift_comp()
Read/write of the shift vector generated by the companion black hole.
Definition: star_bhns.C:442
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:448
void relax_bhns(const Star_bhns &star_prev, double relax_ent, double relax_met, int mer, int fmer_met)
Performs a relaxation on ent , lapse_auto , shift_auto , confo_auto .
Definition: star_bhns.C:594
Star_bhns(Map &mp_i, int nzet_i, const Eos &eos_i, bool irrot_i)
Standard constructor.
Definition: star_bhns.C:73
Scalar & set_confo_auto()
Read/write of the conformal factor generated by the neutron star.
Definition: star_bhns.C:449
Scalar lapconf_comp
Lapconf function generated by the companion black hole.
Definition: star_bhns.h:116
Vector & set_shift_auto()
Read/write of the shift vector generated by the neutron star.
Definition: star_bhns.C:435
Vector u_euler
Fluid 3-velocity with respect to the Eulerian observer.
Definition: star.h:207
const Scalar & deriv(int i) const
Returns of *this , where .
Definition: scalar_deriv.C:359
double ray_eq_pi() const
Coordinate radius at , [r_unit].
Definition: star_global.C:186
Scalar taij_quad_auto
Part of the scalar generated by .
Definition: star_bhns.h:187
Vector shift_comp
Shift vector generated by the companion black hole.
Definition: star_bhns.h:141
Scalar gam0
Lorentz factor between the co-orbiting observer and the Eulerian one.
Definition: star_bhns.h:107
Vector ssjm1_wshift
Effective source at the previous step for the resolution of the vector Poisson equation for by means...
Definition: star_bhns.h:220
Scalar pot_centri
Centrifugal potential.
Definition: star_bhns.h:110
double ray_pole() const
Coordinate radius at [r_unit].
Definition: star_global.C:278
Scalar lapconf_auto
Lapconf function generated by the star.
Definition: star_bhns.h:113
bool irrotational
true for an irrotational star, false for a corotating one
Definition: star_bhns.h:72
virtual void del_deriv() const
Deletes all the derived quantities.
Definition: star.C:298
Scalar & set_lapconf_comp()
Read/write of the lapconf function generated by the companion black hole.
Definition: star_bhns.C:428
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
Vector d_confo_auto
Derivative of the conformal factor generated by the star .
Definition: star_bhns.h:168
Scalar confo_auto
Conformal factor generated by the star.
Definition: star_bhns.h:157
Scalar confo_comp
Conformal factor generated by the companion black hole.
Definition: star_bhns.h:160
double ray_eq() const
Coordinate radius at , [r_unit].
Definition: star_global.C:108
Sym_tensor taij_auto
Part of the extrinsic curvature tensor generated by shift_auto , lapse_auto , and confo_auto ...
Definition: star_bhns.h:182
Scalar ssjm1_lapconf
Effective source at the previous step for the resolution of the Poisson equation for lapconf_auto ...
Definition: star_bhns.h:197
Scalar & set_confo_comp()
Read/write of the conformal factor generated by the companion black hole.
Definition: star_bhns.C:456
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
void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition: star_bhns.C:352
Scalar & set_lapconf_auto()
Read/write of the lapconf function generated by the neutron star.
Definition: star_bhns.C:421
Scalar & set(int)
Read/write access to a component.
Definition: vector.C:296
void equation_of_state()
Computes the proper baryon and energy density, as well as pressure from the enthalpy.
Definition: star.C:462
Scalar ssjm1_confo
Effective source at the previous step for the resolution of the Poisson equation for confo_auto ...
Definition: star_bhns.h:202
virtual void sauve(FILE *) const
Save in a file.
Definition: star_bhns.C:470
virtual void annule(int l_min, int l_max)
Sets the Tensor to zero in several domains.
Definition: tensor.C:671
Scalar lapse_auto
Lapse function generated by the "star".
Definition: star_bhns.h:122
Scalar psi0
Scalar potential of the non-translational part of fluid 4-velocity (in the irrotational case) ...
Definition: star_bhns.h:77
Valeur va
The numerical value of the Cmp.
Definition: cmp.h:461
Tensor d_shift_comp
Derivative of the shift vector generated by the companion black hole.
Definition: star_bhns.h:154