diff --git a/src/main.cpp b/src/main.cpp index 28ef2a46e73bff2f49dbb52ed047e89c39d24367..4fea916dc4a45f3c6ad4b795031d1a1e3cd30e39 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 f63ac8498724f7c4037a8668b901234886fb7b54..b0fd382868987170161a72d2160a04b3e01b778f 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 3acadd1b5a63aaa1a1b446289d31c9166ca3b2c0..b4a794d0f5c254ba72014096d40cccc70ab97736 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();