LORENE
star_bhns.h
1/*
2 * Definition of Lorene class Star_bhns
3 *
4 */
5
6/*
7 * Copyright (c) 2005-2007 Keisuke Taniguchi
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#ifndef __STAR_BHNS_H_
27#define __STAR_BHNS_H_
28
29/*
30 * $Id: star_bhns.h,v 1.3 2014/10/13 08:52:36 j_novak Exp $
31 * $Log: star_bhns.h,v $
32 * Revision 1.3 2014/10/13 08:52:36 j_novak
33 * Lorene classes and functions now belong to the namespace Lorene.
34 *
35 * Revision 1.2 2008/05/15 18:55:55 k_taniguchi
36 * Change of some parameters.
37 *
38 * Revision 1.1 2007/06/22 01:04:35 k_taniguchi
39 * *** empty log message ***
40 *
41 *
42 * $Header: /cvsroot/Lorene/C++/Include/star_bhns.h,v 1.3 2014/10/13 08:52:36 j_novak Exp $
43 *
44 */
45
46// Lorene headers
47#include "star.h"
48
49namespace Lorene {
50
51// External classes which appear in the declaration of class Star_bhns:
52class Hole_bhns ;
53
59class Star_bhns : public Star {
60
61 // Data :
62 // -----
63 protected:
68
73
78
83
89
94
100
103
108
111
114
117
120
122 Scalar lapse_auto ; // = lapconf_auto / confo_tot
123
126
131
136
139
142
145
150
155
158
161
164
169
174
176 Scalar psi4 ; // psi4 = pow(confo_tot, 4.)
177
183
188
193
198
203
211
221
222 // Derived data
223 // ------------
224 protected:
225 mutable double* p_mass_b_bhns ;
226 mutable double* p_mass_g_bhns ;
227
228 // Constructors - Destructor
229 // -------------------------
230 public:
239 Star_bhns(Map& mp_i, int nzet_i, const Eos& eos_i, bool irrot_i) ;
240
241 Star_bhns(const Star_bhns& ) ;
242
249 Star_bhns(Map& mp_i, const Eos& eos_i, FILE* fich) ;
250
251 virtual ~Star_bhns() ;
252
253
254 // Memory management
255 // -----------------
256 protected:
258 virtual void del_deriv() const ;
259
261 void set_der_0x0() const ;
262
263 // Mutators / assignment
264 // ---------------------
265 public:
267 void operator=(const Star_bhns&) ;
268
271
274
279
282
287
290
295
296 // Accessors
297 // ---------
298 public:
302 bool is_irrotational() const {return irrotational; } ;
303
305 const Scalar& get_psi0() const {return psi0; } ;
306
310 const Vector& get_d_psi() const {return d_psi; } ;
311
316 const Vector& get_wit_w() const {return wit_w; } ;
317
321 const Scalar& get_loggam() const {return loggam; } ;
322
327 const Vector& get_bsn() const {return bsn; } ;
328
330 const Scalar& get_gam() const {return gam; } ;
331
333 const Scalar& get_gam0() const {return gam0; } ;
334
336 const Scalar& get_pot_centri() const {return pot_centri; } ;
337
339 const Scalar& get_lapconf_auto() const {return lapconf_auto; } ;
340
344 const Scalar& get_lapconf_comp() const {return lapconf_comp; } ;
345
347 const Scalar& get_lapconf_tot() const {return lapconf_tot; } ;
348
349 // Returns the part of the lapse function generated by the star
350 const Scalar& get_lapse_auto() const {return lapse_auto; } ;
351
353 const Scalar& get_lapse_tot() const {return lapse_tot; } ;
354
356 const Vector& get_d_lapconf_auto() const {return d_lapconf_auto; } ;
357
361 const Vector& get_d_lapconf_comp() const {return d_lapconf_comp; } ;
362
364 const Vector& get_shift_auto() const {return shift_auto; } ;
365
369 const Vector& get_shift_comp() const {return shift_comp; } ;
370
372 const Vector& get_shift_tot() const {return shift_tot; } ;
373
375 const Tensor& get_d_shift_auto() const {return d_shift_auto; } ;
376
380 const Tensor& get_d_shift_comp() const {return d_shift_comp; } ;
381
383 const Scalar& get_confo_auto() const {return confo_auto; } ;
384
388 const Scalar& get_confo_comp() const {return confo_comp; } ;
389
391 const Scalar& get_confo_tot() const {return confo_tot; } ;
392
396 const Vector& get_d_confo_auto() const {return d_confo_auto; } ;
397
401 const Vector& get_d_confo_comp() const {return d_confo_comp; } ;
402
404 const Scalar& get_psi4() const {return psi4; } ;
405
409 const Sym_tensor& get_taij_auto() const {return taij_auto; } ;
410
415 const Scalar& get_taij_quad_auto() const {return taij_quad_auto; } ;
416
417 // Outputs
418 // -------
419 public:
420 virtual void sauve(FILE *) const ;
421
422 protected:
424 virtual ostream& operator>>(ostream& ) const ;
425
426 // Global quantities
427 // -----------------
428 public:
430 virtual double mass_b() const ;
431
432 virtual double mass_b_bhns(bool kerrschild, const double& mass_bh,
433 const double& sepa) const ;
434
436 virtual double mass_g() const ;
437
438 virtual double mass_g_bhns() const ;
439
440 // Computational routines
441 // ----------------------
442 public:
461 void hydro_euler_bhns(bool kerrschild, const double& mass_bh,
462 const double& sepa) ;
463
480 void update_metric_bhns(const Hole_bhns& hole,
481 const Star_bhns& star_prev,
482 double relax) ;
483
490 void update_met_der_comp_bhns(const Hole_bhns& hole) ;
491
508 void kinema_bhns(bool kerrschild, const double& mass_bh,
509 const double& sepa, double omega,
510 double x_rot, double y_rot) ;
511
513 void fait_d_psi_bhns() ;
514
518 void extr_curv_bhns() ;
519
548 void equilibrium_bhns(double ent_c, const double& mass_bh,
549 const double& sepa, bool kerrschild,
550 int mer, int mermax_ns, int mermax_potvit,
551 int mermax_poisson, int filter_r, int filter_r_s,
552 int filter_p_s, double relax_poisson,
553 double relax_potvit, double thres_adapt,
554 double resize_ns,
555 const Tbl& fact_resize, Tbl& diff) ;
556
574 double velo_pot_bhns(const double& mass_bh, const double& sepa,
575 bool kerrschild,
576 int mermax, double precis, double relax) ;
577
584 double chi_rp(double radius, double phi) ;
585
591 double radius_p(double phi) ;
592
596 double phi_min() ;
597
603 double phi_local_min(double phi_ini) ;
604
616 void relax_bhns(const Star_bhns& star_prev, double relax_ent,
617 double relax_met, int mer, int fmer_met) ;
618
624 void equil_spher_bhns(double ent_c, double precis) ;
625
626 friend class Bin_bhns ;
627
628};
629
630}
631#endif
Class for computing a black hole - neutron star binary system with comparable mass ()
Definition bin_bhns.h:63
Equation of state base class.
Definition eos.h:190
Class for black holes in black hole-neutron star binaries.
Definition hole_bhns.h:65
Affine radial mapping.
Definition map.h:2027
Base class for coordinate mappings.
Definition map.h:670
Flat metric for tensor calculation.
Definition metric.h:261
Tensor field of valence 0 (or component of a tensorial field).
Definition scalar.h:387
Class for stars in black hole-neutron star binaries.
Definition star_bhns.h:59
Metric_flat flat
Flat metric defined on the mapping (Spherical components with respect to the mapping of the star ).
Definition star_bhns.h:192
Scalar & set_confo_comp()
Read/write of the conformal factor generated by the companion black hole.
Definition star_bhns.C:456
const Vector & get_d_psi() const
Returns the covariant derivative of the velocity potential (Spherical components with respect to the ...
Definition star_bhns.h:310
virtual void sauve(FILE *) const
Save in a file.
Definition star_bhns.C:470
Vector shift_tot
Total shift vector.
Definition star_bhns.h:144
void equil_spher_bhns(double ent_c, double precis)
Computes a spherical configuration.
Sym_tensor taij_auto
Part of the extrinsic curvature tensor generated by shift_auto , lapse_auto , and confo_auto .
Definition star_bhns.h:182
const Scalar & get_lapconf_tot() const
Returns the total lapconf function.
Definition star_bhns.h:347
Scalar lapconf_auto
Lapconf function generated by the star.
Definition star_bhns.h:113
Scalar confo_comp
Conformal factor generated by the companion black hole.
Definition star_bhns.h:160
Vector wit_w
Spatial projection of the fluid 3-velocity with respect to the co-orbiting observer.
Definition star_bhns.h:88
const Scalar & get_lapse_tot() const
Returns the total lapse function.
Definition star_bhns.h:353
Tensor d_shift_auto
Derivative of the shift vector generated by the star .
Definition star_bhns.h:149
double velo_pot_bhns(const double &mass_bh, const double &sepa, bool kerrschild, int mermax, double precis, double relax)
Computes the non-translational part of the velocity scalar potential by solving the continuity equat...
Vector bsn
3-vector shift, divided by N , of the rotating coordinates, .
Definition star_bhns.h:99
Scalar & set_confo_auto()
Read/write of the conformal factor generated by the neutron star.
Definition star_bhns.C:449
const Tensor & get_d_shift_comp() const
Returns the derivative of the shift vector generated by the companion black hole.
Definition star_bhns.h:380
Vector & set_shift_auto()
Read/write of the shift vector generated by the neutron star.
Definition star_bhns.C:435
const Scalar & get_psi0() const
Returns the non-translational part of the velocity potential.
Definition star_bhns.h:305
Scalar & set_lapconf_auto()
Read/write of the lapconf function generated by the neutron star.
Definition star_bhns.C:421
const Scalar & get_pot_centri() const
Returns the centrifugal potential.
Definition star_bhns.h:336
const Scalar & get_gam() const
Returns the Lorentz factor gam.
Definition star_bhns.h:330
Scalar & set_pot_centri()
Read/write the centrifugal potential.
Definition star_bhns.C:414
const Vector & get_wit_w() const
Returns the spatial projection of the fluid 3-velocity with respect to the co-orbiting observer.
Definition star_bhns.h:316
void equilibrium_bhns(double ent_c, const double &mass_bh, const double &sepa, bool kerrschild, int mer, int mermax_ns, int mermax_potvit, int mermax_poisson, int filter_r, int filter_r_s, int filter_p_s, double relax_poisson, double relax_potvit, double thres_adapt, double resize_ns, const Tbl &fact_resize, Tbl &diff)
Computes an equilibrium configuration.
void fait_d_psi_bhns()
Computes the gradient of the total velocity potential .
Definition star_bhns.C:535
const Scalar & get_loggam() const
Returns the logarithm of the Lorentz factor between the fluid and the co-orbiting observer.
Definition star_bhns.h:321
Scalar lapse_tot
Total lapse function.
Definition star_bhns.h:125
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
virtual void del_deriv() const
Deletes all the derived quantities.
Definition star_bhns.C:341
Vector d_confo_comp
Derivative of the conformal factor generated by the companion black hole.
Definition star_bhns.h:173
void operator=(const Star_bhns &)
Assignment to another Star_bhns.
Definition star_bhns.C:368
void kinema_bhns(bool kerrschild, const double &mass_bh, const double &sepa, double omega, double x_rot, double y_rot)
Computes the quantities bsn and pot_centri .
Scalar confo_tot
Total conformal factor.
Definition star_bhns.h:163
Scalar psi4
Fourth power of the total conformal factor.
Definition star_bhns.h:176
bool is_irrotational() const
Returns true for an irrotational motion, false for a corotating one.
Definition star_bhns.h:302
Vector d_lapconf_auto
Derivative of the lapconf function generated by the star .
Definition star_bhns.h:130
Vector d_psi
Gradient of (in the irrotational case) (Spherical components with respect to the mapping of the star...
Definition star_bhns.h:82
Vector d_confo_auto
Derivative of the conformal factor generated by the star .
Definition star_bhns.h:168
double radius_p(double phi)
Radius of the star to the direction of and .
const Vector & get_d_lapconf_auto() const
Returns the derivative of the lapse function generated by the star.
Definition star_bhns.h:356
Scalar & set_lapconf_comp()
Read/write of the lapconf function generated by the companion black hole.
Definition star_bhns.C:428
Scalar gam0
Lorentz factor between the co-orbiting observer and the Eulerian one.
Definition star_bhns.h:107
const Scalar & get_confo_tot() const
Returns the total conformal factor.
Definition star_bhns.h:391
const Scalar & get_gam0() const
Returns the Lorentz factor gam0.
Definition star_bhns.h:333
void hydro_euler_bhns(bool kerrschild, const double &mass_bh, const double &sepa)
Computes the hydrodynamical quantities relative to the Eulerian observer from those in the fluid fram...
const Vector & get_bsn() const
Returns the shift vector, divided by N , of the rotating coordinates, .
Definition star_bhns.h:327
Vector shift_comp
Shift vector generated by the companion black hole.
Definition star_bhns.h:141
const Vector & get_shift_tot() const
Returns the part total shift vector.
Definition star_bhns.h:372
const Scalar & get_confo_auto() const
Returns the part of the conformal factor generated by the star.
Definition star_bhns.h:383
Scalar lapconf_tot
Total lapconf function.
Definition star_bhns.h:119
const Scalar & get_confo_comp() const
Returns the part of the conformal factor generated by the companion black hole.
Definition star_bhns.h:388
const Sym_tensor & get_taij_auto() const
Returns the part of the extrinsic curvature tensor generated by the neutron star part.
Definition star_bhns.h:409
Vector d_lapconf_comp
Derivative of the lapconf function generated by the companion black hole.
Definition star_bhns.h:135
const Vector & get_d_confo_comp() const
Returns the derivative of the conformal factor generated by the companion black hole.
Definition star_bhns.h:401
Vector & set_shift_comp()
Read/write of the shift vector generated by the companion black hole.
Definition star_bhns.C:442
Scalar gam
Lorentz factor between the fluid and the co-orbiting observer.
Definition star_bhns.h:102
Scalar ssjm1_lapconf
Effective source at the previous step for the resolution of the Poisson equation for lapconf_auto .
Definition star_bhns.h:197
bool irrotational
true for an irrotational star, false for a corotating one
Definition star_bhns.h:72
Scalar pot_centri
Centrifugal potential.
Definition star_bhns.h:110
const Vector & get_shift_auto() const
Returns the part of the shift vector generated by the star.
Definition star_bhns.h:364
double chi_rp(double radius, double phi)
Sensitive indicator of the mass-shedding to the direction of , , .
const Scalar & get_lapconf_comp() const
Returns the part of the lapconf function generated by the companion black hole.
Definition star_bhns.h:344
void extr_curv_bhns()
Computes taij_auto , taij_quad_auto from shift_auto , lapse_auto , confo_auto .
void set_der_0x0() const
Sets to 0x0 all the pointers on derived quantities.
Definition star_bhns.C:352
double * p_mass_b_bhns
Baryon mass.
Definition star_bhns.h:225
Scalar lapse_auto
Lapse function generated by the "star".
Definition star_bhns.h:122
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
Scalar psi0
Scalar potential of the non-translational part of fluid 4-velocity (in the irrotational case)
Definition star_bhns.h:77
void update_met_der_comp_bhns(const Hole_bhns &hole)
Computes derivative of metric quantities from the companion black hole.
virtual double mass_g() const
Gravitational mass.
virtual double mass_b() const
Baryon mass.
Scalar lapconf_comp
Lapconf function generated by the companion black hole.
Definition star_bhns.h:116
const Vector & get_shift_comp() const
Returns the part of the shift vector generated by the companion black hole.
Definition star_bhns.h:369
Scalar taij_quad_auto
Part of the scalar generated by .
Definition star_bhns.h:187
double phi_local_min(double phi_ini)
Azimuthal angle when the indicator of the mass-shedding takes its local minimum.
const Tensor & get_d_shift_auto() const
Returns the derivative of the shift vector generated by the star.
Definition star_bhns.h:375
virtual ostream & operator>>(ostream &) const
Operator >> (virtual function called by the operator <<).
Definition star_bhns.C:495
double phi_min()
Azimuthal angle when the indicator of the mass-shedding takes its minimum chi_min.
void update_metric_bhns(const Hole_bhns &hole, const Star_bhns &star_prev, double relax)
Computes metric coefficients from known potentials with relaxation when the companion is a black hole...
const Scalar & get_lapconf_auto() const
Returns the part of the lapconf function generated by the star.
Definition star_bhns.h:339
Vector shift_auto
Shift vector generated by the star.
Definition star_bhns.h:138
const Scalar & get_psi4() const
Returns the fourth power of the total conformal factor.
Definition star_bhns.h:404
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
const Scalar & get_taij_quad_auto() const
Returns the part of the scalar generated by .
Definition star_bhns.h:415
Scalar confo_auto
Conformal factor generated by the star.
Definition star_bhns.h:157
Scalar ssjm1_confo
Effective source at the previous step for the resolution of the Poisson equation for confo_auto .
Definition star_bhns.h:202
const Vector & get_d_lapconf_comp() const
Returns the derivative of the lapse function generated by the companion black hole.
Definition star_bhns.h:361
Map_af mp_aff
Affine mapping for solving poisson's equations of metric quantities.
Definition star_bhns.h:67
Tensor d_shift_comp
Derivative of the shift vector generated by the companion black hole.
Definition star_bhns.h:154
Scalar loggam
Logarithm of the Lorentz factor between the fluid and the co-orbiting observer.
Definition star_bhns.h:93
const Vector & get_d_confo_auto() const
Returns the derivative of the conformal factor generated by the star.
Definition star_bhns.h:396
virtual ~Star_bhns()
Destructor.
Definition star_bhns.C:329
double * p_mass_g_bhns
Gravitational mass.
Definition star_bhns.h:226
Base class for stars.
Definition star.h:175
Class intended to describe valence-2 symmetric tensors.
Definition sym_tensor.h:223
Basic array class.
Definition tbl.h:161
Tensor handling.
Definition tensor.h:288
Tensor field of valence 1.
Definition vector.h:188
Lorene prototypes.
Definition app_hor.h:64