diff --git a/src/main.cpp b/src/main.cpp index 90feecfff0bcef20d897cf777deabf334f644b61..00086299a545661dd3c9ace26fbc7563cbee055c 100644 --- a/src/main.cpp +++ b/src/main.cpp @@ -144,7 +144,7 @@ int main(int argc, char *argv[]) const Kokkos::View<const double*> Vj = mesh_data.Vj(); const Kokkos::View<const Rd**> Cjr = mesh_data.Cjr(); - const double tmax=0.8; + const double tmax=0.2; double t=0.; int itermax=std::numeric_limits<int>::max(); @@ -320,7 +320,7 @@ int main(int argc, char *argv[]) // NAVIER-STOKES SANS SPLITTING - double dt = 0.9*finite_volumes_diffusion.diffusion_dt(rhoj, kj,nuj, cj, nuL, nuR, kL,kR); + double dt = 0.05*finite_volumes_diffusion.diffusion_dt(rhoj, kj,nuj, cj, nuL, nuR, kL,kR); if (t+dt > tmax) { dt = tmax-t; } diff --git a/src/scheme/FiniteVolumesEulerUnknowns.hpp b/src/scheme/FiniteVolumesEulerUnknowns.hpp index 2b78943cdf57d9719d22932f528a58ea04d0e299..81b26c05ffc635c34391b3a8b51a12a04a2d9110 100644 --- a/src/scheme/FiniteVolumesEulerUnknowns.hpp +++ b/src/scheme/FiniteVolumesEulerUnknowns.hpp @@ -12,6 +12,7 @@ public: typedef TinyVector<dimension> Rd; private: + const MeshDataType& m_mesh_data; const MeshType& m_mesh; @@ -26,6 +27,7 @@ private: Kokkos::View<double*> m_mj; Kokkos::View<double*> m_inv_mj; Kokkos::View<double*> m_kj; + Kokkos::View<double*> m_PTj; Kokkos::View<Rd*> m_uL; Kokkos::View<Rd*> m_uR; @@ -143,7 +145,17 @@ public: { return m_kj; } + + Kokkos::View<double*> PTj() + { + return m_PTj; + } + const Kokkos::View<const double*> PTj() const + { + return m_PTj; + } + Kokkos::View<Rd*> uL() { return m_uL; @@ -283,6 +295,8 @@ public: { const Kokkos::View<const Rd*> xj = m_mesh_data.xj(); + const Kokkos::View<const double*> Vj = m_mesh_data.Vj(); + Kokkos::parallel_for(m_mesh.numberOfCells(), KOKKOS_LAMBDA(const int& j){ // TAC /* @@ -328,7 +342,6 @@ public: m_Ej[j] = m_ej[j]+0.5*(m_uj[j],m_uj[j]); }); - const Kokkos::View<const double*> Vj = m_mesh_data.Vj(); Kokkos::parallel_for(m_mesh.numberOfCells(), KOKKOS_LAMBDA(const int& j){ m_mj[j] = m_rhoj[j] * Vj[j]; }); @@ -386,7 +399,11 @@ public: */ }); - + + Kokkos::parallel_for(m_mesh.numberOfCells(), KOKKOS_LAMBDA(const int& j){ + m_PTj(j) = 0; + }); + // Conditions aux bords de Dirichlet sur u et k m_uL[0] = zero; @@ -524,7 +541,8 @@ public: m_nuR("nuR", 1), m_entropy("entropy", m_mesh.numberOfCells()), m_entropy_new("entropy_new", m_mesh.numberOfCells()), - m_S0("S0", m_mesh.numberOfCells()) + m_S0("S0", m_mesh.numberOfCells()), + m_PTj("PTj",m_mesh.numberOfCells()) { ; } diff --git a/src/scheme/NoSplitting.hpp b/src/scheme/NoSplitting.hpp index f046a3d8664ae65fd7d5c437c20537a6032db128..039b535f2648cafb28574c971fabbf1edd82c2f1 100644 --- a/src/scheme/NoSplitting.hpp +++ b/src/scheme/NoSplitting.hpp @@ -71,7 +71,7 @@ private: } }; - + /* // ---- KOKKOS_INLINE_FUNCTION @@ -100,7 +100,7 @@ private: } // ---- - + */ KOKKOS_INLINE_FUNCTION const Kokkos::View<const double*> @@ -154,7 +154,7 @@ private: computeBr(const Kokkos::View<const Rdd**>& Ajr, const Kokkos::View<const Rd**>& Cjr, const Kokkos::View<const Rd*>& uj, - const Kokkos::View<const double*>& pj, + const Kokkos::View<const double*>& PTj, const double t) { const Kokkos::View<const unsigned int**>& node_cells = m_connectivity.nodeCells(); const Kokkos::View<const unsigned short**>& node_cell_local_node = m_connectivity.nodeCellLocalNode(); @@ -167,7 +167,7 @@ private: 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); - br += Ajr(J,R)*uj(J) + pj(J)*Cjr(J,R); + br += Ajr(J,R)*uj(J) + PTj(J)*Cjr(J,R); } }); @@ -212,14 +212,14 @@ private: const Kokkos::View<const Rd*>& ur, const Kokkos::View<const Rd**>& Cjr, const Kokkos::View<const Rd*>& uj, - const Kokkos::View<const double*>& pj) { + const Kokkos::View<const double*>& PTj) { const Kokkos::View<const unsigned int**>& cell_nodes = m_connectivity.cellNodes(); 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_Fjr(j,r) = Ajr(j,r)*(uj(j)-ur(cell_nodes(j,r)))+pj(j)*Cjr(j,r); + m_Fjr(j,r) = Ajr(j,r)*(uj(j)-ur(cell_nodes(j,r)))+PTj(j)*Cjr(j,r); } }); @@ -251,7 +251,7 @@ private: const Kokkos::View<const double*>& kj, const Kokkos::View<const double*>& rhoj, const Kokkos::View<const Rd*>& uj, - const Kokkos::View<const double*>& pj, + const Kokkos::View<const double*>& PTj, const Kokkos::View<const double*>& cj, const Kokkos::View<const double*>& Vj, const Kokkos::View<const Rd**>& Cjr, @@ -259,15 +259,13 @@ private: const Kokkos::View<const double*> rhocj = computeRhoCj(rhoj, cj); const Kokkos::View<const Rdd**> Ajr = computeAjr(rhocj, Cjr); - const Kokkos::View<const double*> PTj = computePj(pj, kj, uj); - const Kokkos::View<const Rdd*> Ar = computeAr(Ajr); const Kokkos::View<const Rd*> br = computeBr(Ajr, Cjr, uj, PTj, t); Kokkos::View<Rd*> ur = m_ur; Kokkos::View<Rd**> Fjr = m_Fjr; ur = computeUr(Ar, br, t); - Fjr = computeFjr(Ajr, ur, Cjr, uj, pj); + Fjr = computeFjr(Ajr, ur, Cjr, uj, PTj); } Kokkos::View<Rd*> m_br; @@ -313,21 +311,32 @@ public: Kokkos::View<double*> cj = unknowns.cj(); Kokkos::View<double*> kj = unknowns.kj(); Kokkos::View<double*> nuj = unknowns.nuj(); + Kokkos::View<double*> PTj = unknowns.PTj(); 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(); Kokkos::View<Rd*> xr = m_mesh.xr(); + const Kokkos::View<const unsigned int**>& cell_nodes = m_connectivity.cellNodes(); + + const Kokkos::View<const unsigned short*> cell_nb_nodes + = m_connectivity.cellNbNodes(); + + Kokkos::parallel_for(m_mesh.numberOfCells(), KOKKOS_LAMBDA(const int& j) { + double sum = 0; + for (int m=0; m<cell_nb_nodes(j); ++m) { + //sum += uj[cell_nodes(j,m)][0]; + sum -= (uj(cell_nodes(j,m)), Cjr(cell_nodes(j,m), m)); + } + PTj(j) = pj(j) - kj(j)*sum/Vj(j); + }); + // Calcule les flux - computeExplicitFluxes(xr, xj, kj, rhoj, uj, pj, cj, Vj, Cjr, t); + computeExplicitFluxes(xr, xj, kj, rhoj, uj, PTj, cj, Vj, Cjr, t); const Kokkos::View<const Rd**> Fjr = m_Fjr; const Kokkos::View<const Rd*> ur = m_ur; - const Kokkos::View<const unsigned short*> cell_nb_nodes - = m_connectivity.cellNbNodes(); - const Kokkos::View<const unsigned int**>& cell_nodes - = m_connectivity.cellNodes(); // Mise a jour de la vitesse et de l'energie totale specifique const Kokkos::View<const double*> inv_mj = unknowns.invMj(); @@ -364,12 +373,12 @@ public: }); // Mise a jour de k - /* + Kokkos::parallel_for(m_mesh.numberOfCells(), KOKKOS_LAMBDA(const int& j) { kj(j) = xj[j][0]; }); - */ + } };