135 #include "metrique.h" 138 #include "graphique.h" 139 #include "utilitaires.h" 147 int mermax_poisson,
double relax_poisson,
148 double relax_potvit,
double thres_adapt,
149 Tbl& diff,
double om) {
166 int i_b = mg->
get_nr(l_b) - 1 ;
176 double& diff_ent = diff.
set(0) ;
177 double& diff_vel_pot = diff.
set(1) ;
178 double& diff_logn = diff.
set(2) ;
179 double& diff_lnq = diff.
set(3) ;
180 double& diff_beta_x = diff.
set(4) ;
181 double& diff_beta_y = diff.
set(5) ;
182 double& diff_beta_z = diff.
set(6) ;
183 double& diff_h11 = diff.
set(7) ;
184 double& diff_h21 = diff.
set(8) ;
185 double& diff_h31 = diff.
set(9) ;
186 double& diff_h22 = diff.
set(10) ;
187 double& diff_h32 = diff.
set(11) ;
188 double& diff_h33 = diff.
set(12) ;
203 int nz_search =
nzet ;
206 double precis_secant = 1.e-14 ;
208 double reg_map = 1. ;
213 par_adapt.
add_int(nitermax, 0) ;
218 par_adapt.
add_int(nz_search, 2) ;
220 par_adapt.
add_int(adapt_flag, 3) ;
240 par_adapt.
add_tbl(ent_limit, 0) ;
264 double precis_poisson = 1.e-16 ;
271 par_logn.
add_int(mermax_poisson, 0) ;
282 par_lnq.
add_int(mermax_poisson, 0) ;
293 par_beta2.
add_int(mermax_poisson, 0) ;
301 for (
int i=0; i<3; i++) {
302 ssjm1wbeta.set(i) =
Cmp(ssjm1_wbeta(i+1)) ;
313 par_h11.
add_int(mermax_poisson, 0) ;
324 par_h21.
add_int(mermax_poisson, 0) ;
335 par_h31.
add_int(mermax_poisson, 0) ;
346 par_h22.
add_int(mermax_poisson, 0) ;
357 par_h32.
add_int(mermax_poisson, 0) ;
368 par_h33.
add_int(mermax_poisson, 0) ;
398 for(
int mer=0 ; mer<mermax ; mer++ ) {
400 cout <<
"-----------------------------------------------" << endl ;
401 cout <<
"step: " << mer << endl ;
402 cout <<
"diff_ent = " << diff_ent << endl ;
432 double alpha_r2 = 0 ;
433 for (
int k=0; k<mg->
get_np(l_b); k++) {
434 for (
int j=0; j<mg->
get_nt(l_b); j++) {
439 double alpha_r2_jk = ( ent_c - ent_b + pot_ext_c - pot_ext_b) /
441 ( logn_auto_b - logn_auto_c ) ;
443 if (alpha_r2_jk > alpha_r2) {
444 alpha_r2 = alpha_r2_jk ;
452 alpha_r =
sqrt(alpha_r2) ;
454 cout <<
"k_b, j_b, alpha_r: " << k_b <<
" " << j_b <<
" " 475 ent = (ent_c + logn_auto_c + pot_ext_c) -
logn_auto - pot_ext ;
477 cout <<
"pot" <<
norme(pot_ext) << endl ;
490 double rap_dent = fabs( dent_eq / dent_pole ) ;
491 cout <<
"| dH/dr_eq / dH/dr_pole | = " << rap_dent << endl ;
493 if ( rap_dent < thres_adapt ) {
495 cout <<
"******* FROZEN MAPPING *********" << endl ;
503 for (
int l=0; l<
nzet; l++) {
506 ent_limit.
set(
nzet-1) = ent_b ;
519 double ray_eqq =
ray_eq() ;
521 double new_rr_out_2 = (separation - ray_eqq) * 0.95 ;
522 double new_rr_out_3 = (separation + ray_eqq_pi) * 1.05 ;
524 double rr_in_1 =
mp.
val_r(1,-1., M_PI/2, 0.) ;
525 double rr_out_1 =
mp.
val_r(1, 1., M_PI/2, 0.) ;
526 double rr_out_2 =
mp.
val_r(2, 1., M_PI/2, 0.) ;
527 double rr_out_3 =
mp.
val_r(3, 1., M_PI/2, 0.) ;
529 mp.
resize(1, 0.5*(new_rr_out_2 + rr_in_1) / rr_out_1) ;
530 mp.
resize(2, new_rr_out_2 / rr_out_2) ;
531 mp.
resize(3, new_rr_out_3 / rr_out_3) ;
533 for (
int dd=4; dd<=nz-2; dd++) {
535 mp.
val_r(dd, 1., M_PI/2, 0.)) ;
540 cout <<
"too small number of domains" << endl ;
553 double ent_s_max = -1 ;
556 for (
int k=0; k<mg->
get_np(l_b); k++) {
557 for (
int j=0; j<mg->
get_nt(l_b); j++) {
559 if (xx > ent_s_max) {
566 cout <<
"Max. abs(enthalpy) at the boundary between domains nzet-1" 567 <<
" and nzet : " << endl ;
568 cout <<
" " << ent_s_max <<
" reached for k = " << k_s_max <<
569 " and j = " << j_s_max << endl ;
644 double lambda_dirac = 0. ;
683 0, dcov_logn_auto, 0,
true) ;
685 source4 = -
contract(
hij, 0, 1, dcovdcov_logn_auto +
688 source5 = -
contract(hdirac, 0, dcov_logn_auto, 0) ;
690 source_tot = source1 + source2 + source3 + source4 + source5 ;
693 cout <<
"moyenne de la source 1 pour logn_auto" << endl ;
694 cout <<
norme(source1/(nr*nt*np)) << endl ;
695 cout <<
"moyenne de la source 2 pour logn_auto" << endl ;
696 cout <<
norme(source2/(nr*nt*np)) << endl ;
697 cout <<
"moyenne de la source 3 pour logn_auto" << endl ;
698 cout <<
norme(source3/(nr*nt*np)) << endl ;
699 cout <<
"moyenne de la source 4 pour logn_auto" << endl ;
700 cout <<
norme(source4/(nr*nt*np)) << endl ;
701 cout <<
"moyenne de la source 5 pour logn_auto" << endl ;
702 cout <<
norme(source5/(nr*nt*np)) << endl ;
703 cout <<
"moyenne de la source pour logn_auto" << endl ;
704 cout <<
norme(source_tot/(nr*nt*np)) << endl ;
712 cout <<
"logn_auto" << endl <<
norme(
logn_auto/(nr*nt*np)) << endl ;
719 "Relative error in the resolution of the equation for logn_auto : " 721 for (
int l=0; l<nz; l++) {
722 cout << tdiff_logn(l) <<
" " ;
725 diff_logn =
max(
abs(tdiff_logn)) ;
738 source3 = -
contract(dcon_lnq, 0, dcov_lnq_auto, 0,
true) ;
741 dcov_phi_auto + dcov_logn_auto, 0,
true) ;
744 dcov_hij_auto, 0, 1, dcov_gtilde, 0, 1), 0, 1) ;
747 dcov_hij_auto, 0, 1, dcov_gtilde, 0, 2), 0, 1) ;
749 source7 = -
contract(
hij, 0, 1, dcovdcov_lnq_auto + dcov_lnq_auto *
752 source8 = - 0.25 *
contract(dcov_hdirac_auto, 0, 1)
753 -
contract(hdirac, 0, dcov_lnq_auto, 0) ;
755 source_tot = source1 + source2 + source3 + source4 + source5 +
756 source6 + source7 + source8 ;
759 cout <<
"moyenne de la source 1 pour lnq_auto" << endl ;
760 cout <<
norme(source1/(nr*nt*np)) << endl ;
761 cout <<
"moyenne de la source 2 pour lnq_auto" << endl ;
762 cout <<
norme(source2/(nr*nt*np)) << endl ;
763 cout <<
"moyenne de la source 3 pour lnq_auto" << endl ;
764 cout <<
norme(source3/(nr*nt*np)) << endl ;
765 cout <<
"moyenne de la source 4 pour lnq_auto" << endl ;
766 cout <<
norme(source4/(nr*nt*np)) << endl ;
767 cout <<
"moyenne de la source 5 pour lnq_auto" << endl ;
768 cout <<
norme(source5/(nr*nt*np)) << endl ;
769 cout <<
"moyenne de la source 6 pour lnq_auto" << endl ;
770 cout <<
norme(source6/(nr*nt*np)) << endl ;
771 cout <<
"moyenne de la source 7 pour lnq_auto" << endl ;
772 cout <<
norme(source7/(nr*nt*np)) << endl ;
773 cout <<
"moyenne de la source 8 pour lnq_auto" << endl ;
774 cout <<
norme(source8/(nr*nt*np)) << endl ;
775 cout <<
"moyenne de la source pour lnq_auto" << endl ;
776 cout <<
norme(source_tot/(nr*nt*np)) << endl ;
785 cout <<
"lnq_auto" << endl <<
norme(
lnq_auto/(nr*nt*np)) << endl ;
792 "Relative error in the resolution of the equation for lnq : " 794 for (
int l=0; l<nz; l++) {
795 cout << tdiff_lnq(l) <<
" " ;
798 diff_lnq =
max(
abs(tdiff_lnq)) ;
815 source1_beta = (4.*qpig) *
nn %
psi4 824 source4_beta = -
contract(
hij, 0, 1, dcovdcov_beta_auto, 1, 2) ;
827 dcovdcov_beta_auto, 0, 1), 0) ;
829 source6_beta.set_etat_zero() ;
830 source6_beta.inc_dzpuis() ;
833 + 2./3. * hdirac * divbeta_auto
834 -
contract(hdirac, 0, dcov_beta_auto, 1) ;
836 source_beta = source1_beta + source2_beta + source3_beta
837 + source4_beta + source5_beta + source6_beta + source7_beta ;
840 cout <<
"moyenne de la source 1 pour beta_auto" << endl ;
841 cout <<
norme(source1_beta(1)/(nr*nt*np)) << endl ;
842 cout <<
norme(source1_beta(2)/(nr*nt*np)) << endl ;
843 cout <<
norme(source1_beta(3)/(nr*nt*np)) << endl ;
844 cout <<
"moyenne de la source 2 pour beta_auto" << endl ;
845 cout <<
norme(source2_beta(1)/(nr*nt*np)) << endl ;
846 cout <<
norme(source2_beta(2)/(nr*nt*np)) << endl ;
847 cout <<
norme(source2_beta(3)/(nr*nt*np)) << endl ;
848 cout <<
"moyenne de la source 3 pour beta_auto" << endl ;
849 cout <<
norme(source3_beta(1)/(nr*nt*np)) << endl ;
850 cout <<
norme(source3_beta(2)/(nr*nt*np)) << endl ;
851 cout <<
norme(source3_beta(3)/(nr*nt*np)) << endl ;
852 cout <<
"moyenne de la source 4 pour beta_auto" << endl ;
853 cout <<
norme(source4_beta(1)/(nr*nt*np)) << endl ;
854 cout <<
norme(source4_beta(2)/(nr*nt*np)) << endl ;
855 cout <<
norme(source4_beta(3)/(nr*nt*np)) << endl ;
856 cout <<
"moyenne de la source 5 pour beta_auto" << endl ;
857 cout <<
norme(source5_beta(1)/(nr*nt*np)) << endl ;
858 cout <<
norme(source5_beta(2)/(nr*nt*np)) << endl ;
859 cout <<
norme(source5_beta(3)/(nr*nt*np)) << endl ;
860 cout <<
"moyenne de la source 6 pour beta_auto" << endl ;
861 cout <<
norme(source6_beta(1)/(nr*nt*np)) << endl ;
862 cout <<
norme(source6_beta(2)/(nr*nt*np)) << endl ;
863 cout <<
norme(source6_beta(3)/(nr*nt*np)) << endl ;
864 cout <<
"moyenne de la source 7 pour beta_auto" << endl ;
865 cout <<
norme(source7_beta(1)/(nr*nt*np)) << endl ;
866 cout <<
norme(source7_beta(2)/(nr*nt*np)) << endl ;
867 cout <<
norme(source7_beta(3)/(nr*nt*np)) << endl ;
868 cout <<
"moyenne de la source pour beta_auto" << endl ;
869 cout <<
norme(source_beta(1)/(nr*nt*np)) << endl ;
870 cout <<
norme(source_beta(2)/(nr*nt*np)) << endl ;
871 cout <<
norme(source_beta(3)/(nr*nt*np)) << endl ;
878 for (
int i=1; i<=3; i++) {
879 if (source_beta(i).get_etat() != ETATZERO)
880 source_beta.
set(i).filtre(4) ;
883 for (
int i=1; i<=3; i++) {
884 if(source_beta(i).dz_nonzero()) {
885 assert( source_beta(i).get_dzpuis() == 4 ) ;
888 (source_beta.set(i)).set_dzpuis(4) ;
892 double lambda = double(1) / double(3) ;
896 for (
int i=0; i<3; i++) {
897 source_p.set(i) =
Cmp(source_beta(i+1)) ;
902 for (
int i=0; i<3 ;i++){
903 vect_auxi.set(i) = 0. ;
912 for (
int i=0; i<3 ;i++)
913 resu_p.set(i).annule_hard() ;
914 resu_p.set_std_base() ;
919 source_p.poisson_vect_oohara(lambda, par_beta2, resu_p, scal_auxi) ;
921 for (
int i=1; i<=3; i++)
925 for (
int i=0; i<3; i++){
926 ssjm1_wbeta.
set(i+1) = ssjm1wbeta(i) ;
944 Tbl tdiff_beta_x =
diffrel(lap_beta(1), source_beta(1)) ;
945 Tbl tdiff_beta_y =
diffrel(lap_beta(2), source_beta(2)) ;
946 Tbl tdiff_beta_z =
diffrel(lap_beta(3), source_beta(3)) ;
949 "Relative error in the resolution of the equation for beta_auto : " 951 cout <<
"x component : " ;
952 for (
int l=0; l<nz; l++) {
953 cout << tdiff_beta_x(l) <<
" " ;
956 cout <<
"y component : " ;
957 for (
int l=0; l<nz; l++) {
958 cout << tdiff_beta_y(l) <<
" " ;
961 cout <<
"z component : " ;
962 for (
int l=0; l<nz; l++) {
963 cout << tdiff_beta_z(l) <<
" " ;
967 diff_beta_x =
max(
abs(tdiff_beta_x)) ;
968 diff_beta_y =
max(
abs(tdiff_beta_y)) ;
969 diff_beta_z =
max(
abs(tdiff_beta_z)) ;
998 source_1 =
contract(dcon_hij, 1, dcov_lnq_auto, 0) ;
1000 source_2 = -
contract(dcon_hij, 2, dcov_lnq_auto, 0)
1011 double r0 =
mp.
val_r(nz-2, 1, 0, 0) ;
1012 double sigma = 1.*r0 ;
1018 ff =
exp( -(rr - r0)*(rr - r0)/sigma/sigma ) ;
1019 for (
int ii=0; ii<nz-1; ii++)
1020 ff.set_domain(ii) = 1. ;
1021 ff.set_outer_boundary(nz-1, 0) ;
1022 ff.std_spectral_base() ;
1036 dcov_omdsdphi.
set(1,1) = 0. ;
1037 dcov_omdsdphi.
set(2,1) = om * ff ;
1038 dcov_omdsdphi.set(3,1) = 0. ;
1039 dcov_omdsdphi.set(2,2) = 0. ;
1040 dcov_omdsdphi.set(3,2) = 0. ;
1041 dcov_omdsdphi.set(3,3) = 0. ;
1042 dcov_omdsdphi.set(1,2) = -om * ff ;
1043 dcov_omdsdphi.set(1,3) = 0. ;
1044 dcov_omdsdphi.set(2,3) = 0. ;
1045 dcov_omdsdphi.std_spectral_base() ;
1048 source_3a.inc_dzpuis() ;
1062 omdsdp.
set(1) = - om * yya * ff ;
1063 omdsdp.set(2) = om * xxa * ff ;
1064 omdsdp.set(3).annule_hard() ;
1067 omdsdp.set(1) = om * yya * ff ;
1068 omdsdp.set(2) = - om * xxa * ff ;
1069 omdsdp.set(3).annule_hard() ;
1072 omdsdp.set(1).set_outer_boundary(nz-1, 0) ;
1073 omdsdp.set(2).set_outer_boundary(nz-1, 0) ;
1074 omdsdp.std_spectral_base() ;
1085 source_5 = dcon_hdirac_auto ;
1101 source_Sij += -
nn / (3.*
psi4) * gtilde_con *
1103 dcov_gtilde, 0, 1), 0, 1)
1105 dcov_gtilde, 0, 2), 0, 1)) ;
1107 source_Sij += - 8.*
nn / (3.*
psi4) * gtilde_con *
1113 source_Sij += tens_temp ;
1115 source_Sij += - 8./(3.*
psi4) *
contract(dcov_phi_auto, 0,
1122 - 0.33333333333333333 *
s_euler * gtilde_con ) ;
1124 source_Sij += - 1./(
psi4*psi2) *
contract(gtilde_con, 1,
1125 contract(gtilde_con, 1, qq*dcovdcov_lnq_auto +
1126 qq*dcov_lnq_auto*dcov_lnq, 0), 1) ;
1129 dcov_hij_auto, 2), 1, dcov_qq, 0) -
1131 hij, 1), 1, dcov_qq, 0) ;
1134 dcov_hij_auto, 2), 0, dcov_qq, 0) ;
1136 source_Sij += 1./(3.*
psi4*psi2)*
contract(gtilde_con, 0, 1,
1137 qq*dcovdcov_lnq_auto + qq*dcov_lnq_auto*dcov_lnq, 0, 1)
1154 source_Rij +=
contract(hdirac, 0, dcov_hij_auto, 2) ;
1160 gtilde_con, 0, dcov_hij_auto, 2), 0, dcov_hij, 2), 1, 3) ;
1163 dcov_hij_auto, 1), 2, gtilde_con, 1), 0, dcov_hij, 1), 0, 3) +
1165 dcov_hij_auto, 1), 0, dcov_hij, 1), 0, 3), 0, gtilde_con, 1) ;
1167 source_Rij += 0.5 *
contract(gtilde_con*gtilde_con, 1, 3,
1168 contract(dcov_hij_auto, 0, 1, dcov_gtilde, 0, 1), 0, 1) ;
1170 source_Rij = source_Rij * 0.5 ;
1172 for(
int i=1; i<=3; i++)
1173 for(
int j=1; j<=i; j++) {
1175 source_tot_hij = source_1(i,j) + source_1(j,i)
1176 + source_2(i,j) + 2.*
psi4/
nn * (
1177 source_4(i,j) - source_Sij(i,j))
1178 - 2.* source_Rij(i,j) +
1179 source_5(i,j) + source_5(j,i) + source_6(i,j) ;
1182 source3 = 2.*
psi4/
nn * (source_3a(i,j) + source_3a(j,i) +
1185 source_tot_hij = source_tot_hij + source3 ;
1189 cout <<
"source_mat" << endl
1191 - 0.33333333333333333 *
s_euler * gtilde_con ))
1192 (i,j))/(nr*nt*np) << endl ;
1193 cout <<
"max source_mat" << endl
1195 - 0.33333333333333333 *
s_euler * gtilde_con ))
1198 cout <<
"source1" << endl
1199 <<
norme(source_1(i,j)/(nr*nt*np)) << endl ;
1200 cout <<
"max source1" << endl
1201 <<
max(source_1(i,j)) << endl ;
1202 cout <<
"source2" << endl
1203 <<
norme(source_2(i,j)/(nr*nt*np)) << endl ;
1204 cout <<
"max source2" << endl
1205 <<
max(source_2(i,j)) << endl ;
1206 cout <<
"source3a" << endl
1207 <<
norme(source_3a(i,j)/(nr*nt*np)) << endl ;
1208 cout <<
"max source3a" << endl
1209 <<
max(source_3a(i,j)) << endl ;
1210 cout <<
"source3b" << endl
1211 <<
norme(source_3b(i,j)/(nr*nt*np)) << endl ;
1212 cout <<
"max source3b" << endl
1213 <<
max(source_3b(i,j)) << endl ;
1214 cout <<
"source4" << endl
1215 <<
norme(source_4(i,j)/(nr*nt*np)) << endl ;
1216 cout <<
"max source4" << endl
1217 <<
max(source_4(i,j)) << endl ;
1218 cout <<
"source5" << endl
1219 <<
norme(source_5(i,j)/(nr*nt*np)) << endl ;
1220 cout <<
"max source5" << endl
1221 <<
max(source_5(i,j)) << endl ;
1222 cout <<
"source6" << endl
1223 <<
norme(source_6(i,j)/(nr*nt*np)) << endl ;
1224 cout <<
"max source6" << endl
1225 <<
max(source_6(i,j)) << endl ;
1226 cout <<
"source_Rij" << endl
1227 <<
norme(source_Rij(i,j)/(nr*nt*np)) << endl ;
1228 cout <<
"max source_Rij" << endl
1229 <<
max(source_Rij(i,j)) << endl ;
1230 cout <<
"source_Sij" << endl
1231 <<
norme(source_Sij(i,j)/(nr*nt*np)) << endl ;
1232 cout <<
"max source_Sij" << endl
1233 <<
max(source_Sij(i,j)) << endl ;
1234 cout <<
"source_tot" << endl
1235 <<
norme(source_tot_hij/(nr*nt*np)) << endl ;
1236 cout <<
"max source_tot" << endl
1237 <<
max(source_tot_hij) << endl ;
1245 source_tot_hij.poisson(par_h11,
hij_auto.
set(1,1)) ;
1249 Tbl tdiff_h11 =
diffrel(laplac, source_tot_hij) ;
1250 cout <<
"Relative error in the resolution of the equation for " 1251 <<
"h11_auto : " << endl ;
1252 for (
int l=0; l<nz; l++) {
1253 cout << tdiff_h11(l) <<
" " ;
1256 diff_h11 =
max(
abs(tdiff_h11)) ;
1261 source_tot_hij.poisson(par_h21,
hij_auto.
set(2,1)) ;
1265 Tbl tdiff_h21 =
diffrel(laplac, source_tot_hij) ;
1268 "Relative error in the resolution of the equation for " 1269 <<
"h21_auto : " << endl ;
1270 for (
int l=0; l<nz; l++) {
1271 cout << tdiff_h21(l) <<
" " ;
1274 diff_h21 =
max(
abs(tdiff_h21)) ;
1279 source_tot_hij.poisson(par_h31,
hij_auto.
set(3,1)) ;
1283 Tbl tdiff_h31 =
diffrel(laplac, source_tot_hij) ;
1286 "Relative error in the resolution of the equation for " 1287 <<
"h31_auto : " << endl ;
1288 for (
int l=0; l<nz; l++) {
1289 cout << tdiff_h31(l) <<
" " ;
1292 diff_h31 =
max(
abs(tdiff_h31)) ;
1297 source_tot_hij.poisson(par_h22,
hij_auto.
set(2,2)) ;
1301 Tbl tdiff_h22 =
diffrel(laplac, source_tot_hij) ;
1304 "Relative error in the resolution of the equation for " 1305 <<
"h22_auto : " << endl ;
1306 for (
int l=0; l<nz; l++) {
1307 cout << tdiff_h22(l) <<
" " ;
1310 diff_h22 =
max(
abs(tdiff_h22)) ;
1315 source_tot_hij.poisson(par_h32,
hij_auto.
set(3,2)) ;
1319 Tbl tdiff_h32 =
diffrel(laplac, source_tot_hij) ;
1322 "Relative error in the resolution of the equation for " 1323 <<
"h32_auto : " << endl ;
1324 for (
int l=0; l<nz; l++) {
1325 cout << tdiff_h32(l) <<
" " ;
1328 diff_h32 =
max(
abs(tdiff_h32)) ;
1333 source_tot_hij.poisson(par_h33,
hij_auto.
set(3,3)) ;
1337 Tbl tdiff_h33 =
diffrel(laplac, source_tot_hij) ;
1340 "Relative error in the resolution of the equation for " 1341 <<
"h33_auto : " << endl ;
1342 for (
int l=0; l<nz; l++) {
1343 cout << tdiff_h33(l) <<
" " ;
1346 diff_h33 =
max(
abs(tdiff_h33)) ;
1351 cout <<
"Tenseur hij auto in cartesian coordinates" << endl ;
1352 for (
int i=1; i<=3; i++)
1353 for (
int j=1; j<=i; j++) {
1354 cout <<
" Comp. " << i <<
" " << j <<
" : " ;
1355 for (
int l=0; l<nz; l++){
1379 diff_ent = diff_ent_tbl(0) ;
1380 for (
int l=1; l<
nzet; l++) {
1381 diff_ent += diff_ent_tbl(l) ;
Vector dcov_phi
Covariant derivative of the logarithm of the conformal factor.
Coord xa
Absolute x coordinate.
bool irrotational
true for an irrotational star, false for a corotating one
virtual const Sym_tensor & con() const
Read-only access to the contravariant representation.
void add_tenseur_mod(Tenseur &ti, int position=0)
Adds the address of a new modifiable Tenseur to the list.
Component of a tensorial field *** DEPRECATED : use class Scalar instead ***.
Cmp exp(const Cmp &)
Exponential.
int & set(int i)
Read/write of a particular element (index i ) (1D case)
Radial mapping of rather general form.
void add_int(const int &n, int position=0)
Adds the address of a new int to the list.
Map & mp
Mapping associated with the star.
Scalar ssjm1_lnq
Effective source at the previous step for the resolution of the Poisson equation for lnq_auto...
Sym_tensor hij_auto
Deviation of the inverse conformal metric from the inverse flat metric generated principally by the ...
int get_np(int l) const
Returns the number of points in the azimuthal direction ( ) in domain no. l.
Cmp sqrt(const Cmp &)
Square root.
void set_std_base()
Set the standard spectal basis of decomposition for each component.
Scalar ssjm1_h31
Effective source at the previous step for the resolution of the Poisson equation for h20_auto...
Scalar ssjm1_h11
Effective source at the previous step for the resolution of the Poisson equation for h00_auto...
Standard units of space, time and mass.
Sym_tensor stress_euler
Spatial part of the stress-energy tensor with respect to the Eulerian observer.
const Mg3d * get_mg() const
Gives the Mg3d on which the mapping is defined.
const Scalar & laplacian(int ced_mult_r=4) const
Returns the Laplacian of *this.
double & set(int i)
Read/write of a particular element (index i) (1D case)
Tensor field of valence 0 (or component of a tensorial field).
Scalar poisson() const
Solves the scalar Poisson equation with *this as a source.
double get_ori_x() const
Returns the x coordinate of the origin.
virtual void adapt(const Cmp &ent, const Param &par, int nbr=0)=0
Adaptation of the mapping to a given scalar field.
Scalar ssjm1_h33
Effective source at the previous step for the resolution of the Poisson equation for h22_auto...
const Tensor_sym & derive_con(const Metric &gam) const
Returns the "contravariant" derivative of this with respect to some metric , by raising the last inde...
Scalar logn_comp
Part of the lapse logarithm (gravitational potential at the Newtonian limit) generated principally by...
Basic integer array class.
Scalar loggam
Logarithm of the Lorentz factor between the fluid and the co-orbiting observer.
virtual void std_spectral_base()
Sets the spectral bases of the Valeur va to the standard ones for a scalar field. ...
Scalar kcar_auto
Part of the scalar generated by beta_auto, i.e.
void annule_hard()
Sets the Cmp to zero in a hard way.
double get_rot_phi() const
Returns the angle between the x –axis and X –axis.
Scalar s_euler
Trace of the stress scalar in the Eulerian frame.
void set(const Map *mp, Mtbl *(*construct)(const Map *))
Semi-constructor from a mapping and a method.
const Vector & derive_con(const Metric &gam) const
Returns the "contravariant" derivative of *this with respect to some metric , by raising the index of...
Tensor field of valence 1.
Scalar logn_auto
Part of the lapse logarithm (gravitational potential at the Newtonian limit) generated principally by...
const Tensor_sym & derive_cov(const Metric &gam) const
Returns the covariant derivative of this with respect to some metric .
void inc_dzpuis()
dzpuis += 1 ;
void smooth(int nzet, Valeur &uuva) const
Changes the function *this as a smooth one when there exists a discontinuity between the nucleus and ...
const Tensor_sym & get_delta() const
Returns the tensor which defines the connection with respect to the flat one: is the difference bet...
Tbl diffrel(const Cmp &a, const Cmp &b)
Relative difference between two Cmp (norme version).
Scalar ssjm1_h21
Effective source at the previous step for the resolution of the Poisson equation for h10_auto...
double val_grid_point(int l, int k, int j, int i) const
Returns the value of the field at a specified grid point.
virtual void dec_dzpuis(int dec=1)
Decreases by dec units the value of dzpuis and changes accordingly the values in the compactified ext...
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
int nzet
Number of domains of *mp occupied by the star.
virtual const Sym_tensor & con() const
Read-only access to the contravariant representation.
virtual double val_r(int l, double xi, double theta, double pphi) const =0
Returns the value of the radial coordinate r for a given in a given domain.
Cmp & set()
Read/write for a scalar (see also operator=(const Cmp&) ).
const Vector & divergence(const Metric &) const
Returns the divergence of this with respect to a Metric .
Tbl norme(const Cmp &)
Sums of the absolute values of all the values of the Cmp in each domain.
Sym_tensor hij
Total deviation of the inverse conformal metric from the inverse flat metric.
Scalar pot_centri
Centrifugal potential.
Scalar press
Fluid pressure.
Metric_flat flat
Flat metric defined on the mapping (Spherical components with respect to the mapping of the star) ...
Scalar lnq_auto
Scalar field generated principally by the star.
void equilibrium(double ent_c, int mermax, int mermax_potvit, int mermax_poisson, double relax_poisson, double relax_potvit, double thres_adapt, Tbl &diff, double om)
Computes an equilibrium configuration.
void add_tbl(const Tbl &ti, int position=0)
Adds the address of a new Tbl to the list.
int get_nzone() const
Returns the number of domains.
virtual void homothetie(double lambda)
Sets a new radial scale.
Tbl max(const Cmp &)
Maximum values of a Cmp in each domain.
const Scalar & divergence(const Metric &) const
The divergence of this with respect to a Metric .
Vector u_euler
Fluid 3-velocity with respect to the Eulerian observer.
Sym_tensor tkij_comp
Part of the extrinsic curvature tensor generated by beta_comp.
double ray_eq_pi() const
Coordinate radius at , [r_unit].
Cmp pow(const Cmp &, int)
Power .
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
Vector dcov_logn
Covariant derivative of the total logarithm of the lapse.
Tenseur contract(const Tenseur &, int id1, int id2)
Self contraction of two indices of a Tenseur .
const Tensor & derive_cov(const Metric &gam) const
Returns the covariant derivative of this with respect to some metric .
Scalar logn
Logarithm of the lapse N .
virtual void resize(int l, double lambda)=0
Rescales the outer boundary of one domain.
double val_point(double r, double theta, double phi) const
Computes the value of the field at an arbitrary point , by means of the spectral expansion.
Metric gtilde
Conformal metric .
virtual void inc_dzpuis(int inc=1)
Increases by inc units the value of dzpuis and changes accordingly the values in the compactified ext...
virtual void hydro_euler()
Computes the hydrodynamical quantities relative to the Eulerian observer from those in the fluid fram...
Vector dcon_logn
Contravariant derivative of the total logarithm of the lapse.
double ray_pole() const
Coordinate radius at [r_unit].
int get_nr(int l) const
Returns the number of points in the radial direction ( ) in domain no. l.
virtual const Sym_tensor & cov() const
Read-only access to the covariant representation.
Coord ya
Absolute y coordinate.
const Tensor & derive_con(const Metric &gam) const
Returns the "contravariant" derivative of this with respect to some metric , by raising the last inde...
double velocity_potential(int mermax, double precis, double relax)
Computes the non-translational part of the velocity scalar potential by solving the continuity equat...
Tbl & set(int l)
Read/write of the value in a given domain.
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.
Scalar & set(const Itbl &ind)
Returns the value of a component (read/write version).
Scalar nn
Lapse function N .
Scalar ssjm1_h22
Effective source at the previous step for the resolution of the Poisson equation for h11_auto...
Cmp abs(const Cmp &)
Absolute value.
Sym_tensor tkij_auto
Part of the extrinsic curvature tensor generated by beta_auto.
Coord za
Absolute z coordinate.
void add_double(const double &x, int position=0)
Adds the the address of a new double to the list.
bool conf_flat
true if the 3-metric is conformally flat, false for a more general metric.
double ray_eq() const
Coordinate radius at , [r_unit].
Scalar ssjm1_h32
Effective source at the previous step for the resolution of the Poisson equation for h21_auto...
const Scalar & dsdr() const
Returns of *this .
Vector beta_auto
Part of the shift vector generated principally by the star (Spherical components with respect to the ...
Scalar kcar_comp
Part of the scalar generated by beta_auto and beta_comp, i.e.
int get_nt(int l) const
Returns the number of points in the co-latitude direction ( ) in domain no. l.
Valeur & set_spectral_va()
Returns va (read/write version)
Scalar ssjm1_khi
Effective source at the previous step for the resolution of the Poisson equation for khi...
Scalar & set(int)
Read/write access to a component.
void set_etat_qcq()
Sets the logical state to ETATQCQ (ordinary state).
void add_cmp_mod(Cmp &ti, int position=0)
Adds the address of a new modifiable Cmp to the list.
void equation_of_state()
Computes the proper baryon and energy density, as well as pressure from the enthalpy.
const Vector & derive_cov(const Metric &gam) const
Returns the gradient (1-form = covariant vector) of *this
Scalar ener_euler
Total energy density in the Eulerian frame.
virtual void dec_dzpuis(int dec=1)
Decreases by dec units the value of dzpuis and changes accordingly the values of the Scalar in the co...
Sym_tensor derive_lie(const Vector &v) const
Computes the Lie derivative of this with respect to some vector field v.
Scalar psi4
Conformal factor .
Class intended to describe valence-2 symmetric tensors.
Tensor handling *** DEPRECATED : use class Tensor instead ***.
void add_int_mod(int &n, int position=0)
Adds the address of a new modifiable int to the list.
Vector dcon_phi
Contravariant derivative of the logarithm of the conformal factor.
Coord r
r coordinate centered on the grid
Scalar ssjm1_logn
Effective source at the previous step for the resolution of the Poisson equation for logn_auto...