From f79d1e8bfebce7da028a8be794ce56debf7f7fa4 Mon Sep 17 00:00:00 2001 From: Fanny CHOPOT <fanny.chopot.ocre@cea.fr> Date: Wed, 23 May 2018 09:11:55 +0200 Subject: [PATCH] correction solution exacte gnuplot et modification CL dans computeUr (Acoustic) --- src/main.cpp | 6 +++--- src/scheme/AcousticSolver.hpp | 29 ++++++++++++++++++++++++--- src/scheme/FiniteVolumesDiffusion.hpp | 2 -- 3 files changed, 29 insertions(+), 8 deletions(-) diff --git a/src/main.cpp b/src/main.cpp index 28ef2a46e..4fea916dc 100644 --- a/src/main.cpp +++ b/src/main.cpp @@ -265,7 +265,7 @@ int main(int argc, char *argv[]) fout.precision(15); for (size_t j=0; j<mesh.numberOfCells(); ++j) { //fout << xj[j][0] << ' ' << rhoj[j] << '\n'; - fout << xj[j][0] << ' ' << rhoj[j] << ' ' << std::sqrt(((4.*((xj[j][0]*xj[j][0])/(h*h)) + 100.-(xj[j][0]*xj[j][0])/(h*h))/100.)/h) << '\n'; + fout << xj[j][0] << ' ' << rhoj[j] << ' ' << std::sqrt((4.*((xj[j][0]*xj[j][0])/(h*h)) + 100.-(xj[j][0]*xj[j][0])/(h*h))/100.)/h << '\n'; } } @@ -279,7 +279,7 @@ int main(int argc, char *argv[]) //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] << '\n'; - fout << xj[j][0] << ' ' << uj[j][0] << ' ' << -(xj[j][0]*0.2)/(1.-0.2*0.2) << '\n'; + fout << xj[j][0] << ' ' << uj[j][0] << ' ' << -(xj[j][0]*0.2)/((50./9.)-0.2*0.2) << '\n'; } } @@ -294,7 +294,7 @@ int main(int argc, char *argv[]) //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] << ' ' << ((xj[j][0]*pi*pi*0.5)*(std::sin(pi*xj[j][0])*std::sin(pi*xj[j][0]) - std::cos(xj[j][0]*pi)*std::cos(pi*xj[j][0])) - pi*0.5*std::sin(pi*xj[j][0])*std::cos(pi*xj[j][0]))*(std::exp(-2.*0.2)-1.) + 2. <<'\n' ; // cas k non constant //fout << xj[j][0] << ' ' << Ej[j] << '\n'; - fout << xj[j][0] << ' ' << ej[j] << ' ' << std::sqrt(((4.*((xj[j][0]*xj[j][0])/(h*h)) + 100.-(xj[j][0]*xj[j][0])/(h*h))/100.)/h)*std::sqrt(((4.*((xj[j][0]*xj[j][0])/(h*h)) + 100.-(xj[j][0]*xj[j][0])/(h*h))/100.)/h) << '\n'; + fout << xj[j][0] << ' ' << ej[j] << ' ' << ( std::sqrt((4.*((xj[j][0]*xj[j][0])/(h*h)) + 100.-(xj[j][0]*xj[j][0])/(h*h))/100.)/h)*( std::sqrt((4.*((xj[j][0]*xj[j][0])/(h*h)) + 100.-(xj[j][0]*xj[j][0])/(h*h))/100.)/h) << '\n'; } } diff --git a/src/scheme/AcousticSolver.hpp b/src/scheme/AcousticSolver.hpp index f63ac8498..b0fd38286 100644 --- a/src/scheme/AcousticSolver.hpp +++ b/src/scheme/AcousticSolver.hpp @@ -147,6 +147,7 @@ private: return m_br; } + /* Kokkos::View<Rd*> // calcule u_r (vitesse au sommet du maillage et flux de vitesse) computeUr(const Kokkos::View<const Rdd*>& Ar, const Kokkos::View<const Rd*>& br) { @@ -160,6 +161,26 @@ private: return m_ur; } + */ + + Kokkos::View<Rd*> // calcule u_r (vitesse au sommet du maillage et flux de vitesse) + computeUr(const Kokkos::View<const Rdd*>& Ar, + const Kokkos::View<const Rd*>& br, + const double& t) { + + inverse(Ar, m_inv_Ar); + const Kokkos::View<const Rdd*> invAr = m_inv_Ar; + + Kokkos::View<Rd*> xr = m_mesh.xr(); + + Kokkos::parallel_for(m_mesh.numberOfNodes(), KOKKOS_LAMBDA(const int& r) { + m_ur[r]=invAr(r)*br(r); + }); + m_ur[0]=zero; + m_ur[m_mesh.numberOfNodes()-1] = (-t/((50./9.)-t*t))*xr[m_mesh.numberOfNodes()-1]; + + return m_ur; + } Kokkos::View<Rd**> // Fonction qui calcule F_jr computeFjr(const Kokkos::View<const Rdd**>& Ajr, @@ -207,7 +228,8 @@ private: const Kokkos::View<const double*>& pj, const Kokkos::View<const double*>& cj, const Kokkos::View<const double*>& Vj, - const Kokkos::View<const Rd**>& Cjr) { + const Kokkos::View<const Rd**>& Cjr, + const double& t) { const Kokkos::View<const double*> rhocj = computeRhoCj(rhoj, cj); const Kokkos::View<const Rdd**> Ajr = computeAjr(rhocj, Cjr); @@ -216,7 +238,8 @@ private: Kokkos::View<Rd*> ur = m_ur; Kokkos::View<Rd**> Fjr = m_Fjr; - ur = computeUr(Ar, br); + //ur = computeUr(Ar, br); + ur = computeUr(Ar, br, t); Fjr = computeFjr(Ajr, ur, Cjr, uj, pj); } @@ -282,7 +305,7 @@ public: Kokkos::View<Rd*> xr = m_mesh.xr(); // Calcule les flux - computeExplicitFluxes(xr, xj, rhoj, uj, pj, cj, Vj, Cjr); + computeExplicitFluxes(xr, xj, rhoj, uj, pj, cj, Vj, Cjr, t); const Kokkos::View<const Rd**> Fjr = m_Fjr; const Kokkos::View<const Rd*> ur = m_ur; diff --git a/src/scheme/FiniteVolumesDiffusion.hpp b/src/scheme/FiniteVolumesDiffusion.hpp index 3acadd1b5..b4a794d0f 100644 --- a/src/scheme/FiniteVolumesDiffusion.hpp +++ b/src/scheme/FiniteVolumesDiffusion.hpp @@ -278,8 +278,6 @@ public: Kokkos::View<double*> kL = unknowns.kL(); Kokkos::View<double*> kR = unknowns.kR(); - uR[0][0] = -t/(1.-t*t); - const Kokkos::View<const Rd*> xj = m_mesh_data.xj(); const Kokkos::View<const double*> Vj = m_mesh_data.Vj(); const Kokkos::View<const Rd**> Cjr = m_mesh_data.Cjr(); -- GitLab