diff --git a/src/main.cpp b/src/main.cpp index 2478cef554b7b62b35d731c26069aa731a296afe..7440d67736e1cc88935e721f02ec5bb6f897da09 100644 --- a/src/main.cpp +++ b/src/main.cpp @@ -142,7 +142,7 @@ int main(int argc, char *argv[]) const Kokkos::View<const double*> Vj = mesh_data.Vj(); const Kokkos::View<const Rd**> Cjr = mesh_data.Cjr(); - const double tmax=0.2; + const double tmax=0.8; double t=0.; int itermax=std::numeric_limits<int>::max(); @@ -169,6 +169,7 @@ int main(int argc, char *argv[]) const Kokkos::View<const Rd*> xj = mesh_data.xj(); int size = 3000; std::vector<std::vector<double>> x(size, std::vector<double>(mesh.numberOfCells())); + std::vector<std::vector<double>> rho_marche(size, std::vector<double>(mesh.numberOfCells())); std::vector<double> tempo(size); int i = 0; @@ -274,7 +275,7 @@ int main(int argc, char *argv[]) while((t<tmax) and (iteration<itermax)) { - /* + // ETAPE 1 DU SPLITTING - EULER double dt_euler = 0.9*acoustic_solver.acoustic_dt(Vj, cj); @@ -287,7 +288,7 @@ int main(int argc, char *argv[]) // ETAPE 2 DU SPLITTING - DIFFUSION - double dt_diff = 0.9*finite_volumes_diffusion.diffusion_dt(rhoj, kj,nuj, cj); + double dt_diff = 0.9*finite_volumes_diffusion.diffusion_dt(rhoj, kj,nuj, cj, nuL, nuR, kL,kR); double t_diff = t-dt_euler; if (dt_euler <= dt_diff) { @@ -295,7 +296,7 @@ int main(int argc, char *argv[]) finite_volumes_diffusion.computeNextStep(t_diff, dt_diff, unknowns); } else { while (t > t_diff) { - dt_diff = 0.9*finite_volumes_diffusion.diffusion_dt(rhoj, kj, nuj,cj); + dt_diff = 0.9*finite_volumes_diffusion.diffusion_dt(rhoj, kj, nuj,cj,nuL, nuR, kL,kR); if (t_diff+dt_diff > t) { dt_diff = t-t_diff; } @@ -303,9 +304,9 @@ int main(int argc, char *argv[]) t_diff += dt_diff; } } - */ - + + /* // DIFFUSION PURE double dt = 0.9*finite_volumes_diffusion.diffusion_dt(rhoj,kj,nuj,cj,nuL,nuR,kL,kR); @@ -314,7 +315,7 @@ int main(int argc, char *argv[]) } finite_volumes_diffusion.computeNextStep(t, dt, unknowns); t += dt; - + */ block_eos.updatePandCFromRhoE(); @@ -545,14 +546,15 @@ int main(int argc, char *argv[]) // ENTROPY TEST //finite_volumes_diffusion.entropie(unknowns); - /* + // STOCKAGE COORDONNES ET TEMPS for (size_t j=0; j<mesh.numberOfCells(); ++j) { x[i][j] = xj[j][0]; + rho_marche[i][j] = rhoj[j]; } tempo[i] = t; i = i + 1; - */ + } @@ -572,21 +574,21 @@ int main(int argc, char *argv[]) fout << ' ' << '\n'; } */ - /* - std::ofstream fout("cara1"); + + std::ofstream fout("cararho"); fout.precision(15); for (int j=0; j<mesh.numberOfCells(); ++j) { if (j%10 == 0) { for (int k=0; k<i; ++k) { if ( (k%10 == 0) or (k == i-1) ) { - fout << tempo[k] << ' ' << x[k][j] << '\n'; + fout << tempo[k] << ' ' << x[k][j] << ' ' << rho_marche[k][j] << '\n'; } } fout << ' ' << '\n'; fout << ' ' << '\n'; } } - */ + /* double error1 = 0.; diff --git a/src/mesh/Mesh.hpp b/src/mesh/Mesh.hpp index 176beae12e81f4faad30969301a82d1ea3e13778..adf4b739b6f8c5852f6affb41488335c72267441 100644 --- a/src/mesh/Mesh.hpp +++ b/src/mesh/Mesh.hpp @@ -76,7 +76,7 @@ public: // pas constant - /* + Mesh(const Connectivity& connectivity) : m_connectivity(connectivity), m_xr("xr", connectivity.numberOfNodes()), @@ -92,7 +92,7 @@ public: m_xr[r][0] = a + r*delta_x; }); } - */ + // pas non constant @@ -128,7 +128,7 @@ public: // pas aleatoire - + /* Mesh(const Connectivity& connectivity) : m_connectivity(connectivity), m_xr("xr", connectivity.numberOfNodes()), @@ -160,7 +160,7 @@ public: std::exit(0); } - + */ diff --git a/src/scheme/FiniteVolumesEulerUnknowns.hpp b/src/scheme/FiniteVolumesEulerUnknowns.hpp index 0e3719b8a4ab773074e395ccc6ca4a5eeb8d0c45..9464275392d87641cc1362a8707f06184f0d60a7 100644 --- a/src/scheme/FiniteVolumesEulerUnknowns.hpp +++ b/src/scheme/FiniteVolumesEulerUnknowns.hpp @@ -276,7 +276,7 @@ public: } - /* + // --- Acoustic Solver --- void initializeSod() @@ -338,7 +338,7 @@ public: //m_kj[j] = 0.; - + /* // Sod // k non regulier @@ -373,12 +373,12 @@ public: m_kj[j]=0. ; } } - + */ // k regulier int n = 1; m_kj[j] = std::exp(1.)*std::exp(-1./(1.-( (xj[j][0]-(0.7+0.1/n)) / (0.1/n) )*( (xj[j][0]-(0.7+0.1/n)) / (0.1/n) ))) * (xj[j][0]>0.7)*(xj[j][0]<0.7+0.1/n) + std::exp(1.)*std::exp(-1./(1.-( (xj[j][0]-(0.9-0.1/n)) / (0.1/n) )*( (xj[j][0]-(0.9-0.1/n)) / (0.1/n) ))) * (xj[j][0]>0.9-0.1/n)*(xj[j][0]<0.9) + (xj[j][0]>0.7+0.1/n)*(xj[j][0]<0.9-0.1/n); - m_kj[j] = 0.014*m_kj[j]; + m_kj[j] = 0.14*m_kj[j]; }); @@ -397,8 +397,8 @@ public: }); Kokkos::parallel_for(m_mesh.numberOfCells(), KOKKOS_LAMBDA(const int& j){ //m_nuj(j) = 0.5*(1.+xj[j][0]); - m_nuj(j) = 0.0005; - + m_nuj(j) = 0.; + /* if (xj[j][0]<0.7) { m_nuj[j]=0.; } else { @@ -408,7 +408,7 @@ public: m_nuj[j]=0. ; } } - + */ }); Kokkos::parallel_for(m_mesh.numberOfCells(), KOKKOS_LAMBDA(const int& j){ @@ -419,14 +419,14 @@ public: m_TL[0] = m_ej[0]; m_TR[0] = m_ej[m_mesh.numberOfCells()-1]; - m_nuL[0] = 0.0005; - m_nuR[0] = 0.0005; + m_nuL[0] = 0.; + m_nuR[0] = 0.; } - */ - + + /* // DIFFUSION PURE void initializeSod() @@ -499,7 +499,7 @@ public: m_nuR[0] = 0.; } - + */ FiniteVolumesEulerUnknowns(const MeshDataType& mesh_data)