29char star_bhns_vel_pot_C[] =
"$Header: /cvsroot/Lorene/C++/Source/Star_bhns/star_bhns_vel_pot.C,v 1.3 2014/10/13 08:53:41 j_novak Exp $" ;
60#include "utilitaires.h"
65Cmp raccord_c1(
const Cmp& uu,
int l1) ;
68 bool kerrschild,
int mermax,
double precis,
99 v_orb.set_etat_qcq() ;
101 for (
int i=1;
i<=3;
i++) {
102 v_orb.set(
i) = www(
i).val_grid_point(0, 0, 0, 0) ;
106 v_orb.std_spectral_base() ;
121 zeta_h.std_spectral_base() ;
131 cout <<
"!!!!! WARNING: Not yet available !!!!!" <<
endl ;
146 for (
int i=1;
i<=3;
i++) {
151 dentmb.std_spectral_base() ;
170 par.add_double(precis, 0) ;
171 par.add_double(relax, 1) ;
205 source.set_spectral_va().ylm_i() ;
209 cout <<
"Check of the resolution of the continuity equation : "
218 for (
int i=1;
i<=3;
i++)
228 for (
int i=1;
i<=3;
i++) {
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Time evolution with partial storage (*** under development ***).
const Base_vect_cart & get_bvect_cart() const
Returns the Cartesian basis associated with the coordinates (x,y,z) of the mapping,...
virtual void poisson_compact(const Cmp &source, const Cmp &aa, const Tenseur &bb, const Param &par, Cmp &psi) const =0
Resolution of the elliptic equation in the case where the stellar interior is covered by a single do...
const Base_vect_spher & get_bvect_spher() const
Returns the orthonormal vectorial basis associated with the coordinates of the mapping.
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
const Metric_flat & flat_met_spher() const
Returns the flat metric associated with the spherical coordinates and with components expressed in th...
Flat metric for tensor calculation.
Tensor field of valence 0 (or component of a tensorial field).
const Vector & derive_cov(const Metric &gam) const
Returns the gradient (1-form = covariant vector) of *this
const Scalar & deriv(int i) const
Returns of *this , where .
const Scalar & laplacian(int ced_mult_r=4) const
Returns the Laplacian of *this.
virtual void std_spectral_base()
Sets the spectral bases of the Valeur va to the standard ones for a scalar field.
int get_etat() const
Returns the logical state ETATNONDEF (undefined), ETATZERO (null) or ETATQCQ (ordinary).
const Vector & derive_con(const Metric &gam) const
Returns the "contravariant" derivative of *this with respect to some metric , by raising the index of...
Metric_flat flat
Flat metric defined on the mapping (Spherical components with respect to the mapping of the star ).
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, .
Vector d_confo_comp
Derivative of the conformal factor generated by the companion black hole.
Scalar confo_tot
Total conformal factor.
Scalar psi4
Fourth power of the total conformal factor.
Vector d_lapconf_auto
Derivative of the lapconf function generated by the star .
Vector d_psi
Gradient of (in the irrotational case) (Spherical components with respect to the mapping of the star...
Vector d_confo_auto
Derivative of the conformal factor generated by the star .
Scalar lapconf_tot
Total lapconf function.
Vector d_lapconf_comp
Derivative of the lapconf function generated by the companion black hole.
Scalar psi0
Scalar potential of the non-translational part of fluid 4-velocity (in the irrotational case)
const Eos & eos
Equation of state of the stellar matter.
Scalar gam_euler
Lorentz factor between the fluid and Eulerian observers.
Map & mp
Mapping associated with the star.
int nzet
Number of domains of *mp occupied by the star.
Tensor handling *** DEPRECATED : use class Tensor instead ***.
Tensor field of valence 1.
virtual void change_triad(const Base_vect &)
Sets a new vectorial basis (triad) of decomposition and modifies the components accordingly.
Scalar & set(int)
Read/write access to a component.
Cmp exp(const Cmp &)
Exponential.
Tbl diffrel(const Cmp &a, const Cmp &b)
Relative difference between two Cmp (norme version).
Tbl norme(const Cmp &)
Sums of the absolute values of all the values of the Cmp in each domain.
Cmp log(const Cmp &)
Neperian logarithm.
virtual void annule(int l_min, int l_max)
Sets the Tensor to zero in several domains.
const Base_vect * get_triad() const
Returns the vectorial basis (triad) on which the components are defined.
Tenseur contract(const Tenseur &, int id1, int id2)
Self contraction of two indices of a Tenseur .
Standard units of space, time and mass.