diff --git a/opm/autodiff/BlackoilAmg.hpp b/opm/autodiff/BlackoilAmg.hpp
index dd0ca952e..083dfbf94 100644
--- a/opm/autodiff/BlackoilAmg.hpp
+++ b/opm/autodiff/BlackoilAmg.hpp
@@ -19,6 +19,7 @@
 #ifndef OPM_AMG_HEADER_INCLUDED
 #define OPM_AMG_HEADER_INCLUDED
 
+#include <ewoms/linear/matrixblock.hh>
 #include <opm/autodiff/ParallelOverlappingILU0.hpp>
 #include <opm/autodiff/CPRPreconditioner.hpp>
 #include <dune/istl/paamg/twolevelmethod.hh>
@@ -168,6 +169,12 @@ struct ScalarType<Dune::MatrixBlock<FieldType, ROWS, COLS> >
     typedef Dune::MatrixBlock<FieldType, 1, 1> value;
 };
 
+template<typename FieldType, int ROWS, int COLS>
+struct ScalarType<Ewoms::MatrixBlock<FieldType, ROWS, COLS> >
+{
+    typedef Ewoms::MatrixBlock<FieldType, 1, 1> value;
+};
+
 template<typename BlockType, typename Allocator>
 struct ScalarType<Dune::BCRSMatrix<BlockType, Allocator> >
 {
diff --git a/opm/autodiff/BlackoilModelEbos.hpp b/opm/autodiff/BlackoilModelEbos.hpp
index 649237f64..9b9de2c86 100644
--- a/opm/autodiff/BlackoilModelEbos.hpp
+++ b/opm/autodiff/BlackoilModelEbos.hpp
@@ -109,6 +109,7 @@ namespace Opm {
         typedef typename GET_PROP_TYPE(TypeTag, Simulator)         Simulator;
         typedef typename GET_PROP_TYPE(TypeTag, Grid)              Grid;
         typedef typename GET_PROP_TYPE(TypeTag, ElementContext)    ElementContext;
+        typedef typename GET_PROP_TYPE(TypeTag, SparseMatrixAdapter) SparseMatrixAdapter;
         typedef typename GET_PROP_TYPE(TypeTag, SolutionVector)    SolutionVector ;
         typedef typename GET_PROP_TYPE(TypeTag, PrimaryVariables)  PrimaryVariables ;
         typedef typename GET_PROP_TYPE(TypeTag, FluidSystem)       FluidSystem;
@@ -126,8 +127,8 @@ namespace Opm {
         static const int temperatureIdx = Indices::temperatureIdx;
 
         typedef Dune::FieldVector<Scalar, numEq >        VectorBlockType;
-        typedef Dune::FieldMatrix<Scalar, numEq, numEq >        MatrixBlockType;
-        typedef Dune::BCRSMatrix <MatrixBlockType>      Mat;
+        typedef typename SparseMatrixAdapter::MatrixBlock MatrixBlockType;
+        typedef typename SparseMatrixAdapter::IstlMatrix Mat;
         typedef Dune::BlockVector<VectorBlockType>      BVector;
 
         typedef ISTLSolverEbos<TypeTag> ISTLSolverType;
@@ -358,13 +359,13 @@ namespace Opm {
             ebosSimulator_.model().linearizer().linearize();
             ebosSimulator_.problem().endIteration();
 
-            auto& ebosJac = ebosSimulator_.model().linearizer().matrix();
+            auto& ebosJac = ebosSimulator_.model().linearizer().jacobian();
             if (param_.matrix_add_well_contributions_) {
-                wellModel().addWellContributions(ebosJac);
+                wellModel().addWellContributions(ebosJac.istlMatrix());
             }
             if ( param_.preconditioner_add_well_contributions_ &&
                                   ! param_.matrix_add_well_contributions_ ) {
-                matrix_for_preconditioner_ .reset(new Mat(ebosJac));
+                matrix_for_preconditioner_ .reset(new Mat(ebosJac.istlMatrix()));
                 wellModel().addWellContributions(*matrix_for_preconditioner_);
             }
 
@@ -494,7 +495,7 @@ namespace Opm {
             // r = [r, r_well], where r is the residual and r_well the well residual.
             // r -= B^T * D^-1 r_well
 
-            auto& ebosJac = ebosSimulator_.model().linearizer().matrix();
+            auto& ebosJac = ebosSimulator_.model().linearizer().jacobian();
             auto& ebosResid = ebosSimulator_.model().linearizer().residual();
 
             wellModel().apply(ebosResid);
@@ -502,13 +503,13 @@ namespace Opm {
             // set initial guess
             x = 0.0;
 
-            const Mat& actual_mat_for_prec = matrix_for_preconditioner_ ? *matrix_for_preconditioner_.get() : ebosJac;
+            const Mat& actual_mat_for_prec = matrix_for_preconditioner_ ? *matrix_for_preconditioner_.get() : ebosJac.istlMatrix();
             // Solve system.
             if( isParallel() )
             {	      
                 typedef WellModelMatrixAdapter< Mat, BVector, BVector, BlackoilWellModel<TypeTag>, true > Operator;
 
-                auto ebosJacIgnoreOverlap = Mat(ebosJac);
+                auto ebosJacIgnoreOverlap = Mat(ebosJac.istlMatrix());
                 //remove ghost rows in local matrix
                 makeOverlapRowsInvalid(ebosJacIgnoreOverlap);
 
@@ -522,7 +523,7 @@ namespace Opm {
             else
             {
                 typedef WellModelMatrixAdapter< Mat, BVector, BVector, BlackoilWellModel<TypeTag>, false > Operator;
-                Operator opA(ebosJac, actual_mat_for_prec, wellModel());
+                Operator opA(ebosJac.istlMatrix(), actual_mat_for_prec, wellModel());
                 istlSolver().solve( opA, x, ebosResid );
             }
         }
diff --git a/opm/autodiff/BlackoilWellModel.hpp b/opm/autodiff/BlackoilWellModel.hpp
index d36b30923..1d6a0c0d7 100644
--- a/opm/autodiff/BlackoilWellModel.hpp
+++ b/opm/autodiff/BlackoilWellModel.hpp
@@ -86,7 +86,7 @@ namespace Opm {
             typedef typename GET_PROP_TYPE(TypeTag, Scalar)              Scalar;
             typedef typename GET_PROP_TYPE(TypeTag, RateVector)          RateVector;
             typedef typename GET_PROP_TYPE(TypeTag, GlobalEqVector)      GlobalEqVector;
-            typedef typename GET_PROP_TYPE(TypeTag, JacobianMatrix)      JacobianMatrix;
+            typedef typename GET_PROP_TYPE(TypeTag, SparseMatrixAdapter) SparseMatrixAdapter;
 
             typedef typename Ewoms::BaseAuxiliaryModule<TypeTag>::NeighborSet NeighborSet;
 
@@ -105,7 +105,7 @@ namespace Opm {
 #else
             typedef Dune::FieldMatrix<Scalar, numEq, numEq > MatrixBlockType;
 #endif
-            typedef Dune::BCRSMatrix <MatrixBlockType> Mat;
+            typedef typename SparseMatrixAdapter::IstlMatrix Mat;
 
             typedef Ewoms::BlackOilPolymerModule<TypeTag> PolymerModule;
 
@@ -129,7 +129,7 @@ namespace Opm {
             void applyInitial()
             {}
 
-            void linearize(JacobianMatrix& mat , GlobalEqVector& res);
+            void linearize(SparseMatrixAdapter& mat , GlobalEqVector& res);
 
             void postSolve(GlobalEqVector& deltaX)
             {
diff --git a/opm/autodiff/BlackoilWellModel_impl.hpp b/opm/autodiff/BlackoilWellModel_impl.hpp
index 3ea869fcc..5882ff8ca 100644
--- a/opm/autodiff/BlackoilWellModel_impl.hpp
+++ b/opm/autodiff/BlackoilWellModel_impl.hpp
@@ -94,7 +94,7 @@ namespace Opm {
     template<typename TypeTag>
     void
     BlackoilWellModel<TypeTag>::
-    linearize(JacobianMatrix& mat , GlobalEqVector& res)
+    linearize(SparseMatrixAdapter& mat , GlobalEqVector& res)
     {
         if (!localWellsActive())
             return;
diff --git a/opm/autodiff/ISTLSolverEbos.hpp b/opm/autodiff/ISTLSolverEbos.hpp
index 8d1aca9b5..050c274ce 100644
--- a/opm/autodiff/ISTLSolverEbos.hpp
+++ b/opm/autodiff/ISTLSolverEbos.hpp
@@ -71,10 +71,12 @@ namespace Opm
     class ISTLSolverEbos : public NewtonIterationBlackoilInterface
     {
         typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
-        typedef typename GET_PROP_TYPE(TypeTag, JacobianMatrix) Matrix;
+        typedef typename GET_PROP_TYPE(TypeTag, SparseMatrixAdapter) SparseMatrixAdapter;
         typedef typename GET_PROP_TYPE(TypeTag, GlobalEqVector) Vector;
         typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
 
+        typedef typename SparseMatrixAdapter::IstlMatrix Matrix;
+
         enum { pressureIndex = Indices::pressureSwitchIdx };
 
     public:
diff --git a/opm/autodiff/NewtonIterationBlackoilInterleaved.cpp b/opm/autodiff/NewtonIterationBlackoilInterleaved.cpp
index afaa217f2..2c22bf6fa 100644
--- a/opm/autodiff/NewtonIterationBlackoilInterleaved.cpp
+++ b/opm/autodiff/NewtonIterationBlackoilInterleaved.cpp
@@ -26,13 +26,13 @@
 // Define making clear that the simulator supports AMG
 #define FLOW_SUPPORT_AMG !defined(HAVE_UMFPACK)
 
-#include <opm/autodiff/DuneMatrix.hpp>
 #include <opm/autodiff/CPRPreconditioner.hpp>
 #include <opm/autodiff/NewtonIterationBlackoilInterleaved.hpp>
 #include <opm/autodiff/NewtonIterationUtilities.hpp>
 #include <opm/autodiff/ParallelRestrictedAdditiveSchwarz.hpp>
 #include <opm/autodiff/ParallelOverlappingILU0.hpp>
 #include <opm/autodiff/AutoDiffHelpers.hpp>
+#include <opm/autodiff/DuneMatrix.hpp>
 #include <opm/common/Exceptions.hpp>
 #include <opm/core/linalg/ParallelIstlInformation.hpp>
 
diff --git a/opm/autodiff/WellConnectionAuxiliaryModule.hpp b/opm/autodiff/WellConnectionAuxiliaryModule.hpp
index 702ecb299..611c65303 100644
--- a/opm/autodiff/WellConnectionAuxiliaryModule.hpp
+++ b/opm/autodiff/WellConnectionAuxiliaryModule.hpp
@@ -36,7 +36,7 @@ class WellConnectionAuxiliaryModule
     : public Ewoms::BaseAuxiliaryModule<TypeTag>
 {
     typedef typename GET_PROP_TYPE(TypeTag, GlobalEqVector) GlobalEqVector;
-    typedef typename GET_PROP_TYPE(TypeTag, JacobianMatrix) JacobianMatrix;
+    typedef typename GET_PROP_TYPE(TypeTag, SparseMatrixAdapter) SparseMatrixAdapter;
 
 public:
 
@@ -116,7 +116,7 @@ public:
     void applyInitial()
     {}
 
-    void linearize(JacobianMatrix& , GlobalEqVector&)
+    void linearize(SparseMatrixAdapter& , GlobalEqVector&)
     {
         // Linearization is done in StandardDenseWells
     }
diff --git a/opm/autodiff/WellInterface.hpp b/opm/autodiff/WellInterface.hpp
index e2b59e4cf..673b80975 100644
--- a/opm/autodiff/WellInterface.hpp
+++ b/opm/autodiff/WellInterface.hpp
@@ -80,6 +80,7 @@ namespace Opm
         typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
         typedef typename GET_PROP_TYPE(TypeTag, IntensiveQuantities) IntensiveQuantities;
         typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
+        typedef typename GET_PROP_TYPE(TypeTag, SparseMatrixAdapter) SparseMatrixAdapter;
         typedef typename GET_PROP_TYPE(TypeTag, RateVector) RateVector;
 
         static const int numEq = Indices::numEq;
@@ -87,7 +88,7 @@ namespace Opm
 
         typedef Dune::FieldVector<Scalar, numEq    > VectorBlockType;
         typedef Dune::FieldMatrix<Scalar, numEq, numEq > MatrixBlockType;
-        typedef Dune::BCRSMatrix <MatrixBlockType> Mat;
+        typedef typename SparseMatrixAdapter::IstlMatrix Mat;
         typedef Dune::BlockVector<VectorBlockType> BVector;
         typedef DenseAd::Evaluation<double, /*size=*/numEq> Eval;