From 961c80506a672a5168c4d3af00302dfad3ce13c6 Mon Sep 17 00:00:00 2001 From: Fanny CHOPOT <fanny.chopot.ocre@cea.fr> Date: Mon, 14 May 2018 15:33:49 +0200 Subject: [PATCH] ajout solution E pour k non constant --- src/main.cpp | 7 ++++--- src/mesh/Mesh.hpp | 8 ++++---- src/scheme/FiniteVolumesDiffusion.hpp | 19 ++++++++++--------- src/scheme/FiniteVolumesEulerUnknowns.hpp | 5 ++--- 4 files changed, 20 insertions(+), 19 deletions(-) diff --git a/src/main.cpp b/src/main.cpp index 71524b990..8dc0a24db 100644 --- a/src/main.cpp +++ b/src/main.cpp @@ -206,8 +206,8 @@ int main(int argc, char *argv[]) std::ofstream fout("comparaison u"); fout.precision(15); 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(-0.2) <<'\n'; // cas k non constant + // 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(-0.2) <<'\n'; // cas k non constant } } @@ -218,7 +218,8 @@ int main(int argc, char *argv[]) std::ofstream fout("comparaison E"); fout.precision(15); for (size_t j=0; j<mesh.numberOfCells(); ++j) { - fout << xj[j][0] << ' ' << Ej[j] << ' ' << (-(std::cos(pi*xj[j][0])*std::cos(pi*xj[j][0]))+(std::sin(pi*xj[j][0])*std::sin(pi*xj[j][0])))*0.5*(std::exp(-4.*pi*pi*0.2)-1.) + 2. <<'\n'; // cas k constant + //fout << xj[j][0] << ' ' << Ej[j] << ' ' << (-(std::cos(pi*xj[j][0])*std::cos(pi*xj[j][0]))+(std::sin(pi*xj[j][0])*std::sin(pi*xj[j][0])))*0.5*(std::exp(-4.*pi*pi*0.2)-1.) + 2. <<'\n'; // cas k constant + fout << xj[j][0] << ' ' << Ej[j] << ' ' << (0.5*pi*pi*xj[j][0]*(-std::cos(pi*xj[j][0])*std::cos(pi*xj[j][0]) + std::sin(pi*xj[j][0])*std::sin(pi*xj[j][0])) - cos(xj[j][0])*sin(xj[j][0])*pi*0.5)*(std::exp(-2.*0.2) - 1.) + 2. <<'\n' ; // cas k non constant } } diff --git a/src/mesh/Mesh.hpp b/src/mesh/Mesh.hpp index e9002afbe..b7dd7da9d 100644 --- a/src/mesh/Mesh.hpp +++ b/src/mesh/Mesh.hpp @@ -50,7 +50,7 @@ public: // pas constant - /* + Mesh(const Connectivity& connectivity) : m_connectivity(connectivity), m_xr("xr", connectivity.numberOfNodes()) @@ -60,10 +60,10 @@ public: m_xr[r][0] = r*delta_x; }); } -*/ - // pas non constant + // pas non constant + /* Mesh(const Connectivity& connectivity) : m_connectivity(connectivity), @@ -85,7 +85,7 @@ public: } }); } - + */ ~Mesh() { diff --git a/src/scheme/FiniteVolumesDiffusion.hpp b/src/scheme/FiniteVolumesDiffusion.hpp index f54b291d6..091870631 100644 --- a/src/scheme/FiniteVolumesDiffusion.hpp +++ b/src/scheme/FiniteVolumesDiffusion.hpp @@ -315,9 +315,10 @@ 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; // test avec k non constant - uj[j] += (dt*inv_mj[j]) * momentum_fluxes; - Ej[j] += (dt*inv_mj[j]) * energy_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; + }); // Calcul de e par la formule e = E-0.5 u^2 @@ -351,8 +352,8 @@ 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(-2.*pi*pi*0.2); // solution exacte cas test k constant - //exacte = std::sin(pi*xj[j][0])*std::exp(-0.2); // solution exacte cas test k non constant + //exacte = std::sin(pi*xj[j][0])*std::exp(-2.*pi*pi*0.2); // 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])*Vj(j); } erreur = std::sqrt(erreur); @@ -369,13 +370,13 @@ 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(-2.*pi*pi*0.2); // k constant - //double exacte = std::sin(pi*xj[0][0])*std::exp(-0.2); // k non constant + //double exacte = std::sin(pi*xj[0][0])*std::exp(-2.*pi*pi*0.2); // 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.2); - //exacte = std::sin(pi*xj[j][0])*std::exp(-0.2); + //exacte = std::sin(pi*xj[j][0])*std::exp(-2.*pi*pi*0.2); + exacte = std::sin(pi*xj[j][0])*std::exp(-0.2); if (std::abs(exacte - uj[j][0]) > erreur) { erreur = std::abs(exacte - uj[j][0]); diff --git a/src/scheme/FiniteVolumesEulerUnknowns.hpp b/src/scheme/FiniteVolumesEulerUnknowns.hpp index 5ea0d3390..71700408b 100644 --- a/src/scheme/FiniteVolumesEulerUnknowns.hpp +++ b/src/scheme/FiniteVolumesEulerUnknowns.hpp @@ -183,7 +183,6 @@ public: return m_kR; } - // --- Acoustic Solver --- @@ -277,8 +276,8 @@ void initializeSod() }); Kokkos::parallel_for(m_mesh.numberOfCells(), KOKKOS_LAMBDA(const int& j){ - m_kj[j] = 2.; // k constant - //m_kj[j] = xj[j][0]; // k non constant, k = x + //m_kj[j] = 2.; // k constant + m_kj[j] = xj[j][0]; // k non constant, k = x }); Kokkos::parallel_for(m_mesh.numberOfCells(), KOKKOS_LAMBDA(const int& j){ -- GitLab