From b6444edb28c6849196ad950b2d2854c3c27d165b Mon Sep 17 00:00:00 2001 From: Stephane Del Pino <stephane.delpino44@gmail.com> Date: Thu, 5 Apr 2018 19:27:03 +0200 Subject: [PATCH] AcousticSolverTest replaces AcousticSolver Began vector notations in AcousticSolverTest --- experimental/AcousticSolver.cpp | 120 ++++++++++++++-------------- experimental/AcousticSolver.hpp | 15 +++- experimental/AcousticSolverTest.cpp | 16 +++- experimental/AcousticSolverTest.hpp | 4 + 4 files changed, 87 insertions(+), 68 deletions(-) diff --git a/experimental/AcousticSolver.cpp b/experimental/AcousticSolver.cpp index 3088ff717..80ed8310e 100644 --- a/experimental/AcousticSolver.cpp +++ b/experimental/AcousticSolver.cpp @@ -1,8 +1,6 @@ #include <AcousticSolver.hpp> #include <rang.hpp> -#include <memory> - #include <BlockPerfectGas.hpp> typedef const double my_double; @@ -48,11 +46,10 @@ const Kokkos::View<const double*> AcousticSolver::computeRhoCj(const Kokkos::View<const double*>& rhoj, const Kokkos::View<const double*>& cj) { - Kokkos::View<double*> rhocj("rho_c", m_nj); Kokkos::parallel_for(m_nj, KOKKOS_LAMBDA(const int& j) { - rhocj[j] = rhoj[j]*cj[j]; + m_rhocj[j] = rhoj[j]*cj[j]; }); - return rhocj; + return m_rhocj; } KOKKOS_INLINE_FUNCTION @@ -60,14 +57,13 @@ const Kokkos::View<const double*[2]> AcousticSolver::computeAjr(const Kokkos::View<const double*>& rhocj, const Kokkos::View<const double*[2]>& Cjr) { - Kokkos::View<double*[2]> Ajr("Ajr", m_nj); Kokkos::parallel_for(m_nj, KOKKOS_LAMBDA(const int& j) { for (int r=0; r<2; ++r) { - Ajr(j,r) = rhocj(j)*Cjr(j,r)*Cjr(j,r); + m_Ajr(j,r) = rhocj(j)*Cjr(j,r)*Cjr(j,r); } }); - return Ajr; + return m_Ajr; } KOKKOS_INLINE_FUNCTION @@ -77,7 +73,6 @@ AcousticSolver::computeAr(const Kokkos::View<const double*[2]>& Ajr, const Kokkos::View<const int*[2]>& node_cell_local_node, const Kokkos::View<const int*>& node_nb_cells) { - Kokkos::View<double*> Ar("Ar", m_nr); Kokkos::parallel_for(m_nr, KOKKOS_LAMBDA(const int& r) { double sum = 0; for (int j=0; j<node_nb_cells(r); ++j) { @@ -85,10 +80,10 @@ AcousticSolver::computeAr(const Kokkos::View<const double*[2]>& Ajr, const int R = node_cell_local_node(r,j); sum += Ajr(J,R); } - Ar(r) = sum; + m_Ar(r) = sum; }); - return Ar; + return m_Ar; } KOKKOS_INLINE_FUNCTION @@ -101,18 +96,17 @@ AcousticSolver::computeBr(const Kokkos::View<const double*[2]>& Ajr, const Kokkos::View<const int*[2]>& node_cell_local_node, const Kokkos::View<const int*>& node_nb_cells) { - Kokkos::View<double*> br("br", m_nr); Kokkos::parallel_for(m_nr, KOKKOS_LAMBDA(const int& r) { double sum = 0; for (int j=0; j<node_nb_cells(r); ++j) { - const int J = node_cells(r,j); - const int R = node_cell_local_node(r,j); - sum += Ajr(J,R)*uj(J) + Cjr(J,R)*pj[J]; + const int J = node_cells(r,j); + const int R = node_cell_local_node(r,j); + sum += Ajr(J,R)*uj(J) + Cjr(J,R)*pj(J); } - br(r) = sum; + m_br(r) = sum; }); - return br; + return m_br; } KOKKOS_INLINE_FUNCTION @@ -120,15 +114,15 @@ Kokkos::View<double*> AcousticSolver::computeUr(const Kokkos::View<const double*>& Ar, const Kokkos::View<const double*>& br) { - const Kokkos::View<const double*> invAr = inverse(Ar); - Kokkos::View<double*> ur("ur", m_nr); + inverse(Ar, m_inv_Ar); + const Kokkos::View<const double*> invAr = m_inv_Ar; Kokkos::parallel_for(m_nr, KOKKOS_LAMBDA(const int& r) { - ur[r]=br(r)*invAr(r); + m_ur[r]=br(r)*invAr(r); }); - ur[0]=0; - ur[m_nr-1]=0; + m_ur[0]=0; + m_ur[m_nr-1]=0; - return ur; + return m_ur; } KOKKOS_INLINE_FUNCTION @@ -140,15 +134,13 @@ AcousticSolver::computeFjr(const Kokkos::View<const double*[2]>& Ajr, const Kokkos::View<const double*>& pj, const Kokkos::View<const int*[2]>& cell_nodes) { - Kokkos::View<double*[2]> Fjr("Fjr", m_nj); - Kokkos::parallel_for(m_nj, KOKKOS_LAMBDA(const int& j) { for (int r=0; r<2; ++r) { - Fjr(j,r) = Ajr(j,r)*(uj(j)-ur(cell_nodes(j,r)))+Cjr(j,r)*pj(j); + m_Fjr(j,r) = Ajr(j,r)*(uj(j)-ur(cell_nodes(j,r)))+Cjr(j,r)*pj(j); } }); - return Fjr; + return m_Fjr; } KOKKOS_INLINE_FUNCTION @@ -156,29 +148,24 @@ double AcousticSolver:: acoustic_dt(const Kokkos::View<const double*>& Vj, const Kokkos::View<const double*>& cj) const { - double dt = std::numeric_limits<double>::max(); - - Kokkos::View<double*> Vj_cj("Vj_cj", m_nj); - Kokkos::parallel_for(m_nj, KOKKOS_LAMBDA(const int& j){ - Vj_cj[j] = Vj[j]/cj[j]; + m_Vj_over_cj[j] = Vj[j]/cj[j]; }); - Kokkos::parallel_reduce(m_nj, ReduceMin(Vj_cj), dt); + double dt = std::numeric_limits<double>::max(); + Kokkos::parallel_reduce(m_nj, ReduceMin(m_Vj_over_cj), dt); return dt; } KOKKOS_INLINE_FUNCTION -const Kokkos::View<const double*> -AcousticSolver::inverse(const Kokkos::View<const double*>& x) const +void +AcousticSolver::inverse(const Kokkos::View<const double*>& x, + Kokkos::View<double*>& inv_x) const { - Kokkos::View<double*> inv_x("inv_x", x.size()); Kokkos::parallel_for(x.size(), KOKKOS_LAMBDA(const int& r) { inv_x(r) = 1./x(r); }); - - return inv_x; } KOKKOS_INLINE_FUNCTION @@ -210,7 +197,15 @@ void AcousticSolver::computeExplicitFluxes(const Kokkos::View<const double*>& xr AcousticSolver::AcousticSolver(const long int& nj) : m_nj(nj), - m_nr(nj+1) + m_nr(nj+1), + m_br("br", m_nr), + m_Ajr("Ajr", m_nj), + m_Ar("Ar", m_nr), + m_inv_Ar("inv_Ar", m_nr), + m_Fjr("Fjr", m_nj), + m_ur("ur", m_nr), + m_rhocj("rho_c", m_nj), + m_Vj_over_cj("Vj_over_cj", m_nj) { Kokkos::View<double*> xj("xj",m_nj); Kokkos::View<double*> rhoj("rhoj",m_nj); @@ -276,20 +271,20 @@ AcousticSolver::AcousticSolver(const long int& nj) }); Kokkos::parallel_for(nj, KOKKOS_LAMBDA(const int& j){ - if (xj[j]<0.5) { - rhoj[j]=1; - } else { - rhoj[j]=0.125; - } - }); + if (xj[j]<0.5) { + rhoj[j]=1; + } else { + rhoj[j]=0.125; + } + }); Kokkos::parallel_for(nj, KOKKOS_LAMBDA(const int& j){ - if (xj[j]<0.5) { - pj[j]=1; - } else { - pj[j]=0.1; - } - }); + if (xj[j]<0.5) { + pj[j]=1; + } else { + pj[j]=0.1; + } + }); Kokkos::parallel_for(nj, KOKKOS_LAMBDA(const int& j){ gammaj[j] = 1.4; @@ -301,13 +296,14 @@ AcousticSolver::AcousticSolver(const long int& nj) Kokkos::parallel_for(nj, KOKKOS_LAMBDA(const int& j){ Ej[j] = ej[j]+0.5*uj[j]*uj[j]; - }); + }); Kokkos::parallel_for(nj, KOKKOS_LAMBDA(const int& j){ - mj[j] = rhoj[j] * Vj[j]; - }); + mj[j] = rhoj[j] * Vj[j]; + }); - const Kokkos::View<const double*> inv_mj=inverse(mj); + Kokkos::View<double*> inv_mj("inv_mj",m_nj); + inverse(mj, inv_mj); const double tmax=0.2; double t=0; @@ -329,14 +325,14 @@ AcousticSolver::AcousticSolver(const long int& nj) dt=tmax-t; t=tmax; } - - Kokkos::View<double*> ur; - Kokkos::View<double*[2]> Fjr; - + computeExplicitFluxes(xr, xj, rhoj, uj, pj, cj, Vj, Cjr, cell_nodes, node_cells, node_nb_cells, node_cell_local_node, - ur, Fjr); + m_ur, m_Fjr); + + const Kokkos::View<const double*[2]> Fjr = m_Fjr; + const Kokkos::View<const double*> ur = m_ur; Kokkos::parallel_for(nj, KOKKOS_LAMBDA(const int& j){ double momentum_fluxes = 0; @@ -359,8 +355,8 @@ AcousticSolver::AcousticSolver(const long int& nj) }); Kokkos::parallel_for(nj, KOKKOS_LAMBDA(const int& j){ - xj[j] = 0.5*(xr[j]+xr[j+1]); - Vj[j] = xr[j+1]-xr[j]; + xj[j] = 0.5*(xr[j]+xr[j+1]); + Vj[j] = xr[j+1]-xr[j]; }); Kokkos::parallel_for(nj, KOKKOS_LAMBDA(const int& j){ diff --git a/experimental/AcousticSolver.hpp b/experimental/AcousticSolver.hpp index d8429d5ac..5ac400a9d 100644 --- a/experimental/AcousticSolver.hpp +++ b/experimental/AcousticSolver.hpp @@ -6,7 +6,6 @@ class AcousticSolver { private: - inline const Kokkos::View<const double*> computeRhoCj(const Kokkos::View<const double*>& rhoj, const Kokkos::View<const double*>& cj); @@ -42,8 +41,8 @@ private: const Kokkos::View<const double*>& pj, const Kokkos::View<const int*[2]>& cell_nodes); - const Kokkos::View<const double*> - inverse(const Kokkos::View<const double*>& x) const; + void inverse(const Kokkos::View<const double*>& x, + Kokkos::View<double*>& inv_x) const; inline double acoustic_dt(const Kokkos::View<const double*>& Vj, const Kokkos::View<const double*>& cj) const; @@ -68,6 +67,16 @@ private: const long int m_nj; const long int m_nr; + Kokkos::View<double*> m_br; + Kokkos::View<double*[2]> m_Ajr; + Kokkos::View<double*> m_Ar; + Kokkos::View<double*> m_inv_Ar; + Kokkos::View<double*[2]> m_Fjr; + Kokkos::View<double*> m_ur; + Kokkos::View<double*> m_rhocj; + Kokkos::View<double*> m_Vj_over_cj; + + public: AcousticSolver(const long int& nj); }; diff --git a/experimental/AcousticSolverTest.cpp b/experimental/AcousticSolverTest.cpp index 5c1e190b8..fd5e9f994 100644 --- a/experimental/AcousticSolverTest.cpp +++ b/experimental/AcousticSolverTest.cpp @@ -1,8 +1,6 @@ #include <AcousticSolverTest.hpp> #include <rang.hpp> -#include <memory> - #include <BlockPerfectGas.hpp> typedef const double my_double; @@ -313,6 +311,18 @@ AcousticSolverTest::AcousticSolverTest(const long int& nj) int itermax=std::numeric_limits<int>::max(); int iteration=0; + Kokkos::View<TinyVector<2>*[2]> Bjr("Cjr",m_nj); + Bjr(0,0)=zero; + Kokkos::parallel_for(nj, KOKKOS_LAMBDA(const int& j) { + Bjr(j,0)[0]=j; + Bjr(j,0)[1]=2; + + Bjr(j,1)=zero; +}); + + for (size_t j=0; j<m_nj; ++j) { + std::cout << Bjr(j,0) << " | " << Bjr(j,1) << '\n'; + } Kokkos::View<double*[2]> Cjr("Cjr",m_nj); Kokkos::parallel_for(nj, KOKKOS_LAMBDA(const int& j) { Cjr(j,0)=-1; @@ -336,7 +346,7 @@ AcousticSolverTest::AcousticSolverTest(const long int& nj) const Kokkos::View<const double*[2]> Fjr = m_Fjr; const Kokkos::View<const double*> ur = m_ur; - Kokkos::parallel_for(nj, KOKKOS_LAMBDA(const int& j){ + Kokkos::parallel_for(nj, KOKKOS_LAMBDA(const int& j) { double momentum_fluxes = 0; double energy_fluxes = 0; for (int R=0; R<2; ++R) { diff --git a/experimental/AcousticSolverTest.hpp b/experimental/AcousticSolverTest.hpp index f12816355..8c6f09d6d 100644 --- a/experimental/AcousticSolverTest.hpp +++ b/experimental/AcousticSolverTest.hpp @@ -3,9 +3,13 @@ #include <Kokkos_Core.hpp> +#include <TinyVector.hpp> + class AcousticSolverTest { private: + constexpr static size_t dimension = 1; + inline const Kokkos::View<const double*> computeRhoCj(const Kokkos::View<const double*>& rhoj, const Kokkos::View<const double*>& cj); -- GitLab