23char bhole_coal_C[] =
"$Header: /cvsroot/Lorene/C++/Source/Bhole_binaire/bhole_coal.C,v 1.8 2014/10/13 08:52:40 j_novak Exp $" ;
129 cout <<
"TROUS STATIQUES : " << endl ;
139 for (
int i=1 ; i<nz_un ; i++)
140 if (diff_un(i) > erreur)
141 erreur = diff_un(i) ;
144 for (
int i=1 ; i<nz_deux ; i++)
145 if (diff_deux(i) > erreur)
146 erreur = diff_deux(i) ;
149 cout <<
"PAS TOTAL : " << conte <<
" DIFFERENCE : " << erreur << endl ;
157void Bhole_binaire::coal (
double precis,
double relax,
int nbre_ome,
double seuil_search,
double m1,
double m2,
const int sortie) {
159 assert (
omega == 0) ;
174 char name_iteration[40] ;
175 char name_correction[40] ;
176 char name_viriel[40] ;
178 char name_linear[40] ;
180 char name_error_m1[40] ;
181 char name_error_m2[40] ;
185 sprintf(name_iteration,
"ite.dat") ;
186 sprintf(name_correction,
"cor.dat") ;
187 sprintf(name_viriel,
"vir.dat") ;
188 sprintf(name_ome,
"ome.dat") ;
189 sprintf(name_linear,
"linear.dat") ;
190 sprintf(name_axe,
"axe.dat") ;
191 sprintf(name_error_m1,
"error_m1.dat") ;
192 sprintf(name_error_m2,
"error_m2.dat") ;
193 sprintf(name_r1,
"r1.dat") ;
194 sprintf(name_r2,
"r2.dat") ;
196 ofstream fiche_iteration(name_iteration) ;
197 fiche_iteration.precision(8) ;
199 ofstream fiche_correction(name_correction) ;
200 fiche_correction.precision(8) ;
202 ofstream fiche_viriel(name_viriel) ;
203 fiche_viriel.precision(8) ;
205 ofstream fiche_ome(name_ome) ;
206 fiche_ome.precision(8) ;
208 ofstream fiche_linear(name_linear) ;
209 fiche_linear.precision(8) ;
211 ofstream fiche_axe(name_axe) ;
212 fiche_axe.precision(8) ;
214 ofstream fiche_error_m1 (name_error_m1) ;
215 fiche_error_m1.precision(8) ;
217 ofstream fiche_error_m2 (name_error_m2) ;
218 fiche_error_m2.precision(8) ;
220 ofstream fiche_r1 (name_r1) ;
221 fiche_r1.precision(8) ;
223 ofstream fiche_r2 (name_r2) ;
224 fiche_r2.precision(8) ;
227 cout <<
"OMEGA AUGMENTE A LA MAIN." << endl ;
229 for (
int pas = 0 ; pas <nbre_ome ; pas ++) {
231 homme += angulaire/nbre_ome ;
245 for (
int i=1 ; i<nz1 ; i++)
246 if (diff_un(i) > erreur)
247 erreur = diff_un(i) ;
250 for (
int i=1 ; i<nz2 ; i++)
251 if (diff_deux(i) > erreur)
252 erreur = diff_deux(i) ;
254 double error_viriel =
viriel() ;
255 double error_linear = linear_momentum_systeme_inf()(1)/scale_linear ;
262 fiche_iteration << conte <<
" " << erreur << endl ;
264 fiche_viriel << conte <<
" " << error_viriel << endl ;
265 fiche_ome << conte <<
" " << homme << endl ;
266 fiche_linear << conte <<
" " << error_linear << endl ;
267 fiche_axe << conte <<
" " << pos_axe << endl ;
268 fiche_error_m1 << conte <<
" " << error_m1 << endl ;
269 fiche_error_m2 << conte <<
" " << error_m2 << endl ;
270 fiche_r1 << conte <<
" " << r1 << endl ;
271 fiche_r2 << conte <<
" " << r2 << endl ;
274 cout <<
"PAS TOTAL : " << conte <<
" DIFFERENCE : " << erreur << endl ;
279 cout <<
"OMEGA VARIABLE" << endl ;
296 for (
int i=1 ; i<nz1 ; i++)
297 if (diff_un(i) > erreur)
298 erreur = diff_un(i) ;
301 for (
int i=1 ; i<nz2 ; i++)
302 if (diff_deux(i) > erreur)
303 erreur = diff_deux(i) ;
305 double error_viriel =
viriel() ;
306 double error_linear = linear_momentum_systeme_inf()(1)/scale_linear ;
313 fiche_iteration << conte <<
" " << erreur << endl ;
315 fiche_viriel << conte <<
" " << error_viriel << endl ;
316 fiche_ome << conte <<
" " <<
omega << endl ;
317 fiche_linear << conte <<
" " << error_linear << endl ;
318 fiche_axe << conte <<
" " << pos_axe << endl ;
319 fiche_error_m1 << conte <<
" " << error_m1 << endl ;
320 fiche_error_m2 << conte <<
" " << error_m2 << endl ;
321 fiche_r1 << conte <<
" " << r1 << endl ;
322 fiche_r2 << conte <<
" " << r2 << endl ;
326 if (erreur <= seuil_search)
329 double scaling_ome =
pow((2-error_viriel)/(2-2*error_viriel), 1.) ;
332 double scaling_axe =
pow((2-error_linear)/(2-2*error_linear), 0.1) ;
333 set_pos_axe (pos_axe*scaling_axe) ;
335 double scaling_r1 =
pow((2-error_m1)/(2-2*error_m1), 0.1) ;
338 double scaling_r2 =
pow((2-error_m2)/(2-2*error_m2), 0.1) ;
342 cout <<
"PAS TOTAL : " << conte <<
" DIFFERENCE : " << erreur << endl ;
348 fiche_iteration.close() ;
349 fiche_correction.close() ;
350 fiche_viriel.close() ;
352 fiche_linear.close() ;
354 fiche_error_m1.close() ;
355 fiche_error_m2.close() ;
double omega
Position of the axis of rotation.
Bhole hole1
Black hole one.
double viriel() const
Computes the viriel error, that is the difference between the ADM and the Komar masses,...
void set_statiques(double precis, double relax)
Initialize the systeme to Misner Lindquist solution, that is solving for N and in the case .
void solve_shift(double precis, double relax)
Solves the equation for the shift, using the Oohara-Nakarmure scheme~:
void set_omega(double ome)
Sets the orbital velocity to ome.
Bhole hole2
Black hole two.
void solve_psi(double precis, double relax)
Solves the equation for the conformal factor~:
void coal(double precis, double relax, int nbre_ome, double search_ome, double m1, double m2, const int sortie=0)
Solves the equation for a particular angular velocity, the systeme being initialized to Misner-Lindqu...
void solve_lapse(double precis, double relax)
Solves the equation for the lapse~:
void fait_tkij()
Calculation af the extrinsic curvature tensor.
void init_bhole_binaire()
Initialisation of the system.
const Tenseur & get_n_auto() const
Returns the part of N generated by the hole.
double regul
Intensity of the correction on the shift vector.
double get_regul() const
Returns the intensity of the regularisation on .
double rayon
Radius of the horizon in LORENE's units.
double area() const
Computes the area of the throat.
Map_af & mp
Affine mapping.
const Tenseur & get_shift_auto() const
Returns the part of generated by the hole.
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
void homothetie_interne(double lambda)
Sets a new radial scale at the bondary between the nucleus and the first shell.
virtual double val_r(int l, double xi, double theta, double pphi) const
Returns the value of the radial coordinate r for a given in a given domain.
double get_ori_x() const
Returns the x coordinate of the origin.
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
int get_nzone() const
Returns the number of domains.
Cmp sqrt(const Cmp &)
Square root.
Cmp pow(const Cmp &, int)
Power .
Tbl diffrelmax(const Cmp &a, const Cmp &b)
Relative difference between two Cmp (max version).