61 int nz =
star1.mp.get_mg()->get_nzone() ;
62 int nr =
star1.mp.get_mg()->get_nr(0);
63 int nt =
star1.mp.get_mg()->get_nt(0);
64 int np =
star1.mp.get_mg()->get_np(0);
71 star1.hij_comp.set_triad(
star1.mp.get_bvect_cart()) ;
76 assert ( *(
star1.hij_comp.get_triad()) == *(comp_hij1.
get_triad())) ;
78 for(
int i=1; i<=3; i++)
79 for(
int j=i; j<=3; j++) {
80 star1.hij_comp.set(i,j).set_etat_qcq() ;
81 star1.hij_comp.set(i,j).import( (comp_hij1)(i,j) ) ;
83 star1.hij_comp.std_spectral_base() ;
85 for(
int i=1; i<=3; i++)
86 for(
int j=i; j<=3; j++)
91 star2.hij_comp.set_triad(
star2.mp.get_bvect_cart()) ;
96 assert ( *(
star2.hij_comp.get_triad()) == *(comp_hij2.
get_triad())) ;
98 for(
int i=1; i<=3; i++)
99 for(
int j=i; j<=3; j++) {
100 star2.hij_comp.set(i,j).set_etat_qcq() ;
101 star2.hij_comp.set(i,j).import( (comp_hij2)(i,j) ) ;
103 star2.hij_comp.std_spectral_base() ;
105 for(
int i=1; i<=3; i++)
106 for(
int j=i; j<=3; j++)
113 cout <<
"Function Binary::dirac_gauge()" << endl ;
119 double precis = 1e-5 ;
120 double precis_poisson = 1e-14 ;
121 double relax_poisson = 1.5 ;
122 int mer_poisson = 4 ;
136 Scalar ssjm1_xi11 (xi1(1)) ;
137 Scalar ssjm1_xi12 (xi1(2)) ;
138 Scalar ssjm1_xi13 (xi1(3)) ;
141 for(
int mer=0; mer<mermax; mer++){
148 double r0_1 =
star1.mp.val_r(nz-2, 1, 0, 0) ;
149 double sigma = 3.*r0_1 ;
152 ff1 =
exp( -(rr1 - r0_1)*(rr1 - r0_1)/sigma/sigma ) ;
153 for (
int ii=0; ii<nz-1; ii++)
166 source_xi1 += source_reg1 ;
170 Cmp ssjm1xi11 (ssjm1_xi11) ;
171 Cmp ssjm1xi12 (ssjm1_xi12) ;
172 Cmp ssjm1xi13 (ssjm1_xi13) ;
179 par_xi11.
add_int(mer_poisson, 0) ;
186 par_xi12.
add_int(mer_poisson, 0) ;
193 par_xi13.
add_int(mer_poisson, 0) ;
203 ssjm1_xi11 = ssjm1xi11 ;
204 ssjm1_xi12 = ssjm1xi12 ;
205 ssjm1_xi13 = ssjm1xi13 ;
213 Tbl tdiff_xi1_x =
diffrel(lap_xi1(1), source_xi1(1)) ;
214 Tbl tdiff_xi1_y =
diffrel(lap_xi1(2), source_xi1(2)) ;
215 Tbl tdiff_xi1_z =
diffrel(lap_xi1(3), source_xi1(3)) ;
218 "Relative error in the resolution of the equation for xi1 : "
220 cout <<
"x component : " ;
221 for (
int l=0; l<nz; l++) {
222 cout << tdiff_xi1_x(l) <<
" " ;
225 cout <<
"y component : " ;
226 for (
int l=0; l<nz; l++) {
227 cout << tdiff_xi1_y(l) <<
" " ;
230 cout <<
"z component : " ;
231 for (
int l=0; l<nz; l++) {
232 cout << tdiff_xi1_z(l) <<
" " ;
239 for (
int i=1 ; i<nz ; i++)
240 if (diff(i) > erreur)
243 cout <<
"Step : " << mer <<
" Difference : " << erreur << endl ;
244 cout <<
"-------------------------------------" << endl ;
260 Scalar ssjm1_xi21 (xi2(1)) ;
261 Scalar ssjm1_xi22 (xi2(2)) ;
262 Scalar ssjm1_xi23 (xi2(3)) ;
265 for(
int mer=0; mer<mermax; mer++){
272 double r0_2 =
star2.mp.val_r(nz-2, 1, 0, 0) ;
273 double sigma = 3.*r0_2 ;
276 ff2 =
exp( -(rr2 - r0_2)*(rr2 - r0_2)/sigma/sigma ) ;
277 for (
int ii=0; ii<nz-1; ii++)
290 source_xi2 += source_reg2 ;
294 Cmp ssjm1xi21 (ssjm1_xi21) ;
295 Cmp ssjm1xi22 (ssjm1_xi22) ;
296 Cmp ssjm1xi23 (ssjm1_xi23) ;
303 par_xi21.
add_int(mer_poisson, 0) ;
310 par_xi22.
add_int(mer_poisson, 0) ;
317 par_xi23.
add_int(mer_poisson, 0) ;
327 ssjm1_xi21 = ssjm1xi21 ;
328 ssjm1_xi22 = ssjm1xi22 ;
329 ssjm1_xi23 = ssjm1xi23 ;
337 Tbl tdiff_xi2_x =
diffrel(lap_xi2(1), source_xi2(1)) ;
338 Tbl tdiff_xi2_y =
diffrel(lap_xi2(2), source_xi2(2)) ;
339 Tbl tdiff_xi2_z =
diffrel(lap_xi2(3), source_xi2(3)) ;
342 "Relative error in the resolution of the equation for xi2 : "
344 cout <<
"x component : " ;
345 for (
int l=0; l<nz; l++) {
346 cout << tdiff_xi2_x(l) <<
" " ;
349 cout <<
"y component : " ;
350 for (
int l=0; l<nz; l++) {
351 cout << tdiff_xi2_y(l) <<
" " ;
354 cout <<
"z component : " ;
355 for (
int l=0; l<nz; l++) {
356 cout << tdiff_xi2_z(l) <<
" " ;
363 for (
int i=1 ; i<nz ; i++)
364 if (diff(i) > erreur)
367 cout <<
"Step : " << mer <<
" Difference : " << erreur << endl ;
368 cout <<
"-------------------------------------" << endl ;
382 guu_dirac1 =
star1.gamma.con().derive_lie(xi1) ;
384 guu_dirac1 = guu_dirac1 +
star1.gamma.con() ;
385 star1.gamma = guu_dirac1 ;
390 gtilde_con1 =
pow(
star1.gamma.determinant(), 1./3.) * guu_dirac1 ;
392 for(
int i=1; i<=3; i++)
393 for(
int j=i; j<=3; j++)
394 hij_dirac1.
set(i,j) = gtilde_con1(i,j) -
star1.flat.con()(i,j) ;
397 star1.gtilde = gtilde_con1 ;
399 star1.psi4.std_spectral_base() ;
401 cout <<
"norme de h_uu avant :" << endl ;
402 for (
int i=1; i<=3; i++)
403 for (
int j=1; j<=i; j++) {
404 cout <<
" Comp. " << i <<
" " << j <<
" : " ;
405 for (
int l=0; l<nz; l++){
406 cout <<
norme(
star1.hij(i,j)/(nr*nt*np))(l) <<
" " ;
412 cout <<
"norme de h_uu en jauge de dirac :" << endl ;
413 for (
int i=1; i<=3; i++)
414 for (
int j=1; j<=i; j++) {
415 cout <<
" Comp. " << i <<
" " << j <<
" : " ;
416 for (
int l=0; l<nz; l++){
417 cout <<
norme(hij_dirac1(i,j)/(nr*nt*np))(l) <<
" " ;
428 cout <<
"For comparaison H^i before computation = " << endl
429 <<
norme(hh_dirac(1))/(nr*nt*np)
431 <<
norme(hh_dirac(2))/(nr*nt*np)
433 <<
norme(hh_dirac(3))/(nr*nt*np)
437 cout <<
"Vector H^i after the computation" << endl ;
438 for (
int i=1; i<=3; i++){
439 cout <<
" Comp. " << i <<
" : " <<
norme(hh_dirac_new(i)
440 /(nr*nt*np)) << endl ;
446 (1 -
star1.decouple) ;
447 star1.hij = hij_dirac1 ;
454 guu_dirac2 =
star2.gamma.con().derive_lie(xi2) ;
456 guu_dirac2 = guu_dirac2 +
star2.gamma.con() ;
457 star2.gamma = guu_dirac2 ;
462 gtilde_con2 =
pow(
star2.gamma.determinant(), 1./3.) * guu_dirac2 ;
464 for(
int i=1; i<=3; i++)
465 for(
int j=i; j<=3; j++)
466 hij_dirac2.
set(i,j) = gtilde_con2(i,j) -
star2.flat.con()(i,j) ;
469 star2.gtilde = gtilde_con2 ;
471 star2.psi4.std_spectral_base() ;
477 (1 -
star2.decouple) ;
478 star2.hij = hij_dirac2 ;