diff --git a/src/scheme/AcousticSolver.hpp b/src/scheme/AcousticSolver.hpp
index aed106a9ba467d230d74627f796650e02232f4e4..d9c5965195256bda75f1c87ab60a4c56e8cc4f64 100644
--- a/src/scheme/AcousticSolver.hpp
+++ b/src/scheme/AcousticSolver.hpp
@@ -81,35 +81,23 @@ private:
   }
 
   KOKKOS_INLINE_FUNCTION
-  const Kokkos::View<const Rdd**>
-  computeAjr(const Kokkos::View<const double*>& rhocj,
-             const Kokkos::View<const double**>& ljr,
-             const Kokkos::View<const Rd**>& njr)
+  void computeAjr(const Kokkos::View<const double*>& rhocj,
+                  const Kokkos::View<const double**>& ljr,
+                  const Kokkos::View<const Rd**>& njr)
   {
     const Kokkos::View<const Rd**> Cjr = m_mesh_data.Cjr();
-
-    Kokkos::parallel_for("old Ajr", m_mesh.numberOfCells(), KOKKOS_LAMBDA(const int& j) {
-        const auto& cell_nodes = m_mesh.connectivity().m_cell_to_node_matrix.rowConst(j);
-        const double& rho_c = rhocj(j);
-        for (size_t r=0; r<cell_nodes.length; ++r) {
-          m_Ajr(j,r) = tensorProduct(rho_c*Cjr(j,r), njr(j,r));
-        }
-      });
-
     Kokkos::parallel_for("new nested Ajr", m_mesh.numberOfCells(), KOKKOS_LAMBDA(const int& j) {
-        const size_t& nb_nodes =m_Ajr_new.numberOfSubValues(j);
+        const size_t& nb_nodes =m_Ajr.numberOfSubValues(j);
         const double& rho_c = rhocj(j);
         for (size_t r=0; r<nb_nodes; ++r) {
-          m_Ajr_new(j,r) = tensorProduct(rho_c*Cjr(j,r), njr(j,r));
+          m_Ajr(j,r) = tensorProduct(rho_c*Cjr(j,r), njr(j,r));
         }
       });
-
-    return m_Ajr;
   }
 
   KOKKOS_INLINE_FUNCTION
   const Kokkos::View<const Rdd*>
-  computeAr(const Kokkos::View<const Rdd**>& Ajr) {
+  computeAr(const NodeByCellData<Rdd>& Ajr) {
     Kokkos::parallel_for(m_mesh.numberOfNodes(), KOKKOS_LAMBDA(const int& r) {
         Rdd sum = zero;
         const auto& node_to_cell = m_connectivity.m_node_to_cell_matrix.rowConst(r);
@@ -127,7 +115,7 @@ private:
 
   KOKKOS_INLINE_FUNCTION
   const Kokkos::View<const Rd*>
-  computeBr(const Kokkos::View<const Rdd**>& Ajr,
+  computeBr(const NodeByCellData<Rdd>& Ajr,
             const Kokkos::View<const Rd**>& Cjr,
             const Kokkos::View<const Rd*>& uj,
             const Kokkos::View<const double*>& pj) {
@@ -178,7 +166,6 @@ private:
 
         Kokkos::parallel_for(symmetry_bc.numberOfNodes(), KOKKOS_LAMBDA(const int& r_number) {
             const int r = symmetry_bc.nodeList()[r_number];
-            // assert(m_connectivity.nodeNbCells()[r] == 1);
 
             m_Ar(r) = P*m_Ar(r)*P + nxn;
             m_br(r) = P*m_br(r);
@@ -202,8 +189,8 @@ private:
     return m_ur;
   }
 
-  Kokkos::View<Rd**>
-  computeFjr(const Kokkos::View<const Rdd**>& Ajr,
+  void
+  computeFjr(const NodeByCellData<Rdd>& Ajr,
              const Kokkos::View<const Rd*>& ur,
              const Kokkos::View<const Rd**>& Cjr,
              const Kokkos::View<const Rd*>& uj,
@@ -216,8 +203,6 @@ private:
           m_Fjr(j,r) = Ajr(j,r)*(uj(j)-ur(cell_nodes(r)))+pj(j)*Cjr(j,r);
         }
       });
-
-    return m_Fjr;
   }
 
   void inverse(const Kokkos::View<const Rdd*>& A,
@@ -240,25 +225,24 @@ private:
                              const Kokkos::View<const double**>& ljr,
                              const Kokkos::View<const Rd**>& njr) {
     const Kokkos::View<const double*> rhocj  = computeRhoCj(rhoj, cj);
-    const Kokkos::View<const Rdd**> Ajr = computeAjr(rhocj, ljr, njr);
+    // const Kokkos::View<const Rdd**> Ajr =
+    computeAjr(rhocj, ljr, njr);
 
-    const Kokkos::View<const Rdd*> Ar = computeAr(Ajr);
-    const Kokkos::View<const Rd*> br = computeBr(Ajr, Cjr, uj, pj);
+    const Kokkos::View<const Rdd*> Ar = computeAr(m_Ajr);
+    const Kokkos::View<const Rd*> br = computeBr(m_Ajr, Cjr, uj, pj);
 
     this->applyBoundaryConditions();
 
     Kokkos::View<Rd*> ur = m_ur;
-    Kokkos::View<Rd**> Fjr = m_Fjr;
     ur  = computeUr(Ar, br);
-    Fjr = computeFjr(Ajr, ur, Cjr, uj, pj);
+    computeFjr(m_Ajr, ur, Cjr, uj, pj);
   }
 
   Kokkos::View<Rd*> m_br;
-  NodeByCellData<Rdd> m_Ajr_new;
-  Kokkos::View<Rdd**> m_Ajr;
+  NodeByCellData<Rdd> m_Ajr;
   Kokkos::View<Rdd*> m_Ar;
   Kokkos::View<Rdd*> m_inv_Ar;
-  Kokkos::View<Rd**> m_Fjr;
+  NodeByCellData<Rd> m_Fjr;
   Kokkos::View<Rd*> m_ur;
   Kokkos::View<double*> m_rhocj;
   Kokkos::View<double*> m_Vj_over_cj;
@@ -272,11 +256,10 @@ public:
       m_connectivity(m_mesh.connectivity()),
       m_boundary_condition_list(bc_list),
       m_br("br", m_mesh.numberOfNodes()),
-      m_Ajr_new(m_connectivity.m_node_id_per_cell_matrix),
-      m_Ajr("Ajr", m_mesh.numberOfCells(), m_connectivity.maxNbNodePerCell()),
+      m_Ajr(m_connectivity.m_node_id_per_cell_matrix),
       m_Ar("Ar", m_mesh.numberOfNodes()),
       m_inv_Ar("inv_Ar", m_mesh.numberOfNodes()),
-      m_Fjr("Fjr", m_mesh.numberOfCells(), m_connectivity.maxNbNodePerCell()),
+      m_Fjr(m_connectivity.m_node_id_per_cell_matrix),
       m_ur("ur", m_mesh.numberOfNodes()),
       m_rhocj("rho_c", m_mesh.numberOfCells()),
       m_Vj_over_cj("Vj_over_cj", m_mesh.numberOfCells())
@@ -328,7 +311,7 @@ public:
 
     computeExplicitFluxes(xr, xj, rhoj, uj, pj, cj, Vj, Cjr, ljr, njr);
 
-    const Kokkos::View<const Rd**> Fjr = m_Fjr;
+    const NodeByCellData<Rd>& Fjr = m_Fjr;
     const Kokkos::View<const Rd*> ur = m_ur;
 
     const Kokkos::View<const double*> inv_mj = unknowns.invMj();