diff --git a/experimental/AcousticSolver.cpp b/experimental/AcousticSolver.cpp
index 3088ff7176c03d46e6e429cb0b503123d0ecd66f..80ed8310e81e6c19e682d889b344953f3de418bf 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 d8429d5acc023acea9ef6930c8de1306f0bc1370..5ac400a9db22b4a96bcbabdf90fdd26d53eaa2ba 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 5c1e190b8f06163cb952db133523067d9b142bca..fd5e9f994ffcc6c9cbfb65614632beed246c2d4f 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 f1281635550cb9323452c5d4cb85e688b6216704..8c6f09d6d26f0226b9c72f6366502e0e500e5b64 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);