24char binhor_coal_C[] =
"$Header: /cvsroot/Lorene/C++/Source/Bin_hor/binhor_coal.C,v 1.15 2014/10/13 08:52:42 j_novak Exp $" ;
97 double lim_nn,
int bound_psi) {
110 cout <<
"Static black holes : " << endl ;
121 for (
int i=1 ; i<nz ; i++)
122 if (diff(i) > erreur)
125 cout <<
"Step : " << conte <<
" Difference : " << erreur << endl ;
134 int nb_it,
int bound_nn,
double lim_nn,
135 int bound_psi,
int bound_beta,
double omega_eff,
137 ostream& fich_iteration, ostream& fich_correction,
138 ostream& fich_viriel, ostream& fich_kss,
139 int step,
int search_mass,
double mass_irr,
144 double precis = 1e-7 ;
147 cout <<
"OMEGA INCREASED STEP BY STEP." << endl ;
149 double inc_homme = (angu_vel - homme)/nb_ome ;
150 for (
int pas = 0 ; pas <nb_ome ; pas ++) {
153 if (omega_eff == alpha*homme ) verif = true ;
158 omega_eff = alpha*homme ;
161 solve_shift (precis, relax, bound_beta, omega_eff) ;
168 if (search_mass == 1 && step >= 30) {
171 double error_m = (mass_irr - mass_area) / mass_irr ;
172 double scaling_r =
pow((2-error_m)/(2-2*error_m), 1.) ;
188 cout <<
"Angular momentum computed at the horizon : " <<
ang_mom_hor()
193 for (
int i=1 ; i<nz ; i++)
194 if (diff(i) > erreur)
207 for (
int k=0 ; k<nnp ; k++)
208 for (
int j=0 ; j<nnt ; j++){
216 fich_iteration << step <<
" " <<
log10(erreur) <<
" " << homme << endl ;
217 fich_correction << step <<
" " <<
log10(
hole1.
regul) <<
" " << homme << endl ;
219 fich_viriel << step <<
" " <<
viriel() <<
" " << homme <<
" " <<
hole1.
omega_hor() - alpha*homme <<
" " << omega_eff << endl ;
220 fich_kss << step <<
" " << max_kss <<
" " << min_kss << endl ;
223 cout <<
"STEP : " << step <<
" DIFFERENCE : " << erreur << endl ;
230 cout <<
"OMEGA FIXED" << endl ;
233 for (
int pas = 0 ; pas <nb_it ; pas ++) {
237 solve_shift (precis, relax, bound_beta, omega_eff) ;
244 if (search_mass == 1 && step >= 30) {
247 double error_m = (mass_irr - mass_area) / mass_irr ;
248 double scaling_r =
pow((2-error_m)/(2-2*error_m), 1.) ;
258 for (
int i=1 ; i<nz ; i++)
259 if (diff(i) > erreur)
272 for (
int k=0 ; k<nnp ; k++)
273 for (
int j=0 ; j<nnt ; j++){
282 fich_iteration << step <<
" " <<
log10(erreur) <<
" " << homme << endl ;
283 fich_correction << step <<
" " <<
log10(
hole1.
regul) <<
" " << homme << endl ;
285 fich_viriel << step <<
" " <<
viriel() <<
" " << homme <<
" " <<
hole1.
omega_hor() - alpha*homme <<
" " << omega_eff << endl ;
286 fich_kss << step <<
" " << max_kss <<
" " << min_kss << endl ;
289 cout <<
"STEP : " << step <<
" DIFFERENCE : " << erreur << endl ;
294 fich_iteration <<
"#----------------------------" << endl ;
295 fich_correction <<
"#-----------------------------" << endl ;
296 fich_viriel <<
"#------------------------------" << endl ;
double viriel() const
Computes the viriel error, that is the difference between the ADM and the Komar masses,...
void solve_psi(double precis, double relax, int bound_psi)
Solves the equation for the conformal factor : The fields are the total values excpet those with subs...
double coal(double ang_vel, double relax, int nb_om, int nb_it, int bound_nn, double lim_nn, int bound_psi, int bound_beta, double omega_eff, double alpha, ostream &fich_iteration, ostream &fich_correction, ostream &fich_viriel, ostream &fich_kss, int step, int search_mass, double mass_irr, const int sortie=0)
Solves the equation for a particular angular velocity, the systeme being initialized to Misner-Lindqu...
Single_hor hole1
Black hole one.
Single_hor hole2
Black hole two.
double ang_mom_hor() const
Calculates the angular momentum of the black hole using the formula at the horizon.
void init_bin_hor()
Initialisation of the system.
void solve_shift(double precis, double relax, int bound_beta, double omega_eff)
Solves the equation for the shift, using the Oohara-Nakarmure scheme : The fields are the total value...
void set_statiques(double precis, double relax, int bound_nn, double lim_nn, int bound_psi)
Initialize the systeme to Misner Lindquist solution, that is solving for N and in the case .
void solve_lapse(double precis, double relax, int bound_nn, double lim_nn)
Solves the equation for the lapse : The fields are the total values except those with subscript ,...
void set_omega(double ome)
Sets the orbital velocity to ome.
void extrinsic_curvature()
Calculation of the extrinsic curvature tensor.
double get_omega() const
Returns the angular velocity.
void homothetie_interne(double lambda)
Sets a new radial scale at the bondary between the nucleus and the first shell.
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...
virtual const Vector & radial_vect() const
Returns the radial vector normal to a spherical slicing and pointing toward spatial infinity.
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
int get_nzone() const
Returns the number of domains.
Tensor field of valence 0 (or component of a tensorial field).
double val_grid_point(int l, int k, int j, int i) const
Returns the value of the field at a specified grid point.
double radius
Radius of the horizon in LORENE's units.
const Metric & get_gam() const
metric
Vector beta_auto
Shift function .
const Sym_tensor & get_k_dd() const
k_dd
Metric_flat ff
3 metric flat
Map_af & mp
Affine mapping.
void init_met_trK()
Sets the 3-metric tilde to the flat metric and gamt_point, trK and trK_point to zero.
double area_hor() const
Area of the horizon.
double regul
Intensity of the correction on the shift vector.
double omega_hor() const
Orbital velocity
Scalar n_auto
Lapse function .
Metric tgam
3 metric tilde
Cmp sqrt(const Cmp &)
Square root.
Cmp log10(const Cmp &)
Basis 10 logarithm.
Cmp pow(const Cmp &, int)
Power .
Tbl diffrelmax(const Cmp &a, const Cmp &b)
Relative difference between two Cmp (max version).
Tenseur contract(const Tenseur &, int id1, int id2)
Self contraction of two indices of a Tenseur .