maillage aleatoire ok, mais le pas n est pas tres equilibre

src/main.cpp

@@ -158,6 +158,11 @@ int main(int argc, char *argv[])
     Kokkos::View<Rd*> uj = unknowns.uj();
     BlockPerfectGas block_eos(rhoj, ej, pj, gammaj, cj);
 
+    Kokkos::View<double*> nuL = unknowns.nuL();
+    Kokkos::View<double*> nuR = unknowns.nuR();
+    Kokkos::View<double*> kL = unknowns.kL();
+    Kokkos::View<double*> kR = unknowns.kR();
+    
     double c = 0.;
     c = finite_volumes_diffusion.conservatif(unknowns);
     
@@ -303,7 +308,7 @@ int main(int argc, char *argv[])
 
       // DIFFUSION PURE
       
-      double dt = 0.9*finite_volumes_diffusion.diffusion_dt(rhoj,kj,nuj,cj);
+      double dt = 0.9*finite_volumes_diffusion.diffusion_dt(rhoj,kj,nuj,cj,nuL,nuR,kL,kR);
       if (t+dt > tmax) {
 	dt = tmax-t;
       }
@@ -540,14 +545,14 @@ int main(int argc, char *argv[])
       // ENTROPY TEST
       //finite_volumes_diffusion.entropie(unknowns);
 
-      
+      /*
       // STOCKAGE COORDONNES ET TEMPS
       for (size_t j=0; j<mesh.numberOfCells(); ++j) {
 	x[i][j] = xj[j][0];
       }
       tempo[i] = t;
       i = i + 1;
-      
+      */
 
     }
 
@@ -567,6 +572,7 @@ int main(int argc, char *argv[])
       fout << ' ' << '\n';
     }
     */
+    /*
     std::ofstream fout("cara1");
     fout.precision(15);
     for (int j=0; j<mesh.numberOfCells(); ++j) {
@@ -580,7 +586,7 @@ int main(int argc, char *argv[])
 	fout << ' ' << '\n';
       }
     }
-    
+   */ 
     
     /*
     double error1 = 0.;

src/mesh/Mesh.hpp

@@ -140,15 +140,16 @@ public:
     m_xmax[0][0] = b;
     const double h = (b-a)/connectivity.numberOfCells();
     m_xr[0][0] = a;
-    m_xr[connectivity.numberOfNodes()-1] = b;
+    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/1.01,h);
     for (int r=1; r<connectivity.numberOfNodes()-1; ++r){
       const double delta_xr = dist(mt);
       m_xr[r][0] = m_xr[r-1][0] + std::abs(delta_xr);
       std::cout << m_xr[r][0] << std::endl;
       } 
+    std::exit(0);
   }
   
 

src/scheme/FiniteVolumesDiffusion.hpp

@@ -117,11 +117,12 @@ private:
     int cell_here = face_cells(0,0);
     int local_face_number_in_cell = face_cell_local_face(0,0);
     m_Fl(0) = -(kL(0) + kj(cell_here))*(1./(2*Vl(0)))*(tensorProduct(uj(cell_here), Cjr(cell_here, local_face_number_in_cell)) - tensorProduct(uL(0), Cjr(cell_here, local_face_number_in_cell)));
+    //m_Fl(0) = -((kL(0)/Vj(0) + kj(cell_here)/Vj(cell_here))/(1./Vj(0) + 1./Vj(cell_here)))*(1./Vl(0))*(tensorProduct(uj(cell_here), Cjr(cell_here, local_face_number_in_cell)) - tensorProduct(uL(0), Cjr(cell_here, local_face_number_in_cell)));
     
     cell_here = face_cells(m_mesh.numberOfFaces()-1,0);
     local_face_number_in_cell = face_cell_local_face(m_mesh.numberOfFaces()-1,0);
     m_Fl(m_mesh.numberOfFaces()-1) = -(kR(0) + kj(cell_here))*(1/(2.*Vl(m_mesh.numberOfFaces()-1)))*(tensorProduct(uj(cell_here), Cjr(cell_here, local_face_number_in_cell)) - tensorProduct(uR(0), Cjr(cell_here, local_face_number_in_cell)));
-    
+    //m_Fl(m_mesh.numberOfFaces()-1) = -((kR(0)/Vj(m_mesh.numberOfCells()-1) + kj(cell_here)/Vj(cell_here))/(1./Vj(cell_here) + 1./Vj(m_mesh.numberOfCells()-1)))*(1./Vl(m_mesh.numberOfFaces()-1))*(tensorProduct(uj(cell_here), Cjr(cell_here, local_face_number_in_cell)) - tensorProduct(uR(0), Cjr(cell_here, local_face_number_in_cell)));
     
     // Kidder
 
@@ -237,10 +238,12 @@ private:
     
     int cell_here = face_cells(0,0);
     m_Bl(0) = (nuL(0) + nuj(cell_here))*(1./(2*Vl(0)))*(Tj(cell_here) - TL(0));
+    //m_Bl(0) = (nuL(0)/Vj(0) + nuj(cell_here)/Vj(cell_here))/(1./Vj(0) + 1./Vj(cell_here))*((Tj(cell_here)-TL(0))/Vl(0)); 
+
     
     cell_here = face_cells(m_mesh.numberOfFaces()-1,0);
     m_Bl(m_mesh.numberOfFaces()-1) = -(nuR(0) + nuj(cell_here))*(1/(2.*Vl(m_mesh.numberOfFaces()-1)))*(Tj(cell_here) - TR(0));
-    
+    //m_Bl(m_mesh.numberOfFaces()-1) = (nuR(0)/Vj(m_mesh.numberOfCells()-1) + nuj(cell_here)/Vj(cell_here))/(1./Vj(m_mesh.numberOfCells()-1) + 1./Vj(cell_here))*((Tj(cell_here)-TR(0))/Vl(m_mesh.numberOfFaces()-1)); 
 
     // Kiddder
     /*
@@ -331,7 +334,12 @@ public:
   double diffusion_dt(const Kokkos::View<const double*>& rhoj,
 		      const Kokkos::View<const double*>& kj,
 		      const Kokkos::View<const double*>& nuj,
-		      const Kokkos::View<const double*>& cj) const {
+		      const Kokkos::View<const double*>& cj,
+		      const Kokkos::View<const double*>& nuL,
+		      const Kokkos::View<const double*>& nuR,
+		      const Kokkos::View<const double*>& kL,
+		      const Kokkos::View<const double*>& kR
+		      ) const {
 
     Kokkos::View<double*> dt_j("dt_j", m_mesh.numberOfCells());
 
@@ -350,6 +358,7 @@ public:
       = m_connectivity.cellNbNodes();
 
     Kokkos::parallel_for(m_mesh.numberOfCells(), KOKKOS_LAMBDA(const int& j){
+	
 	double minVl = std::numeric_limits<double>::max();
 	for (int ll=0; ll<cell_nb_faces(j); ++ll) {
 	  minVl = std::min(minVl, Vl(cell_faces(j, ll)));
@@ -358,10 +367,32 @@ public:
 	double sum = 0.;
 	double sum1 = 0.;
 
+	double vol = 0.;
+	
+	for (int m = 0; m < cell_nb_nodes(j); ++m) {
+	   if (cell_nodes(j,m) == 0) {
+	    sum += kL(0)/Vj(0);
+	    sum1 += nuL(0)/Vj(0);
+	    vol += 1./Vj(0);
+	  } else if (cell_nodes(j,m) == m_mesh.numberOfCells()) {
+	    sum += kR(0)/Vj(m_mesh.numberOfCells()-1);
+	    sum1 = nuR(0)/Vj(m_mesh.numberOfCells()-1);
+	    vol += 1./Vj(m_mesh.numberOfCells()-1);
+	    } else {
+	    sum += kj(cell_nodes(j,m))/Vj(cell_nodes(j,m));
+	    sum1 += nuj(cell_nodes(j,m))/Vj(cell_nodes(j,m));
+	    vol += 1./Vj(cell_nodes(j,m));
+	    }
+	}	
+	sum = sum/vol;
+	sum1 = sum1/vol;
+      
+	/*
 	for (int m = 0; m < cell_nb_nodes(j); ++m) {
 	  sum += kj(cell_nodes(j,m));
 	  sum1 += nuj(cell_nodes(j,m));
 	}
+	*/
 
 	if (sum > sum1) {
 	  if (sum == 0.) {