From bd306980a37ced1d55717ca458854108baed9700 Mon Sep 17 00:00:00 2001 From: Stephane Del Pino <stephane.delpino44@gmail.com> Date: Mon, 23 Apr 2018 14:10:06 +0200 Subject: [PATCH] Fixed acoustic scheme for d>1 (introducing ljr=|Cjr| and njr = 1/ljr Cjr) --- src/mesh/MeshData.hpp | 22 +++++++++++++++++++++- src/scheme/AcousticSolver.hpp | 15 ++++++++++----- 2 files changed, 31 insertions(+), 6 deletions(-) diff --git a/src/mesh/MeshData.hpp b/src/mesh/MeshData.hpp index bed5cf97e..7fc8f8c80 100644 --- a/src/mesh/MeshData.hpp +++ b/src/mesh/MeshData.hpp @@ -20,6 +20,8 @@ public: private: const MeshType& m_mesh; Kokkos::View<Rd**> m_Cjr; + Kokkos::View<double**> m_ljr; + Kokkos::View<Rd**> m_njr; Kokkos::View<Rd*> m_xj; Kokkos::View<double*> m_Vj; @@ -74,7 +76,7 @@ private: KOKKOS_INLINE_FUNCTION void _updateCjr() { if(dimension == 1) { - // Cjr are constant overtime + // Cjr/njr/ljr are constant overtime } static_assert(dimension==1, "only 1d is implemented"); } @@ -90,6 +92,16 @@ public: return m_Cjr; } + const Kokkos::View<const double**> ljr() const + { + return m_ljr; + } + + const Kokkos::View<const Rd**> njr() const + { + return m_njr; + } + const Kokkos::View<const Rd*> xj() const { return m_xj; @@ -110,6 +122,8 @@ public: MeshData(const MeshType& mesh) : m_mesh(mesh), m_Cjr("Cjr", mesh.numberOfCells(), mesh.connectivity().maxNbNodePerCell()), + m_ljr("ljr", mesh.numberOfCells(), mesh.connectivity().maxNbNodePerCell()), + m_njr("njr", mesh.numberOfCells(), mesh.connectivity().maxNbNodePerCell()), m_xj("xj", mesh.numberOfCells()), m_Vj("Vj", mesh.numberOfCells()) { @@ -119,6 +133,12 @@ public: m_Cjr(j,0)=-1; m_Cjr(j,1)= 1; }); + // in 1d njr=Cjr + m_njr=m_Cjr; + Kokkos::parallel_for(m_mesh.numberOfCells(), KOKKOS_LAMBDA(const int& j) { + m_ljr(j,0)= 1; + m_ljr(j,1)= 1; + }); } this->updateAllData(); } diff --git a/src/scheme/AcousticSolver.hpp b/src/scheme/AcousticSolver.hpp index 1889f9f83..6f6d9a864 100644 --- a/src/scheme/AcousticSolver.hpp +++ b/src/scheme/AcousticSolver.hpp @@ -79,13 +79,14 @@ private: KOKKOS_INLINE_FUNCTION const Kokkos::View<const Rdd**> computeAjr(const Kokkos::View<const double*>& rhocj, - const Kokkos::View<const Rd**>& Cjr) { + const Kokkos::View<const double**>& ljr, + const Kokkos::View<const Rd**>& njr) { const Kokkos::View<const unsigned short*> cell_nb_nodes = m_connectivity.cellNbNodes(); Kokkos::parallel_for(m_mesh.numberOfCells(), KOKKOS_LAMBDA(const int& j) { for (int r=0; r<cell_nb_nodes[j]; ++r) { - m_Ajr(j,r) = tensorProduct(rhocj(j)*Cjr(j,r), Cjr(j,r)); + m_Ajr(j,r) = tensorProduct((rhocj(j)*ljr(j,r))*njr(j,r), njr(j,r)); } }); @@ -190,9 +191,11 @@ 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 Kokkos::View<const double**>& ljr, + const Kokkos::View<const Rd**>& njr) { const Kokkos::View<const double*> rhocj = computeRhoCj(rhoj, cj); - const Kokkos::View<const Rdd**> Ajr = computeAjr(rhocj, Cjr); + const Kokkos::View<const Rdd**> Ajr = computeAjr(rhocj, ljr, njr); const Kokkos::View<const Rdd*> Ar = computeAr(Ajr); const Kokkos::View<const Rd*> br = computeBr(Ajr, Cjr, uj, pj); @@ -259,9 +262,11 @@ public: 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(); + const Kokkos::View<const double**> ljr = m_mesh_data.ljr(); + const Kokkos::View<const Rd**> njr = m_mesh_data.njr(); Kokkos::View<Rd*> xr = m_mesh.xr(); - computeExplicitFluxes(xr, xj, rhoj, uj, pj, cj, Vj, Cjr); + computeExplicitFluxes(xr, xj, rhoj, uj, pj, cj, Vj, Cjr, ljr, njr); const Kokkos::View<const Rd**> Fjr = m_Fjr; const Kokkos::View<const Rd*> ur = m_ur; -- GitLab