LORENE
et_bin_nsbh_equilibrium.C
1 /*
2  * Method of class Etoile to compute a static spherical configuration
3  * of a neutron star in a NS-BH binary system.
4  *
5  * (see file etoile.h for documentation).
6  *
7  */
8 
9 /*
10  * Copyright (c) 2003 Keisuke Taniguchi
11  *
12  * This file is part of LORENE.
13  *
14  * LORENE is free software; you can redistribute it and/or modify
15  * it under the terms of the GNU General Public License version 2
16  * as published by the Free Software Foundation.
17  *
18  * LORENE is distributed in the hope that it will be useful,
19  * but WITHOUT ANY WARRANTY; without even the implied warranty of
20  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
21  * GNU General Public License for more details.
22  *
23  * You should have received a copy of the GNU General Public License
24  * along with LORENE; if not, write to the Free Software
25  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
26  *
27  */
28 
29 
30 
31 /*
32  * $Id: et_bin_nsbh_equilibrium.C,v 1.14 2016/12/05 16:17:53 j_novak Exp $
33  * $Log: et_bin_nsbh_equilibrium.C,v $
34  * Revision 1.14 2016/12/05 16:17:53 j_novak
35  * Suppression of some global variables (file names, loch, ...) to prevent redefinitions
36  *
37  * Revision 1.13 2014/10/13 08:52:56 j_novak
38  * Lorene classes and functions now belong to the namespace Lorene.
39  *
40  * Revision 1.12 2014/10/06 15:13:08 j_novak
41  * Modified #include directives to use c++ syntax.
42  *
43  * Revision 1.11 2008/09/26 08:38:45 p_grandclement
44  * get rid of desaliasing
45  *
46  * Revision 1.10 2006/09/05 13:39:45 p_grandclement
47  * update of the bin_ns_bh project
48  *
49  * Revision 1.9 2006/06/01 12:47:53 p_grandclement
50  * update of the Bin_ns_bh project
51  *
52  * Revision 1.8 2006/04/25 07:21:58 p_grandclement
53  * Various changes for the NS_BH project
54  *
55  * Revision 1.7 2006/03/30 07:33:47 p_grandclement
56  * *** empty log message ***
57  *
58  * Revision 1.6 2005/10/18 13:12:33 p_grandclement
59  * update of the mixted binary codes
60  *
61  * Revision 1.5 2005/08/29 15:10:17 p_grandclement
62  * Addition of things needed :
63  * 1) For BBH with different masses
64  * 2) Provisory files for the mixted binaries (Bh and NS) : THIS IS NOT
65  * WORKING YET !!!
66  *
67  * Revision 1.4 2004/03/25 10:29:04 j_novak
68  * All LORENE's units are now defined in the namespace Unites (in file unites.h).
69  *
70  * Revision 1.3 2003/10/24 12:34:06 k_taniguchi
71  * Change the notation as it should be
72  *
73  * Revision 1.2 2003/10/21 11:49:33 k_taniguchi
74  * Change the class from Etoile_bin to sub-class Et_bin_nsbh.
75  *
76  * Revision 1.1 2003/10/20 15:01:55 k_taniguchi
77  * Computation of an equilibrium configuration of a neutron star
78  * in a NS-BH binary system.
79  *
80  *
81  * $Header: /cvsroot/Lorene/C++/Source/Etoile/et_bin_nsbh_equilibrium.C,v 1.14 2016/12/05 16:17:53 j_novak Exp $
82  *
83  */
84 
85 // Headers C
86 #include <cmath>
87 
88 // Headers Lorene
89 #include "etoile.h"
90 #include "map.h"
91 #include "nbr_spx.h"
92 #include "et_bin_nsbh.h"
93 #include "param.h"
94 
95 #include "graphique.h"
96 #include "utilitaires.h"
97 #include "unites.h"
98 
99 namespace Lorene {
100 void Et_bin_nsbh::equilibrium_nsbh(bool adapt, double ent_c, int& niter, int mermax,
101  int mermax_poisson, double relax_poisson,
102  int mermax_potvit, double relax_potvit,
103  Tbl& diff) {
104 
105  // Fundamental constants and units
106  // -------------------------------
107  using namespace Unites ;
108 
109  // Initializations
110  // --------------
111 
112  const Mg3d* mg = mp.get_mg() ;
113  int nz = mg->get_nzone() ; // total number of domains
114 
115  // The following is required to initialize mp_prev as a Map_et:
116  Map_et& mp_et = dynamic_cast<Map_et&>(mp) ;
117 
118  // Error indicators
119  // ----------------
120  double& diff_ent = diff.set(0) ;
121  double& diff_vel_pot = diff.set(1) ;
122  double& diff_lapse = diff.set(2) ;
123  double& diff_confpsi = diff.set(3) ;
124  double& diff_shift_x = diff.set(4) ;
125  double& diff_shift_y = diff.set(5) ;
126  double& diff_shift_z = diff.set(6) ;
127 
128 
129  // Parameters for the function Map_et::poisson for n_auto
130  // ------------------------------------------------------
131  double precis_poisson = 1.e-16 ;
132 
133  Param par_poisson1 ;
134  par_poisson1.add_int(mermax_poisson, 0) ; // maximum number of iterations
135  par_poisson1.add_double(relax_poisson, 0) ; // relaxation parameter
136  par_poisson1.add_double(precis_poisson, 1) ; // required precision
137  par_poisson1.add_int_mod(niter, 0) ; // number of iterations actually used
138  par_poisson1.add_cmp_mod( ssjm1_lapse ) ;
139 
140  // Parameters for the function Map_et::poisson for confpsi_auto
141  // ------------------------------------------------------------
142 
143  Param par_poisson2 ;
144  par_poisson2.add_int(mermax_poisson, 0) ; // maximum number of iterations
145  par_poisson2.add_double(relax_poisson, 0) ; // relaxation parameter
146  par_poisson2.add_double(precis_poisson, 1) ; // required precision
147  par_poisson2.add_int_mod(niter, 0) ; // number of iterations actually used
148  par_poisson2.add_cmp_mod( ssjm1_confpsi ) ;
149 
150  // Parameters for the function Tenseur::poisson_vect
151  // -------------------------------------------------
152 
153  Param par_poisson_vect ;
154  par_poisson_vect.add_int(mermax_poisson, 0) ;
155  // maximum number of iterations
156  par_poisson_vect.add_double(relax_poisson, 0) ; // relaxation parameter
157  par_poisson_vect.add_double(precis_poisson, 1) ; // required precision
158  par_poisson_vect.add_cmp_mod( ssjm1_khi ) ;
159  par_poisson_vect.add_tenseur_mod( ssjm1_wshift ) ;
160  par_poisson_vect.add_int_mod(niter, 0) ;
161 
162  // Parameters for the adaptation
163  Param par_adapt ;
164  int nitermax = 100 ;
165  int niter_adapt ;
166  int adapt_flag = (adapt) ? 1 : 0 ;
167  int nz_search = nzet + 1 ;
168  double precis_secant = 1.e-14 ;
169  double alpha_r ;
170  double reg_map = 1. ;
171  int k_b ;
172  int j_b ;
173  Tbl ent_limit(nzet) ;
174 
175  par_adapt.add_int(nitermax, 0) ;
176  par_adapt.add_int(nzet, 1) ;
177  par_adapt.add_int(nz_search, 2) ;
178  par_adapt.add_int(adapt_flag, 3) ;
179  par_adapt.add_int(j_b, 4) ;
180  par_adapt.add_int(k_b, 5) ;
181  par_adapt.add_int_mod(niter_adapt, 0) ;
182  par_adapt.add_double(precis_secant, 0) ;
183  par_adapt.add_double(reg_map, 1) ;
184  par_adapt.add_double(alpha_r, 2) ;
185  par_adapt.add_tbl(ent_limit, 0) ;
186 
187  // External potential
188  // See Eq (99) from Gourgoulhon et al. (2001)
189  // -----------------------------------------
190 
191 
192  Tenseur ent_jm1 = ent ; // Enthalpy at previous step
193  Tenseur source(mp) ; // source term in the equation for logn_auto
194  // and beta_auto
195  Tenseur source_shift(mp, 1, CON, ref_triad) ; // source term in the
196  // equation for shift_auto
197 
198  //=========================================================================
199  // Start of iteration
200  //=========================================================================
201 
202  for(int mer=0 ; mer<mermax ; mer++ ) {
203 
204  cout << "-----------------------------------------------" << endl ;
205  cout << "step: " << mer << endl ;
206  cout << "diff_ent = " << diff_ent << endl ;
207  //-----------------------------------------------------
208  // Resolution of the elliptic equation for the velocity
209  // scalar potential
210  //-----------------------------------------------------
211 
212  if (irrotational) {
213  diff_vel_pot = velocity_potential(mermax_potvit, precis_poisson,
214  relax_potvit) ;
215  }
216 
217  // Equation de la surface
218  //if (adapt) {
219 
220  // Rescaling of the radius : (Be carefull !)
221  int nt = mg->get_nt(nzet-1) ;
222  int np = mg->get_np(nzet-1) ;
223  int nr = mg->get_nr(nzet-1) ;
224 
225  // valeurs au centre
226  double hc = exp(ent_c) ;
227  double gamma_c = exp(loggam())(0,0,0,0) ;
228  double gamma_0_c = exp(-pot_centri())(0,0,0,0) ;
229  double n_auto_c = n_auto()(0,0,0,0) ;
230  double n_comp_c = n_comp()(0,0,0,0) ;
231 
232  double alpha_square = 0 ;
233  double constante = 0;
234  for (int k=0; k<np; k++) {
235  for (int j=0; j<nt; j++) {
236 
237  // valeurs au bord
238  double gamma_b = exp(loggam())(nzet-1,k,j,nr-1) ;
239  double gamma_0_b = exp(-pot_centri())(nzet-1,k,j,nr-1) ;
240  double n_auto_b = n_auto()(nzet-1,k,j,nr-1) ;
241  double n_comp_b = n_comp()(nzet-1,k,j,nr-1) ;
242 
243  // Les solutions :
244  double alpha_square_courant = (gamma_0_c*gamma_b*n_comp_b - hc*gamma_c*gamma_0_b*n_comp_c) /
245  (hc*gamma_c*gamma_0_b*n_auto_c-gamma_0_c*gamma_b*n_auto_b) ;
246  double constante_courant = gamma_b*(n_comp_b+alpha_square_courant*n_auto_b)/gamma_0_b ;
247 
248  if (alpha_square_courant > alpha_square) {
249  alpha_square = alpha_square_courant ;
250  k_b = k ;
251  j_b = j ;
252  constante = constante_courant ;
253  }
254  }
255  }
256 
257  alpha_r = sqrt(alpha_square) ;
258  cout << "Adaptation : " << k_b << " " << j_b << " " << alpha_r << endl ;
259 
260  // Le potentiel :
261  Tenseur potentiel (constante*exp(-loggam-pot_centri)/(n_auto*alpha_square+n_comp)) ;
262  potentiel.set_std_base() ;
263  for (int l=nzet+1 ; l<nz ; l++)
264  potentiel.set().va.set(l) = 1 ;
265 
266  Map_et mp_prev = mp_et ;
267  ent = log(potentiel) ;
268  ent.set_std_base() ;
269  ent().va.smooth(nzet, (ent.set().va)) ;
270 
271  ent_limit.set_etat_qcq() ;
272  for (int l=0; l<nzet; l++) { // loop on domains inside the star
273  ent_limit.set(l) = ent()(l, k_b, j_b, nr-1) ;
274  }
275 
276  // On adapte :
277  mp.adapt(ent(), par_adapt, 4) ;
278  mp_prev.homothetie(alpha_r) ;
279 
280  for (int l=nzet ; l<nz-1 ; l++)
281  mp.resize(l, 1./alpha_r) ;
282  mp.reevaluate_symy (&mp_prev, nzet, ent.set()) ;
283 
284 
285 
286  // Equation of state
287  //----------------------------------------------------
288  equation_of_state() ; // computes new values for nbar (n), ener (e)
289  // and press (p) from the new ent (H)
290 
291  //---------------------------------------------------------
292  // Matter source terms in the gravitational field equations
293  //---------------------------------------------------------
294  hydro_euler() ; // computes new values for ener_euler (E),
295  // s_euler (S) and u_euler (U^i)
296 
297 
298  //-------------------------------------------------
299  // Relative change in enthalpy
300  //-------------------------------------------------
301 
302  Tbl diff_ent_tbl = diffrel( ent(), ent_jm1() ) ;
303  diff_ent = diff_ent_tbl(0) ;
304  for (int l=1; l<nzet; l++) {
305  diff_ent += diff_ent_tbl(l) ;
306  }
307  diff_ent /= nzet ;
308 
309  ent_jm1 = ent ;
310 
311 
312  //--------------------------------------------------------
313  // Poisson equation for n_auto
314  //--------------------------------------------------------
315 
316  // Source
317  // See Eq (50) from Gourgoulhon et al. (2001)
318  // ------------------------------------------
319 
320  Tenseur confpsi_q = pow(confpsi, 4.) ;
321  Tenseur confpsi_c = pow(confpsi, 5.) ;
322 
323  if (relativistic) {
324  Tenseur tmp = flat_scalar_prod(tkij_tot, tkij_auto) ;
325  Tenseur kk (mp) ;
326  kk = 0 ;
327  Tenseur tmp2(mp) ;
328  tmp2.set_etat_qcq() ;
329  for (int i=0 ; i<3 ; i++) {
330  tmp2.set() = tmp(i, i) ;
331  kk = kk + tmp2 ;
332  }
333 
334  source = qpig * nnn * confpsi_q * (ener_euler + s_euler)
335  + nnn * confpsi_q * kk
337  confpsi ;
338  }
339  else {
340  cout <<
341  "WARNING : Et_bin_nsbh is for the relativistic calculation"
342  << endl ;
343  abort() ;
344  }
345 
346  source.set_std_base() ;
347 
348  // Resolution of the Poisson equation
349  // ----------------------------------
350  Cmp n_auto_old (n_auto()) ;
351  source().poisson(par_poisson1, n_auto.set()) ;
352  n_auto.set() = n_auto() + 0.5 ;
353 
354  // Difference pas précédent
355  // -----------------------------------------------------
356 
357  Tbl tdiff_lapse = diffrel(n_auto(), n_auto_old) ;
358  cout <<
359  "Relative difference on n_auto : "
360  << endl ;
361  for (int l=0; l<nz; l++) {
362  cout << tdiff_lapse(l) << " " ;
363  }
364  cout << endl ;
365  diff_lapse = max(abs(tdiff_lapse)) ;
366 
367  if (relativistic) {
368 
369 
370  //--------------------------------------------------------
371  // Poisson equation for confpsi_auto
372  //--------------------------------------------------------
373 
374  // Source
375  // See Eq (51) from Gourgoulhon et al. (2001)
376  // ------------------------------------------
377 
378  Tenseur tmp = flat_scalar_prod(tkij_tot, tkij_auto) ;
379  Tenseur kk (mp) ;
380  kk = 0 ;
381  Tenseur tmp2(mp) ;
382  tmp2.set_etat_qcq() ;
383  for (int i=0 ; i<3 ; i++) {
384  tmp2.set() = tmp(i, i) ;
385  kk = kk + tmp2 ;
386  }
387 
388  source = -0.5 * qpig * confpsi_c * ener_euler
389  - 0.125 * confpsi_c * kk ;
390 
391  source.set_std_base() ;
392 
393  // Resolution of the Poisson equation
394  // ----------------------------------
395  Cmp psi_old (confpsi_auto()) ;
396  source().poisson(par_poisson2, confpsi_auto.set()) ;
397  confpsi_auto.set() = confpsi_auto() + 0.5 ;
398 
399 
400  // Check: has the Poisson equation been correctly solved ?
401  // -----------------------------------------------------
402 
403  Tbl tdiff_confpsi = diffrel(confpsi_auto(), psi_old) ;
404  cout <<
405  "Relative difference on confpsi_auto : "
406  << endl ;
407  for (int l=0; l<nz; l++) {
408  cout << tdiff_confpsi(l) << " " ;
409  }
410  cout << endl ;
411  diff_confpsi = max(abs(tdiff_confpsi)) ;
412 
413  //--------------------------------------------------------
414  // Vector Poisson equation for shift_auto
415  //--------------------------------------------------------
416 
417  // Source
418  // See Eq (52) from Gourgoulhon et al. (2001)
419  // ------
420  Tenseur vtmp = d_n_auto -6. * nnn * d_confpsi_auto / confpsi ;
421  source_shift = 4.*qpig * nnn *confpsi_q *(ener_euler + press)
422  * u_euler ;
423  if (tkij_tot.get_etat() != ETATZERO)
424  source_shift = source_shift + 2.* flat_scalar_prod(tkij_tot, vtmp) ;
425  source_shift.set_std_base() ;
426  // Resolution of the Poisson equation
427  // ----------------------------------
428  // Filter for the source of shift vector
429  for (int i=0 ; i<3 ; i++)
430  if (source_shift(i).get_etat() != ETATZERO)
431  source_shift.set(i).va.coef_i() ;
432 
433 for (int i=0; i<3; i++)
434  if ((source_shift(i).get_etat() != ETATZERO) && (source_shift(i).va.c->t[nz-1]->get_etat() != ETATZERO))
435  source_shift.set(i).filtre(4) ;
436  for (int i=0; i<3; i++) {
437  if(source_shift(i).dz_nonzero()) {
438  assert( source_shift(i).get_dzpuis() == 4 ) ;
439  }
440  else{
441  (source_shift.set(i)).set_dzpuis(4) ;
442  }
443  }
444  //##
445  // source_shift.dec2_dzpuis() ; // dzpuis 4 -> 2
446 
447  double lambda_shift = double(1) / double(3) ;
448  // ON DOIT CHANGER DE TRIADE
449  source_shift.change_triad(mp.get_bvect_cart()) ;
450  Tenseur shift_old (shift_auto) ;
451  source_shift.poisson_vect_oohara(lambda_shift, par_poisson_vect,
454 
455  // Check: has the equation for shift_auto been correctly solved ?
456  // --------------------------------------------------------------
457 
458 
459 
460  Tbl tdiff_shift_x = diffrel(shift_auto(0), shift_old(0)) ;
461  Tbl tdiff_shift_y = diffrel(shift_auto(1), shift_old(1)) ;
462  Tbl tdiff_shift_z = diffrel(shift_auto(2), shift_old(2)) ;
463 
464  cout <<
465  "Relative difference on shift_auto : "
466  << endl ;
467  cout << "x component : " ;
468  for (int l=0; l<nz; l++) {
469  cout << tdiff_shift_x(l) << " " ;
470  }
471  cout << endl ;
472  cout << "y component : " ;
473  for (int l=0; l<nz; l++) {
474  cout << tdiff_shift_y(l) << " " ;
475  }
476  cout << endl ;
477  cout << "z component : " ;
478  for (int l=0; l<nz; l++) {
479  cout << tdiff_shift_z(l) << " " ;
480  }
481  cout << endl ;
482 
483  diff_shift_x = max(abs(tdiff_shift_x)) ;
484  diff_shift_y = max(abs(tdiff_shift_y)) ;
485  diff_shift_z = max(abs(tdiff_shift_z)) ;
486  } // End of relativistic equations
487 
488 
489  } // End of main loop
490 
491  //=========================================================================
492  // End of iteration
493  //=========================================================================
494 }
495 
496 // Truc pourri
497 void Et_bin_nsbh::equilibrium_nsbh (double, int, int, double,
498  int, double, double, const Tbl&, Tbl&) {
499 
500  cout << "Not implemented !" << endl ;
501  abort() ;
502 }
503 }
Cmp log(const Cmp &)
Neperian logarithm.
Definition: cmp_math.C:299
const Base_vect & ref_triad
Reference triad ("absolute frame"), with respect to which the components of all the member Tenseur &#39;s...
Definition: etoile.h:831
Cmp exp(const Cmp &)
Exponential.
Definition: cmp_math.C:273
Cmp sqrt(const Cmp &)
Square root.
Definition: cmp_math.C:223
void set_std_base()
Set the standard spectal basis of decomposition for each component.
Definition: tenseur.C:1186
Tenseur pot_centri
Centrifugal potential.
Definition: etoile.h:956
Lorene prototypes.
Definition: app_hor.h:67
Standard units of space, time and mass.
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
Definition: map.h:783
Tenseur nnn
Total lapse function.
Definition: etoile.h:512
Tenseur s_euler
Trace of the stress tensor in the Eulerian frame.
Definition: etoile.h:471
Tenseur flat_scalar_prod(const Tenseur &t1, const Tenseur &t2)
Scalar product of two Tenseur when the metric is : performs the contraction of the last index of t1 w...
virtual void adapt(const Cmp &ent, const Param &par, int nbr=0)=0
Adaptation of the mapping to a given scalar field.
Tenseur press
Fluid pressure.
Definition: etoile.h:464
Tenseur n_auto
Part of the lapse { N} generated principaly by the star.
Definition: et_bin_nsbh.h:85
Tbl diffrel(const Cmp &a, const Cmp &b)
Relative difference between two Cmp (norme version).
Definition: cmp_math.C:507
Tenseur confpsi
Total conformal factor $$.
Definition: et_bin_nsbh.h:101
Tenseur d_confpsi_comp
Gradient of { confpsi_comp} (Cartesian components with respect to { ref_triad})
Definition: et_bin_nsbh.h:119
Tenseur shift_auto
Part of the shift vector generated principaly by the star.
Definition: etoile.h:892
bool irrotational
true for an irrotational star, false for a corotating one
Definition: etoile.h:825
Cmp & set()
Read/write for a scalar (see also operator=(const Cmp&) ).
Definition: tenseur.C:840
Tenseur u_euler
Fluid 3-velocity with respect to the Eulerian observer.
Definition: etoile.h:477
void change_triad(const Base_vect &new_triad)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
Definition: tenseur.C:684
Map & mp
Mapping associated with the star.
Definition: etoile.h:432
int get_nzone() const
Returns the number of domains.
Definition: grilles.h:465
virtual void equation_of_state()
Computes the proper baryon and energy density, as well as pressure from the enthalpy.
Definition: etoile.C:569
Tbl max(const Cmp &)
Maximum values of a Cmp in each domain.
Definition: cmp_math.C:438
int get_etat() const
Returns the logical state.
Definition: tenseur.h:710
Tenseur khi_shift
Scalar used in the decomposition of shift_auto , following Shibata&#39;s prescription [Prog...
Definition: etoile.h:921
Cmp pow(const Cmp &, int)
Power .
Definition: cmp_math.C:351
virtual void resize(int l, double lambda)=0
Rescales the outer boundary of one domain.
Tenseur d_confpsi_auto
Gradient of { confpsi_auto} (Cartesian components with respect to { ref_triad})
Definition: et_bin_nsbh.h:114
Tenseur d_n_auto
Gradient of { n_auto} (Cartesian components with respect to { ref_triad})
Definition: et_bin_nsbh.h:93
int nzet
Number of domains of *mp occupied by the star.
Definition: etoile.h:435
bool relativistic
Indicator of relativity: true for a relativistic star, false for a Newtonian one. ...
Definition: etoile.h:440
Tenseur confpsi_auto
Part of the conformal factor $$ generated principaly by the star.
Definition: et_bin_nsbh.h:104
virtual void reevaluate_symy(const Map *mp_prev, int nzet, Cmp &uu) const =0
Recomputes the values of a Cmp at the collocation points after a change in the mapping.
const Base_vect_cart & get_bvect_cart() const
Returns the Cartesian basis associated with the coordinates (x,y,z) of the mapping, i.e.
Definition: map.h:809
Cmp ssjm1_confpsi
Effective source at the previous step for the resolution of the Poisson equation for { confpsi_auto} ...
Definition: et_bin_nsbh.h:164
Tenseur loggam
Logarithm of the Lorentz factor between the fluid and the co-orbiting observer.
Definition: etoile.h:852
Tenseur ent
Log-enthalpy (relativistic case) or specific enthalpy (Newtonian case)
Definition: etoile.h:460
Cmp abs(const Cmp &)
Absolute value.
Definition: cmp_math.C:413
Cmp ssjm1_lapse
Effective source at the previous step for the resolution of the Poisson equation for { n_auto} by mea...
Definition: et_bin_nsbh.h:158
double velocity_potential(int mermax, double precis, double relax)
Computes the non-translational part of the velocity scalar potential by solving the continuity equat...
Tenseur_sym tkij_tot
Total extrinsic curvature tensor $K^{ij}$ generated by { shift_auto} and { shift_comp}.
Definition: et_bin_nsbh.h:152
Basic array class.
Definition: tbl.h:164
Tenseur_sym tkij_auto
Part of the extrinsic curvature tensor $K^{ij}$ generated by { shift_auto}.
Definition: et_bin_nsbh.h:146
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Definition: tenseur.C:652
Tenseur ener_euler
Total energy density in the Eulerian frame.
Definition: etoile.h:468
Tenseur n_comp
Part of the lapse { N} generated principaly by the companion star.
Definition: et_bin_nsbh.h:88
virtual void equilibrium_nsbh(double ent_c, int mermax, int mermax_poisson, double relax_poisson, int mermax_potvit, double relax_potvit, double thres_adapt, const Tbl &fact, Tbl &diff)
Computes an equilibrium configuration in a NS-BH binary system.
Valeur va
The numerical value of the Cmp.
Definition: cmp.h:464
Cmp ssjm1_khi
Effective source at the previous step for the resolution of the Poisson equation for the scalar by m...
Definition: etoile.h:976
virtual void hydro_euler()
Computes the hydrodynamical quantities relative to the Eulerian observer from those in the fluid fram...
Definition: et_bin_hydro.C:109
Tenseur ssjm1_wshift
Effective source at the previous step for the resolution of the vector Poisson equation for by means...
Definition: etoile.h:986