diff --git a/src/main.cpp b/src/main.cpp index 443dcc37f8c31decf8b7f90fa76ebd0906224577..9d4682375c52ba917d97f40bf2e797e4c32c5058 100644 --- a/src/main.cpp +++ b/src/main.cpp @@ -170,6 +170,7 @@ int main(int argc, char *argv[]) t += dt; ++iteration; + std::cout << "temps t : " << t << std::endl; } @@ -197,7 +198,7 @@ int main(int argc, char *argv[]) double pi = 4.*std::atan(1.); std::ofstream fout("comparaison u"); for (size_t j=0; j<mesh.numberOfCells(); ++j) { - //fout << xj[j][0] << ' ' << uj[j][0] << ' ' << std::sin(pi*xj[j][0])*std::exp(-2.*pi*pi*0.2) <<'\n'; //cas k constant + //fout << xj[j][0] << ' ' << uj[j][0] << ' ' << std::sin(pi*xj[j][0])*std::exp(-2.*pi*pi*0.001) <<'\n'; //cas k constant fout << xj[j][0] << ' ' << uj[j][0] << ' ' << std::sin(pi*xj[j][0])*std::exp(-0.2) <<'\n'; // cas k non constant } } diff --git a/src/mesh/Mesh.hpp b/src/mesh/Mesh.hpp index 238a2adb4058c9a564f8f88bd344eb4d27ac8217..1061d0afd2fe1fe02edfe77358a032ee9f91ff67 100644 --- a/src/mesh/Mesh.hpp +++ b/src/mesh/Mesh.hpp @@ -50,15 +50,15 @@ public: // pas constant - Mesh(const Connectivity& connectivity) - : m_connectivity(connectivity), - m_xr("xr", connectivity.numberOfNodes()) - { - const double delta_x = 1./connectivity.numberOfCells(); - Kokkos::parallel_for(connectivity.numberOfNodes(), KOKKOS_LAMBDA(const int& r){ - m_xr[r][0] = r*delta_x; - }); - } + // Mesh(const Connectivity& connectivity) + //: m_connectivity(connectivity), + // m_xr("xr", connectivity.numberOfNodes()) + //{ + // const double delta_x = 1./connectivity.numberOfCells(); + // Kokkos::parallel_for(connectivity.numberOfNodes(), KOKKOS_LAMBDA(const int& r){ + // m_xr[r][0] = r*delta_x; + // }); + //} // pas non constant @@ -77,6 +77,18 @@ public: // }); //} + // pas non constant avec fonction x^2 + + Mesh(const Connectivity& connectivity) + : m_connectivity(connectivity), + m_xr("xr", connectivity.numberOfNodes()) + { + const double delta_x = 1./connectivity.numberOfCells(); + Kokkos::parallel_for(connectivity.numberOfNodes(), KOKKOS_LAMBDA(const int& r){ + m_xr[r][0] = (r*delta_x)*(r*delta_x); + }); + } + ~Mesh() { ; diff --git a/src/scheme/FiniteVolumesDiffusion.hpp b/src/scheme/FiniteVolumesDiffusion.hpp index 554879eac359b5962a0fdb246540ea0cab8a6178..adc46f0417c5a4f3dba6c2e4fc7c2b2eb8661b13 100644 --- a/src/scheme/FiniteVolumesDiffusion.hpp +++ b/src/scheme/FiniteVolumesDiffusion.hpp @@ -214,16 +214,17 @@ public: for (int ll=0; ll<cell_nb_faces(j); ++ll) { minVl = std::min(minVl, Vl(cell_faces(j, ll))); } + //k=2 => (kj(j+1) + 2*kj(j) + kj(j-1)) = 8 - // dt_j[j]= 0.5*rhoj(j)*Vj(j)*(2./8.)*minVl; + //dt_j[j]= 0.5*rhoj(j)*Vj(j)*(2./8.)*minVl; - // k=x + // k pas forcement constant double sum = 0.; for (int m = 0; m < cell_nb_nodes(j); ++m) { sum += kj(cell_nodes(j,m)); } - dt_j[j]= 0.5*rhoj(j)*Vj(j)*(1./(2.*kj(j) + sum))*minVl; + dt_j[j]= 0.5*rhoj(j)*Vj(j)*(1./(2.*kj(j) + sum))*minVl; }); @@ -279,7 +280,7 @@ public: momentum_fluxes += Fl(l)*Cjr(j,R); energy_fluxes += (Gl(l), Cjr(j,R)); } - uj[j] += std::exp(-t)*(dt*inv_mj[j])*Vj(j)*Sj(j) + (dt*inv_mj[j]) * momentum_fluxes; + uj[j] += std::exp(-t)*(dt*inv_mj[j])*Vj(j)*Sj(j) + (dt*inv_mj[j]) * momentum_fluxes; // test avec k non constant //uj[j] += (dt*inv_mj[j]) * momentum_fluxes; Ej[j] += (dt*inv_mj[j]) * energy_fluxes; }); @@ -313,13 +314,17 @@ public: double erreur = 0.; double exacte = 0.; for (size_t j=0; j<m_mesh.numberOfCells(); ++j) { - exacte = std::sin(pi*xj[j][0])*std::exp(-0.2); // solution exacte + //exacte = std::sin(pi*xj[j][0])*std::exp(-2.*pi*pi*0.001); // solution exacte cas test k constant + exacte = std::sin(pi*xj[j][0])*std::exp(-0.2); // solution exacte cas test k non constant erreur += (exacte - uj[j][0])*(exacte - uj[j][0]); } erreur = std::sqrt(erreur); return erreur; } + // Calcul erreur entre solution analytique et solution numerique en norme L infini (max) + // (quand la solution exacte est connue) + double error_Linf(UnknownsType& unknowns) { Kokkos::View<Rd*> uj = unknowns.uj(); @@ -327,10 +332,12 @@ public: const Kokkos::View<const Rd*> xj = m_mesh_data.xj(); double pi = 4.*std::atan(1.); - double exacte = std::sin(pi*xj[0][0])*std::exp(-0.2); + //double exacte = std::sin(pi*xj[0][0])*std::exp(-2.*pi*pi*0.001); // k constant + double exacte = std::sin(pi*xj[0][0])*std::exp(-0.2); // k non constant double erreur = std::abs(exacte - uj[0][0]); for (size_t j=1; j<m_mesh.numberOfCells(); ++j) { + //exacte = std::sin(pi*xj[j][0])*std::exp(-2.*pi*pi*0.001); exacte = std::sin(pi*xj[j][0])*std::exp(-0.2); if (std::abs(exacte - uj[j][0]) > erreur) { diff --git a/src/scheme/FiniteVolumesEulerUnknowns.hpp b/src/scheme/FiniteVolumesEulerUnknowns.hpp index 1ef763027290bac95df07d4d87f166259aa4188d..5f8e2fceefe406e0a7effe78c12cc57d8255ae7c 100644 --- a/src/scheme/FiniteVolumesEulerUnknowns.hpp +++ b/src/scheme/FiniteVolumesEulerUnknowns.hpp @@ -228,7 +228,7 @@ void initializeSod() }); Kokkos::parallel_for(m_mesh.numberOfCells(), KOKKOS_LAMBDA(const int& j){ - //m_kj[j] = 2; // k constant + //m_kj[j] = 2.; // k constant m_kj[j] = xj[j][0]; // k non constant, k = x });