diff --git a/src/VDF/Operateurs/Op_Diff_Dift/Op_Diff_Dift_base/Op_Dift_VDF_Elem_base.cpp b/src/VDF/Operateurs/Op_Diff_Dift/Op_Diff_Dift_base/Op_Dift_VDF_Elem_base.cpp
index c15b7e1dca..af945f2b39 100644
--- a/src/VDF/Operateurs/Op_Diff_Dift/Op_Diff_Dift_base/Op_Dift_VDF_Elem_base.cpp
+++ b/src/VDF/Operateurs/Op_Diff_Dift/Op_Diff_Dift_base/Op_Dift_VDF_Elem_base.cpp
@@ -32,41 +32,88 @@ double Op_Dift_VDF_Elem_base::calculer_dt_stab_elem() const
 {
   double dt_stab, coef = -1.e10;
   const Domaine_VDF& domaine_VDF = iter_->domaine();
-  const IntTab& elem_faces = domaine_VDF.elem_faces();
-  const DoubleVect& alpha_t = diffusivite_turbulente().valeurs();
+  const IntTab& tab_elem_faces = domaine_VDF.elem_faces();
+  [[maybe_unused]] CIntTabView elem_faces = tab_elem_faces.view_ro();
+  const DoubleVect& tab_alpha_t = diffusivite_turbulente().valeurs();
+  [[maybe_unused]] CDoubleArrView alpha_t = tab_alpha_t.view_ro();
   bool is_concentration = (equation().que_suis_je().debute_par("Convection_Diffusion_Concentration") || equation().que_suis_je().debute_par("Convection_Diffusion_Espece"));
 
-  ArrOfInt numfa(2*dimension);
-  ToDo_Kokkos("critical");
-  for (int elem = 0; elem < domaine_VDF.nb_elem(); elem++)
-    {
-      // choix du facteur
-      double rcp = 1.;
-      if (!is_concentration)
-        {
-          const int Ccp = sub_type(Champ_Uniforme, mon_equation->milieu().capacite_calorifique());
-          const int Cr = sub_type(Champ_Uniforme, mon_equation->milieu().masse_volumique());
-          const DoubleTab& tab_Cp = mon_equation->milieu().capacite_calorifique().valeurs();
-          const DoubleTab& tab_r = mon_equation->milieu().masse_volumique().valeurs();
-          rcp = tab_r(Cr ? 0 : elem, 0) * tab_Cp(Ccp ? 0 : elem, 0);
-        }
+  int nb_elem = domaine_VDF.nb_elem(); 
+  const int l_dim  = Objet_U::dimension;
 
-      double moy = 0.;
-      for (int i = 0; i < 2 * dimension; i++) numfa[i] = elem_faces(elem, i);
+  Kokkos::View<double*> alpha_view("alpha_view", nb_elem);
+  auto tab_alpha_view = Kokkos::create_mirror_view(alpha_view);
+  for (int elem = 0; elem < nb_elem; elem++)
+  {
+      tab_alpha_view(elem) = alpha_(elem);
+  }
+  Kokkos::deep_copy(alpha_view, tab_alpha_view);
 
-      // XXX : E Saikali j'ai corrige pour multi inco parce que c'etait 1/dx et pas 1/dx^2 ... donc attention si ecart !
-      // c'etait comme ca : for (int d = 0; d < dimension; d++) moy += 1. / (domaine_VDF.dist_face(numfa[d], numfa[dimension + d], d));
-      for (int d = 0; d < dimension; d++)
+  const DoubleTab& tab_xv = domaine_VDF.xv();
+  CDoubleTabView xv = tab_xv.view_ro();
+
+  // choix du facteur
+  static size_t run = 0;
+  printf("DEBUG: run = %lu\n", run);
+  if (!is_concentration)
+    {
+      const int Ccp = sub_type(Champ_Uniforme, mon_equation->milieu().capacite_calorifique());
+      const int Cr = sub_type(Champ_Uniforme, mon_equation->milieu().masse_volumique());
+      const DoubleTab& tab_Cp = mon_equation->milieu().capacite_calorifique().valeurs();
+      CDoubleTabView Cp = tab_Cp.view_ro();
+      const DoubleTab& tab_r = mon_equation->milieu().masse_volumique().valeurs();
+      CDoubleTabView r = tab_r.view_ro();
+      Kokkos::parallel_reduce(start_gpu_timer(__KERNEL_NAME__), range_1D(0, nb_elem), KOKKOS_LAMBDA (const int& elem, double& l_coef)
+        {
+          double rcp = r(Cr ? 0 : elem, 0) * Cp(Ccp ? 0 : elem, 0);
+          if (elem == nb_elem-1) Kokkos::printf("DEBUG: thread = %d, rcp = %f\n", elem, rcp);
+
+          double moy = 0.;
+          for (int d = 0; d < l_dim; d++)
+            {
+              const double hd = xv(tab_elem_faces(elem, l_dim + d), d) - xv(tab_elem_faces(elem, d), d);
+              moy += 1. / (hd * hd);
+            }
+          if (elem == nb_elem-1) Kokkos::printf("DEBUG: thread = %d, moy = %f\n", elem, moy);
+
+          const double alpha_local = (alpha_view(elem) + alpha_t(elem)) / rcp * moy;
+          if (elem == nb_elem-1) Kokkos::printf("DEBUG: thread = %d, alpha_local = %f\n",elem, alpha_local);
+
+          if (elem == nb_elem-1) Kokkos::printf("DEBUG: thread = %d, l_coef = %f\n", elem, l_coef);
+          l_coef = Kokkos::fmax(l_coef, alpha_local);
+          if (elem == nb_elem-1) Kokkos::printf("DEBUG: thread = %d, l_coef = %f\n", elem, l_coef);
+        }, Kokkos::Max<double>(coef));
+        end_gpu_timer(__KERNEL_NAME__);
+    }
+  else
+    {
+      for (int elem = 0; elem < nb_elem; elem++)
         {
-          const double hd = domaine_VDF.dist_face(numfa[d], numfa[dimension + d], d);
-          moy += 1. / (hd * hd);
+          double rcp = 1.;
+          if (elem == nb_elem-1) printf("DEBUG: thread = %d, rcp = %f\n", elem, rcp);
+
+          double moy = 0.;
+          for (int d = 0; d < l_dim; d++)
+            {
+              const double hd = xv(tab_elem_faces(elem, l_dim + d), d) - xv(tab_elem_faces(elem, d), d);
+              moy += 1. / (hd * hd);
+            }
+          if (elem == nb_elem-1) printf("DEBUG: thread = %d, moy = %f\n", elem, moy);
+
+          const double alpha_local = (tab_alpha_view(elem) + tab_alpha_t(elem)) / rcp * moy;
+          if (elem == nb_elem-1) printf("DEBUG: thread = %d, alpha_local = %f\n",elem, alpha_local);
+          if (elem == nb_elem-1) printf("DEBUG: thread = %d, coef = %f\n", elem, coef);
+          coef = std::max(coef, alpha_local);
+          if (elem == nb_elem-1) printf("DEBUG: thread = %d, coef = %f\n", elem, coef);
         }
-      const double alpha_local = (alpha_(elem) + alpha_t(elem)) / rcp * moy;
-      coef = std::max(coef, alpha_local);
     }
+  run++;
 
+  printf("DEBUG: coef = %f\n", coef);
   coef = Process::mp_max(coef);
+  printf("DEBUG: coef = %f\n", coef);
   dt_stab = 1. / (2. * (coef + DMINFLOAT));
+  printf("DEBUG: dt_stab = %f\n", dt_stab);
   return dt_stab;
 }