diff --git a/src/TNL/Matrices/CSR.h b/src/TNL/Matrices/CSR.h
index e42fca6b73990d05e1d5312d61b7022ecc6f2b49..5ad4327a6310ed1bfaa7993771a4387e27c9dae5 100644
--- a/src/TNL/Matrices/CSR.h
+++ b/src/TNL/Matrices/CSR.h
@@ -41,6 +41,7 @@ class CSR : public Sparse< Real, Device, Index >
    typedef CSR< Real, Devices::Cuda, Index > CudaType;
    typedef Sparse< Real, Device, Index > BaseType;
    typedef typename BaseType::MatrixRow MatrixRow;
+   typedef SparseRow< const RealType, const IndexType > ConstMatrixRow;
 
 
    enum SPMVCudaKernel { scalar, vector, hybrid };
@@ -125,7 +126,7 @@ class CSR : public Sparse< Real, Device, Index >
    MatrixRow getRow( const IndexType rowIndex );
 
    __cuda_callable__
-   const MatrixRow getRow( const IndexType rowIndex ) const;
+   ConstMatrixRow getRow( const IndexType rowIndex ) const;
 
    template< typename Vector >
    __cuda_callable__
diff --git a/src/TNL/Matrices/CSR_impl.h b/src/TNL/Matrices/CSR_impl.h
index 0b100505394fb92c5cad027a80bf7988046c64fc..e021dc52c466788bbc39b6552e71420fb128be6d 100644
--- a/src/TNL/Matrices/CSR_impl.h
+++ b/src/TNL/Matrices/CSR_impl.h
@@ -406,16 +406,16 @@ template< typename Real,
           typename Device,
           typename Index >
 __cuda_callable__
-const typename CSR< Real, Device, Index >::MatrixRow
+typename CSR< Real, Device, Index >::ConstMatrixRow
 CSR< Real, Device, Index >::
 getRow( const IndexType rowIndex ) const
 {
    const IndexType rowOffset = this->rowPointers[ rowIndex ];
    const IndexType rowLength = this->rowPointers[ rowIndex + 1 ] - rowOffset;
-   return MatrixRow( &this->columnIndexes[ rowOffset ],
-                     &this->values[ rowOffset ],
-                     rowLength,
-                     1 );
+   return ConstMatrixRow( &this->columnIndexes[ rowOffset ],
+                          &this->values[ rowOffset ],
+                          rowLength,
+                          1 );
 }
 
 template< typename Real,
diff --git a/src/TNL/Matrices/Ellpack.h b/src/TNL/Matrices/Ellpack.h
index f7d0b9cca03beb9c60e5d5c860e76789d1f97c07..64bb18123d1011cb0678c62bcf57936d98617aed 100644
--- a/src/TNL/Matrices/Ellpack.h
+++ b/src/TNL/Matrices/Ellpack.h
@@ -35,6 +35,7 @@ class Ellpack : public Sparse< Real, Device, Index >
    typedef Ellpack< Real, Devices::Cuda, Index > CudaType;
    typedef Sparse< Real, Device, Index > BaseType;
    typedef typename BaseType::MatrixRow MatrixRow;
+   typedef SparseRow< const RealType, const IndexType > ConstMatrixRow;
 
    Ellpack();
 
@@ -128,7 +129,7 @@ class Ellpack : public Sparse< Real, Device, Index >
    MatrixRow getRow( const IndexType rowIndex );
 
    __cuda_callable__
-   const MatrixRow getRow( const IndexType rowIndex ) const;
+   ConstMatrixRow getRow( const IndexType rowIndex ) const;
 
    template< typename Vector >
    __cuda_callable__
diff --git a/src/TNL/Matrices/Ellpack_impl.h b/src/TNL/Matrices/Ellpack_impl.h
index 60c4c18a61121e858c250c336ff8723656202803..c183b4777368e49a30b9b8ba5e6efa97e52c2591 100644
--- a/src/TNL/Matrices/Ellpack_impl.h
+++ b/src/TNL/Matrices/Ellpack_impl.h
@@ -455,16 +455,16 @@ template< typename Real,
           typename Device,
           typename Index >
 __cuda_callable__
-const typename Ellpack< Real, Device, Index >::MatrixRow
+typename Ellpack< Real, Device, Index >::ConstMatrixRow
 Ellpack< Real, Device, Index >::
 getRow( const IndexType rowIndex ) const
 {
    //printf( "this->rowLengths = %d this = %p \n", this->rowLengths, this );
    IndexType rowBegin = DeviceDependentCode::getRowBegin( *this, rowIndex );
-   return MatrixRow( &this->columnIndexes[ rowBegin ],
-                     &this->values[ rowBegin ],
-                     this->rowLengths,
-                     DeviceDependentCode::getElementStep( *this ) );
+   return ConstMatrixRow( &this->columnIndexes[ rowBegin ],
+                          &this->values[ rowBegin ],
+                          this->rowLengths,
+                          DeviceDependentCode::getElementStep( *this ) );
 }
 
 template< typename Real,