diff --git a/src/main.cpp b/src/main.cpp
index d4e0015ffb173771f5e34360709afe3f84412cd6..c86067c92979b7ec72491ff1ec08cc3ca2bb1847 100644
--- a/src/main.cpp
+++ b/src/main.cpp
@@ -284,7 +284,7 @@ int main(int argc, char *argv[])
// NAVIER-STOKES AVEC SPLITTING
/*
// Etape 1 du splitting - Euler
- double dt_euler = 0.9*acoustic_solver.acoustic_dt(Vj, cj)
+ double dt_euler = 0.9*acoustic_solver.acoustic_dt(Vj, cj);
if (t+dt_euler > tmax) {
dt_euler = tmax-t;
}
@@ -319,14 +319,14 @@ 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);
if (t+dt > tmax) {
dt = tmax-t;
}
no_splitting.computeNextStep(t,dt, unknowns);
t += dt;
-
+
block_eos.updatePandCFromRhoE();
diff --git a/src/mesh/Mesh.hpp b/src/mesh/Mesh.hpp
index 4313b3f73e4d2ef0a37feb0312c041f67c859fd9..007c679443e626e49b3aac3c4b5fa9d9af4a61a0 100644
--- a/src/mesh/Mesh.hpp
+++ b/src/mesh/Mesh.hpp
@@ -75,7 +75,7 @@ public:
// pas constant
-
+ /*
Mesh(const Connectivity& connectivity)
: m_connectivity(connectivity),
m_xr("xr", connectivity.numberOfNodes()),
@@ -91,7 +91,7 @@ public:
m_xr[r][0] = a + r*delta_x;
});
}
-
+ */
// pas non constant
/*
@@ -123,9 +123,9 @@ public:
}
*/
-
+
// pas aleatoire
- /*
+
Mesh(const Connectivity& connectivity)
: m_connectivity(connectivity),
m_xr("xr", connectivity.numberOfNodes()),
@@ -141,24 +141,27 @@ public:
m_xr[connectivity.numberOfNodes()-1][0] = b;
std::random_device rd;
std::mt19937 mt(rd());
- std::uniform_real_distribution<double> dist(-h/2.1,h/2.1);
+ std::uniform_real_distribution<double> dist(-h/2.1,h/2.1);
- // creation fichier pour tracer h en fonction de x
- std::ofstream fout("alea");
-
for (int r=1; r<connectivity.numberOfNodes()-1; ++r){
const double delta_xr = dist(mt);
m_xr[r][0] = r*h + delta_xr;
- fout << m_xr[r][0] << ' ' << m_xr[r][0]-m_xr[r-1][0] << '\n';
- //std::cout << m_xr[r][0] << std::endl;
- }
-
- fout.close();
-
- std::exit(0);
+ }
+
+ // creation fichier pour tracer h en fonction de x
+ //std::ofstream fout("alea");
+ //for (int r=1; r<connectivity.numberOfNodes()-1; ++r){
+ //const double delta_xr = dist(mt);
+ //m_xr[r][0] = r*h + delta_xr;
+ //fout << m_xr[r][0] << ' ' << m_xr[r][0]-m_xr[r-1][0] << '\n';
+ //std::cout << m_xr[r][0] << std::endl;
+ //}
+ //fout.close();
+ //std::exit(0);
+
}
- */
+
~Mesh()
{
diff --git a/src/scheme/FiniteVolumesEulerUnknowns.hpp b/src/scheme/FiniteVolumesEulerUnknowns.hpp
index 3cbf6e542e560e3d9041e3b28b52d934b416ae0d..04a1d1afc13f62c77abb3348c49dadbad442431f 100644
--- a/src/scheme/FiniteVolumesEulerUnknowns.hpp
+++ b/src/scheme/FiniteVolumesEulerUnknowns.hpp
@@ -353,7 +353,7 @@ public:
Kokkos::parallel_for(m_mesh.numberOfCells(), KOKKOS_LAMBDA(const int& j){
// Differents k (xi)
//m_kj[j] = xj[j][0];
- m_kj[j] = 0.5;
+ m_kj[j] = 0.005;
// TAC
@@ -401,15 +401,15 @@ public:
});
Kokkos::parallel_for(m_mesh.numberOfCells(), KOKKOS_LAMBDA(const int& j){
- m_PTj(j) = 0;
+ m_PTj(j) = m_pj(j);
});
// Conditions aux bords de Dirichlet sur u et k
m_uL[0] = zero;
m_uR[0] = zero;
- m_kL[0] = 0.5;
- m_kR[0] = 0.5;
+ m_kL[0] = 0.005;
+ m_kR[0] = 0.005;
Kokkos::parallel_for(m_mesh.numberOfCells(), KOKKOS_LAMBDA(const int& j){
diff --git a/src/scheme/NoSplitting.hpp b/src/scheme/NoSplitting.hpp
index 59ed0203d57d7abcd9a3506491b375c44ea8760e..b971ef440318acca9377d9d3a95a763e2a531fd6 100644
--- a/src/scheme/NoSplitting.hpp
+++ b/src/scheme/NoSplitting.hpp
@@ -314,9 +314,12 @@ public:
Kokkos::View<double*> PTj = unknowns.PTj();
Kokkos::View<double*> kL = unknowns.kL();
Kokkos::View<double*> kR = unknowns.kR();
+ Kokkos::View<Rd*> uL = unknowns.uL();
+ Kokkos::View<Rd*> uR = unknowns.uR();
const Kokkos::View<const Rd*> xj = m_mesh_data.xj();
const Kokkos::View<const double*> Vj = m_mesh_data.Vj();
+ const Kokkos::View<const double*> Vl = m_mesh_data.Vl();
const Kokkos::View<const Rd**> Cjr = m_mesh_data.Cjr();
Kokkos::View<Rd*> xr = m_mesh.xr();
@@ -324,23 +327,6 @@ public:
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));
- }
-
- if (j == 0) {
- PTj(j) = pj(j) - kL(0)*sum/Vj(j);
- } else if (j == m_mesh.numberOfCells()-1) {
- PTj(j) = pj(j) - kR(0)*sum/Vj(j);
- } else {
- PTj(j) = pj(j) - kj(j)*sum/Vj(j);
- }
-
- });
// Calcule les flux
computeExplicitFluxes(xr, xj, kj, rhoj, uj, PTj, cj, Vj, Cjr, t);
@@ -363,7 +349,9 @@ public:
// ajout second membre pour kidder (k cst)
Ej[j] -= (dt*inv_mj[j])*Vj(j)*((kj(j)*t*t)/(((50./9.)-t*t)*((50./9.)-t*t)));
-
+ // ajout second membre pour kidder (k = x)
+ //uj[j][0] += (dt*inv_mj[j])*Vj(j)*(t/((50./9.)-t*t));
+ //Ej[j] -= (dt*inv_mj[j])*Vj(j)*((2.*xj[j][0]*t*t)/(((50./9.)-t*t)*((50./9.)-t*t)));
});
// Calcul de e par la formule e = E-0.5 u^2
@@ -386,11 +374,33 @@ public:
rhoj[j] = mj[j]/Vj[j];
});
+ // Calcul de PT
+ Kokkos::parallel_for(m_mesh.numberOfCells(), KOKKOS_LAMBDA(const int& j) {
+ double sum = 0;
+ double sum1 = 0;
+ for (int m=0; m<cell_nb_nodes(j); ++m) {
+ sum += (uj(cell_nodes(j,m)), Cjr(cell_nodes(j,m), m));
+ sum1 += Vj(cell_nodes(j,m));
+ }
+
+ if (j == 0) {
+ //PTj(j) = pj(j) - kj(j)*(uj[j][0]-uL[0][0])/Vl(0);
+ PTj(j) = pj(j) + kj(j)*(t/((50./9.)-t*t));
+ } else if (j == m_mesh.numberOfCells()-1) {
+ PTj(j) = pj(j) + kj(j)*(t/((50./9.)-t*t));
+ //PTj(j) = pj(j) - kj(j)*(uR[0][0]-uj[j][0])/Vl(m_mesh.numberOfFaces()-1);
+ } else {
+ //PTj(j) = pj(j) - kj(j)*sum/Vl(j);
+ PTj(j) = pj(j) - kj(j)*2.*sum/sum1;
+ }
+
+ });
+
// Mise a jour de k
/*
Kokkos::parallel_for(m_mesh.numberOfCells(), KOKKOS_LAMBDA(const int& j) {
kj(j) = xj[j][0];
-
+
});
*/
}