diff --git a/examples/r_KLU_GLU.cpp b/examples/r_KLU_GLU.cpp
index a59441d8..9b1094e9 100644
--- a/examples/r_KLU_GLU.cpp
+++ b/examples/r_KLU_GLU.cpp
@@ -105,12 +105,12 @@ int main(int argc, char *argv[])
     rhs_file.close();
 
     //Now convert to CSR.
-    if (i < 1) { 
-      matrix_handler->coo2csr(A_coo, A,  "cpu");
+    if (i < 1) {
+      A->updateFromCoo(A_coo, ReSolve::memory::HOST);
       vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
       vec_rhs->setDataUpdated(ReSolve::memory::HOST);
     } else { 
-      matrix_handler->coo2csr(A_coo, A, "cuda");
+      A->updateFromCoo(A_coo, ReSolve::memory::DEVICE);
       vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
     }
     std::cout<<"COO to CSR completed. Expanded NNZ: "<< A->getNnzExpanded()<<std::endl;
diff --git a/examples/r_KLU_GLU_matrix_values_update.cpp b/examples/r_KLU_GLU_matrix_values_update.cpp
index 3dde51eb..b4f003c6 100644
--- a/examples/r_KLU_GLU_matrix_values_update.cpp
+++ b/examples/r_KLU_GLU_matrix_values_update.cpp
@@ -116,11 +116,11 @@ int main(int argc, char *argv[])
 
       //Now convert to CSR.
       if (i < 1) { 
-        matrix_handler->coo2csr(A_coo, A,  "cpu");
+        A->updateFromCoo(A_coo, ReSolve::memory::HOST);
         vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
         vec_rhs->setDataUpdated(ReSolve::memory::HOST);
       } else { 
-        matrix_handler->coo2csr(A_coo, A, "cuda");
+        A->updateFromCoo(A_coo, ReSolve::memory::DEVICE);
         vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
       }
       std::cout<<"COO to CSR completed. Expanded NNZ: "<< A->getNnzExpanded()<<std::endl;
diff --git a/examples/r_KLU_KLU.cpp b/examples/r_KLU_KLU.cpp
index bf4b9179..91450413 100644
--- a/examples/r_KLU_KLU.cpp
+++ b/examples/r_KLU_KLU.cpp
@@ -92,8 +92,7 @@ int main(int argc, char *argv[])
       vec_rhs = new vector_type(A->getNumRows());
       vec_x = new vector_type(A->getNumRows());
       vec_r = new vector_type(A->getNumRows());
-    }
-    else {
+    } else {
       ReSolve::io::readAndUpdateMatrix(mat_file, A_coo);
       ReSolve::io::readAndUpdateRhs(rhs_file, &rhs);
     }
@@ -107,11 +106,11 @@ int main(int argc, char *argv[])
 
     //Now convert to CSR.
     if (i < 2) { 
-      matrix_handler->coo2csr(A_coo, A, "cpu");
+      A->updateFromCoo(A_coo, ReSolve::memory::HOST);
       vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
       vec_rhs->setDataUpdated(ReSolve::memory::HOST);
     } else { 
-      matrix_handler->coo2csr(A_coo, A, "cpu");
+      A->updateFromCoo(A_coo, ReSolve::memory::HOST);
       vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
     }
     std::cout<<"COO to CSR completed. Expanded NNZ: "<< A->getNnzExpanded()<<std::endl;
diff --git a/examples/r_KLU_KLU_standalone.cpp b/examples/r_KLU_KLU_standalone.cpp
index 3821f36d..a2aa6d81 100644
--- a/examples/r_KLU_KLU_standalone.cpp
+++ b/examples/r_KLU_KLU_standalone.cpp
@@ -82,7 +82,7 @@ int main(int argc, char *argv[])
   rhs_file.close();
 
   //Now convert to CSR.
-  matrix_handler->coo2csr(A_coo, A, "cpu");
+  A->updateFromCoo(A_coo, ReSolve::memory::HOST);
   vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
   vec_rhs->setDataUpdated(ReSolve::memory::HOST);
   std::cout << "COO to CSR completed. Expanded NNZ: " << A->getNnzExpanded() << std::endl;
diff --git a/examples/r_KLU_rf.cpp b/examples/r_KLU_rf.cpp
index 2e99c1e7..1dbd4413 100644
--- a/examples/r_KLU_rf.cpp
+++ b/examples/r_KLU_rf.cpp
@@ -106,11 +106,11 @@ int main(int argc, char *argv[] )
 
     //Now convert to CSR.
     if (i < 2) { 
-      matrix_handler->coo2csr(A_coo, A, "cpu");
+      A->updateFromCoo(A_coo, ReSolve::memory::HOST);
       vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
       vec_rhs->setDataUpdated(ReSolve::memory::HOST);
     } else { 
-      matrix_handler->coo2csr(A_coo, A, "cuda");
+      A->updateFromCoo(A_coo, ReSolve::memory::DEVICE);
       vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
     }
     std::cout<<"COO to CSR completed. Expanded NNZ: "<< A->getNnzExpanded()<<std::endl;
diff --git a/examples/r_KLU_rf_FGMRES.cpp b/examples/r_KLU_rf_FGMRES.cpp
index 95857596..79ab7465 100644
--- a/examples/r_KLU_rf_FGMRES.cpp
+++ b/examples/r_KLU_rf_FGMRES.cpp
@@ -110,11 +110,11 @@ int main(int argc, char *argv[])
 
     //Now convert to CSR.
     if (i < 2) { 
-      matrix_handler->coo2csr(A_coo, A, "cpu");
+      A->updateFromCoo(A_coo, ReSolve::memory::HOST);
       vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
       vec_rhs->setDataUpdated(ReSolve::memory::HOST);
     } else { 
-      matrix_handler->coo2csr(A_coo,A, "cuda");
+      A->updateFromCoo(A_coo, ReSolve::memory::DEVICE);
       vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
     }
     std::cout<<"COO to CSR completed. Expanded NNZ: "<< A->getNnzExpanded()<<std::endl;
diff --git a/examples/r_KLU_rf_FGMRES_reuse_factorization.cpp b/examples/r_KLU_rf_FGMRES_reuse_factorization.cpp
index 04586e53..4610be4a 100644
--- a/examples/r_KLU_rf_FGMRES_reuse_factorization.cpp
+++ b/examples/r_KLU_rf_FGMRES_reuse_factorization.cpp
@@ -112,11 +112,11 @@ int main(int argc, char *argv[])
 
     //Now convert to CSR.
     if (i < 2) { 
-      matrix_handler->coo2csr(A_coo,A, "cpu");
+      A->updateFromCoo(A_coo, ReSolve::memory::HOST);
       vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
       vec_rhs->setDataUpdated(ReSolve::memory::HOST);
     } else { 
-      matrix_handler->coo2csr(A_coo, A, "cuda");
+      A->updateFromCoo(A_coo, ReSolve::memory::DEVICE);
       vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
     }
     std::cout<<"COO to CSR completed. Expanded NNZ: "<< A->getNnzExpanded()<<std::endl;
diff --git a/examples/r_KLU_rocSolverRf_FGMRES.cpp b/examples/r_KLU_rocSolverRf_FGMRES.cpp
index 4650b661..84a7bcf5 100644
--- a/examples/r_KLU_rocSolverRf_FGMRES.cpp
+++ b/examples/r_KLU_rocSolverRf_FGMRES.cpp
@@ -99,8 +99,7 @@ int main(int argc, char *argv[])
       vec_x->allocate(ReSolve::memory::HOST);//for KLU
       vec_x->allocate(ReSolve::memory::DEVICE);
       vec_r = new vector_type(A->getNumRows());
-    }
-    else {
+    } else {
       ReSolve::io::readAndUpdateMatrix(mat_file, A_coo);
       ReSolve::io::readAndUpdateRhs(rhs_file, &rhs);
     }
@@ -110,11 +109,11 @@ int main(int argc, char *argv[])
 
     //Now convert to CSR.
     if (i < 2) { 
-      matrix_handler->coo2csr(A_coo, A, "cpu");
+      A->updateFromCoo(A_coo, ReSolve::memory::HOST);
       vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
       vec_rhs->setDataUpdated(ReSolve::memory::HOST);
     } else { 
-      matrix_handler->coo2csr(A_coo,A, "hip");
+      A->updateFromCoo(A_coo, ReSolve::memory::DEVICE);
       vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
     }
     std::cout<<"COO to CSR completed. Expanded NNZ: "<< A->getNnzExpanded()<<std::endl;
diff --git a/examples/r_KLU_rocsolverrf.cpp b/examples/r_KLU_rocsolverrf.cpp
index 73e39131..89f74120 100644
--- a/examples/r_KLU_rocsolverrf.cpp
+++ b/examples/r_KLU_rocsolverrf.cpp
@@ -95,8 +95,7 @@ int main(int argc, char *argv[] )
       vec_rhs = new vector_type(A->getNumRows());
       vec_x = new vector_type(A->getNumRows());
       vec_r = new vector_type(A->getNumRows());
-    }
-    else {
+    } else {
       ReSolve::io::readAndUpdateMatrix(mat_file, A_coo);
       ReSolve::io::readAndUpdateRhs(rhs_file, &rhs);
     }
@@ -106,11 +105,11 @@ int main(int argc, char *argv[] )
 
     //Now convert to CSR.
     if (i < 2) { 
-      matrix_handler->coo2csr(A_coo, A, "cpu");
+      A->updateFromCoo(A_coo, ReSolve::memory::HOST);
       vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
       vec_rhs->setDataUpdated(ReSolve::memory::HOST);
     } else { 
-      matrix_handler->coo2csr(A_coo, A, "hip");
+      A->updateFromCoo(A_coo, ReSolve::memory::DEVICE);
       vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
     }
     std::cout<<"COO to CSR completed. Expanded NNZ: "<< A->getNnzExpanded()<<std::endl;
diff --git a/examples/r_KLU_rocsolverrf_redo_factorization.cpp b/examples/r_KLU_rocsolverrf_redo_factorization.cpp
index 232c8e69..54f28e8f 100644
--- a/examples/r_KLU_rocsolverrf_redo_factorization.cpp
+++ b/examples/r_KLU_rocsolverrf_redo_factorization.cpp
@@ -98,8 +98,7 @@ int main(int argc, char *argv[] )
       vec_rhs = new vector_type(A->getNumRows());
       vec_x = new vector_type(A->getNumRows());
       vec_r = new vector_type(A->getNumRows());
-    }
-    else {
+    } else {
       ReSolve::io::readAndUpdateMatrix(mat_file, A_coo);
       ReSolve::io::readAndUpdateRhs(rhs_file, &rhs);
     }
@@ -109,11 +108,11 @@ int main(int argc, char *argv[] )
 
     //Now convert to CSR.
     if (i < 2) { 
-      matrix_handler->coo2csr(A_coo, A, "cpu");
+      A->updateFromCoo(A_coo, ReSolve::memory::HOST);
       vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
       vec_rhs->setDataUpdated(ReSolve::memory::HOST);
     } else { 
-      matrix_handler->coo2csr(A_coo, A, "hip");
+      A->updateFromCoo(A_coo, ReSolve::memory::DEVICE);
       vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
     }
     std::cout<<"COO to CSR completed. Expanded NNZ: "<< A->getNnzExpanded()<<std::endl;
diff --git a/examples/r_SysSolver.cpp b/examples/r_SysSolver.cpp
index 88ca7324..bc170351 100644
--- a/examples/r_SysSolver.cpp
+++ b/examples/r_SysSolver.cpp
@@ -107,19 +107,16 @@ int main(int argc, char *argv[])
 
     //Now convert to CSR.
     if (i < 2) { 
-      matrix_handler->coo2csr(A_coo, A, "cpu");
+      A->updateFromCoo(A_coo, ReSolve::memory::HOST);
       vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
       vec_rhs->setDataUpdated(ReSolve::memory::HOST);
     } else { 
-      matrix_handler->coo2csr(A_coo, A, "cpu");
+      A->updateFromCoo(A_coo, ReSolve::memory::HOST);
       vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
     }
     std::cout<<"COO to CSR completed. Expanded NNZ: "<< A->getNnzExpanded()<<std::endl;
     //Now call direct solver
     solver->setMatrix(A);
-    // if (i == 0) {
-    //   solver->setupParameters(1, 0.1, false);
-    // }
     int status;
     if (i < 2){
       // solver->setup(A);
diff --git a/examples/r_SysSolverCuda.cpp b/examples/r_SysSolverCuda.cpp
index 43a2a31d..abffc586 100644
--- a/examples/r_SysSolverCuda.cpp
+++ b/examples/r_SysSolverCuda.cpp
@@ -106,11 +106,11 @@ int main(int argc, char *argv[])
 
     //Now convert to CSR.
     if (i < 1) { 
-      matrix_handler->coo2csr(A_coo, A,  "cpu");
+      A->updateFromCoo(A_coo, ReSolve::memory::HOST);
       vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
       vec_rhs->setDataUpdated(ReSolve::memory::HOST);
     } else { 
-      matrix_handler->coo2csr(A_coo, A, "cuda");
+      A->updateFromCoo(A_coo, ReSolve::memory::DEVICE);
       vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
     }
     std::cout<<"COO to CSR completed. Expanded NNZ: "<< A->getNnzExpanded()<<std::endl;
diff --git a/examples/r_SysSolverHip.cpp b/examples/r_SysSolverHip.cpp
index dc067a59..9a1be1fc 100644
--- a/examples/r_SysSolverHip.cpp
+++ b/examples/r_SysSolverHip.cpp
@@ -106,11 +106,11 @@ int main(int argc, char *argv[] )
 
     //Now convert to CSR.
     if (i < 2) { 
-      matrix_handler->coo2csr(A_coo, A, "cpu");
+      A->updateFromCoo(A_coo, ReSolve::memory::HOST);
       vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
       vec_rhs->setDataUpdated(ReSolve::memory::HOST);
     } else { 
-      matrix_handler->coo2csr(A_coo, A, "hip");
+      A->updateFromCoo(A_coo, ReSolve::memory::DEVICE);
       vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
     }
     std::cout<<"COO to CSR completed. Expanded NNZ: "<< A->getNnzExpanded()<<std::endl;
diff --git a/examples/r_SysSolverHipRefine.cpp b/examples/r_SysSolverHipRefine.cpp
index be46e182..bdc415fb 100644
--- a/examples/r_SysSolverHipRefine.cpp
+++ b/examples/r_SysSolverHipRefine.cpp
@@ -108,11 +108,11 @@ int main(int argc, char *argv[])
 
     //Now convert to CSR.
     if (i < 2) { 
-      matrix_handler->coo2csr(A_coo, A, "cpu");
+      A->updateFromCoo(A_coo, ReSolve::memory::HOST);
       vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
       vec_rhs->setDataUpdated(ReSolve::memory::HOST);
     } else { 
-      matrix_handler->coo2csr(A_coo, A, "hip");
+      A->updateFromCoo(A_coo, ReSolve::memory::DEVICE);
       vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
     }
     std::cout<<"COO to CSR completed. Expanded NNZ: "<< A->getNnzExpanded()<<std::endl;
diff --git a/examples/r_randGMRES.cpp b/examples/r_randGMRES.cpp
index b3f85cc9..a491ca2d 100644
--- a/examples/r_randGMRES.cpp
+++ b/examples/r_randGMRES.cpp
@@ -82,7 +82,7 @@ int main(int argc, char *argv[])
   mat_file.close();
   rhs_file.close();
 
-  matrix_handler->coo2csr(A_coo,A, "hip");
+  A->updateFromCoo(A_coo, ReSolve::memory::DEVICE);
   vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
   //Now call direct solver
   real_type norm_b;
diff --git a/examples/r_randGMRES_CUDA.cpp b/examples/r_randGMRES_CUDA.cpp
index 814ca2d4..9e451395 100644
--- a/examples/r_randGMRES_CUDA.cpp
+++ b/examples/r_randGMRES_CUDA.cpp
@@ -82,7 +82,7 @@ int main(int argc, char *argv[])
   mat_file.close();
   rhs_file.close();
 
-  matrix_handler->coo2csr(A_coo,A, "cuda");
+  A->updateFromCoo(A_coo, ReSolve::memory::DEVICE);
   vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
   //Now call direct solver
   real_type norm_b;
diff --git a/resolve/matrix/Csr.cpp b/resolve/matrix/Csr.cpp
index 0c08b641..914d853c 100644
--- a/resolve/matrix/Csr.cpp
+++ b/resolve/matrix/Csr.cpp
@@ -1,9 +1,16 @@
 #include <cstring>  // <-- includes memcpy
+#include <algorithm>
+#include <cassert>
 
 #include "Csr.hpp"
+#include "Coo.hpp"
+#include <resolve/utilities/misc/IndexValuePair.hpp>
+#include <resolve/utilities/logger/Logger.hpp>
 
 namespace ReSolve 
 {
+  using out = io::Logger;
+
   matrix::Csr::Csr()
   {
   }
@@ -20,6 +27,16 @@ namespace ReSolve
   {
   }
 
+  matrix::Csr::Csr(matrix::Coo* A_coo, memory::MemorySpace memspace)
+    : Sparse(A_coo->getNumRows(),
+             A_coo->getNumColumns(),
+             A_coo->getNnz(),
+             A_coo->symmetric(),
+             A_coo->expanded())
+  {
+    coo2csr(A_coo, memspace);
+  }
+
   matrix::Csr::~Csr()
   {
   }
@@ -234,7 +251,176 @@ namespace ReSolve
       default:
         return -1;
     } // switch
-}
+  }
+
+  int matrix::Csr::updateFromCoo(matrix::Coo* A_coo, memory::MemorySpace memspaceOut)
+  {
+    assert(n_            == A_coo->getNumRows());
+    assert(m_            == A_coo->getNumColumns());
+    assert(nnz_          == A_coo->getNnz());
+    assert(is_symmetric_ == A_coo->symmetric()); // <- Do we need to check for this?
+
+    return coo2csr(A_coo, memspaceOut);
+  }
+
+
+  int matrix::Csr::coo2csr(matrix::Coo* A_coo, memory::MemorySpace memspace)
+  {
+    //count nnzs first
+    index_type nnz_unpacked = 0;
+    index_type nnz = A_coo->getNnz();
+    index_type n = A_coo->getNumRows();
+    bool symmetric = A_coo->symmetric();
+    bool expanded = A_coo->expanded();
+
+    index_type* nnz_counts =  new index_type[n];
+    std::fill_n(nnz_counts, n, 0);
+    index_type* coo_rows = A_coo->getRowData(memory::HOST);
+    index_type* coo_cols = A_coo->getColData(memory::HOST);
+    real_type* coo_vals  = A_coo->getValues( memory::HOST);
+
+    index_type* diag_control = new index_type[n]; //for DEDUPLICATION of the diagonal
+    std::fill_n(diag_control, n, 0);
+    index_type nnz_unpacked_no_duplicates = 0;
+    index_type nnz_no_duplicates = nnz;
+
+
+    //maybe check if they exist?
+    for (index_type i = 0; i < nnz; ++i)
+    {
+      nnz_counts[coo_rows[i]]++;
+      nnz_unpacked++;
+      nnz_unpacked_no_duplicates++;
+      if ((coo_rows[i] != coo_cols[i])&& (symmetric) && (!expanded))
+      {
+        nnz_counts[coo_cols[i]]++;
+        nnz_unpacked++;
+        nnz_unpacked_no_duplicates++;
+      }
+      if (coo_rows[i] == coo_cols[i]){
+        if (diag_control[coo_rows[i]] > 0) {
+         //duplicate
+          nnz_unpacked_no_duplicates--;
+          nnz_no_duplicates--;
+        }
+        diag_control[coo_rows[i]]++;
+      }
+    }
+    this->setExpanded(true);
+    this->setNnzExpanded(nnz_unpacked_no_duplicates);
+    index_type* csr_ia = new index_type[n+1];
+    std::fill_n(csr_ia, n + 1, 0);
+    index_type* csr_ja = new index_type[nnz_unpacked];
+    real_type* csr_a = new real_type[nnz_unpacked];
+    index_type* nnz_shifts = new index_type[n];
+    std::fill_n(nnz_shifts, n , 0);
+
+    IndexValuePair* tmp = new IndexValuePair[nnz_unpacked]; 
+
+    csr_ia[0] = 0;
+
+    for (index_type i = 1; i < n + 1; ++i){
+      csr_ia[i] = csr_ia[i - 1] + nnz_counts[i - 1] - (diag_control[i-1] - 1);
+    }
+
+    int r, start;
+
+
+    for (index_type i = 0; i < nnz; ++i){
+      //which row
+      r = coo_rows[i];
+      start = csr_ia[r];
+
+      if ((start + nnz_shifts[r]) > nnz_unpacked) {
+        out::warning() << "index out of bounds (case 1) start: " << start << "nnz_shifts[" << r << "] = " << nnz_shifts[r] << std::endl;
+      }
+      if ((r == coo_cols[i]) && (diag_control[r] > 1)) {//diagonal, and there are duplicates
+        bool already_there = false;  
+        for (index_type j = start; j < start + nnz_shifts[r]; ++j)
+        {
+          index_type c = tmp[j].getIdx();
+          if (c == r) {
+            real_type val = tmp[j].getValue();
+            val += coo_vals[i];
+            tmp[j].setValue(val);
+            already_there = true;
+            out::warning() << " duplicate found, row " << c << " adding in place " << j << " current value: " << val << std::endl;
+          }  
+        }  
+        if (!already_there){ // first time this duplicates appears
+
+          tmp[start + nnz_shifts[r]].setIdx(coo_cols[i]);
+          tmp[start + nnz_shifts[r]].setValue(coo_vals[i]);
+
+          nnz_shifts[r]++;
+        }
+      } else {//not diagonal
+        tmp[start + nnz_shifts[r]].setIdx(coo_cols[i]);
+        tmp[start + nnz_shifts[r]].setValue(coo_vals[i]);
+        nnz_shifts[r]++;
+
+        if ((coo_rows[i] != coo_cols[i]) && (symmetric == 1))
+        {
+          r = coo_cols[i];
+          start = csr_ia[r];
+
+          if ((start + nnz_shifts[r]) > nnz_unpacked)
+            out::warning() << "index out of bounds (case 2) start: " << start << "nnz_shifts[" << r << "] = " << nnz_shifts[r] << std::endl;
+          tmp[start + nnz_shifts[r]].setIdx(coo_rows[i]);
+          tmp[start + nnz_shifts[r]].setValue(coo_vals[i]);
+          nnz_shifts[r]++;
+        }
+      }
+    }
+    //now sort whatever is inside rows
+
+    for (int i = 0; i < n; ++i)
+    {
+
+      //now sorting (and adding 1)
+      int colStart = csr_ia[i];
+      int colEnd = csr_ia[i + 1];
+      int length = colEnd - colStart;
+      std::sort(&tmp[colStart],&tmp[colStart] + length);
+    }
+
+    for (index_type i = 0; i < nnz_unpacked; ++i)
+    {
+      csr_ja[i] = tmp[i].getIdx();
+      csr_a[i] = tmp[i].getValue();
+    }
+#if 0
+    for (int i = 0; i<n; ++i){
+      printf("Row: %d \n", i);
+      for (int j = csr_ia[i]; j<csr_ia[i+1]; ++j){
+        printf("(%d %16.16f) ", csr_ja[j], csr_a[j]);
+      }
+      printf("\n");
+    }
+#endif
+    this->setNnz(nnz_no_duplicates);
+    this->updateData(csr_ia, csr_ja, csr_a, memory::HOST, memspace);
+    // if (memspace == "cpu"){
+    //   this->updateData(csr_ia, csr_ja, csr_a, memory::HOST, memory::HOST);
+    // } else {
+    //   if (memspace == "cuda"){      
+    //     this->updateData(csr_ia, csr_ja, csr_a, memory::HOST, memory::DEVICE);
+    //   } else if (memspace == "hip"){      
+    //     this->updateData(csr_ia, csr_ja, csr_a, memory::HOST, memory::DEVICE);
+    //   } else {
+    //     //display error
+    //   }
+    // }
+    delete [] nnz_counts;
+    delete [] tmp;
+    delete [] nnz_shifts;
+    delete [] csr_ia;
+    delete [] csr_ja;
+    delete [] csr_a;
+    delete [] diag_control; 
+
+    return 0;
+  }
 
 } // namespace ReSolve 
 
diff --git a/resolve/matrix/Csr.hpp b/resolve/matrix/Csr.hpp
index a5d8f682..7b3e3930 100644
--- a/resolve/matrix/Csr.hpp
+++ b/resolve/matrix/Csr.hpp
@@ -3,6 +3,9 @@
 
 namespace ReSolve { namespace matrix {
 
+  // Forward declaration of Coo
+  class Coo;
+
   class Csr : public Sparse
   {
     public:
@@ -15,6 +18,8 @@ namespace ReSolve { namespace matrix {
           index_type nnz,
           bool symmetric,
           bool expanded);
+      
+      Csr(matrix::Coo* mat, memory::MemorySpace memspace);
 
       ~Csr();
 
@@ -23,13 +28,18 @@ namespace ReSolve { namespace matrix {
       virtual real_type*  getValues( memory::MemorySpace memspace); 
 
       virtual int updateData(index_type* row_data, index_type* col_data, real_type* val_data, memory::MemorySpace memspaceIn, memory::MemorySpace memspaceOut); 
-      virtual int updateData(index_type* row_data, index_type* col_data, real_type* val_data, index_type new_nnz, memory::MemorySpace memspaceIn, memory::MemorySpace memspaceOut); 
+      virtual int updateData(index_type* row_data, index_type* col_data, real_type* val_data, index_type new_nnz, memory::MemorySpace memspaceIn, memory::MemorySpace memspaceOut);
 
       virtual int allocateMatrixData(memory::MemorySpace memspace); 
 
       virtual void print() {return;}
 
       virtual int copyData(memory::MemorySpace memspaceOut);
+
+      int updateFromCoo(matrix::Coo* mat, memory::MemorySpace memspaceOut);
+
+    private:
+      int coo2csr(matrix::Coo* mat, memory::MemorySpace memspace);
   };
 
 }} // namespace ReSolve::matrix
diff --git a/tests/functionality/testKLU.cpp b/tests/functionality/testKLU.cpp
index ad8410cd..40b7b140 100644
--- a/tests/functionality/testKLU.cpp
+++ b/tests/functionality/testKLU.cpp
@@ -50,11 +50,7 @@ int main(int argc, char *argv[])
     return -1;
   }
   ReSolve::matrix::Coo* A_coo = ReSolve::io::readMatrixFromFile(mat1);
-  ReSolve::matrix::Csr* A = new ReSolve::matrix::Csr(A_coo->getNumRows(),
-                                                     A_coo->getNumColumns(),
-                                                     A_coo->getNnz(),
-                                                     A_coo->symmetric(),
-                                                     A_coo->expanded());
+  ReSolve::matrix::Csr* A = new ReSolve::matrix::Csr(A_coo, ReSolve::memory::HOST);
   mat1.close();
 
   // Read first rhs vector
@@ -72,7 +68,6 @@ int main(int argc, char *argv[])
   rhs1_file.close();
 
   // Convert first matrix to CSR format
-  matrix_handler->coo2csr(A_coo, A, "cpu");
   vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
   vec_rhs->setDataUpdated(ReSolve::memory::HOST);
 
@@ -163,7 +158,7 @@ int main(int argc, char *argv[])
   ReSolve::io::readAndUpdateRhs(rhs2_file, &rhs);
   rhs2_file.close();
 
-  matrix_handler->coo2csr(A_coo, A, "cpu");
+  A->updateFromCoo(A_coo, ReSolve::memory::HOST);
   vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
 
   // and solve it too
diff --git a/tests/functionality/testKLU_GLU.cpp b/tests/functionality/testKLU_GLU.cpp
index cb14c2dc..275737f0 100644
--- a/tests/functionality/testKLU_GLU.cpp
+++ b/tests/functionality/testKLU_GLU.cpp
@@ -56,11 +56,7 @@ int main(int argc, char *argv[])
     return -1;
   }
   ReSolve::matrix::Coo* A_coo = ReSolve::io::readMatrixFromFile(mat1);
-  ReSolve::matrix::Csr* A = new ReSolve::matrix::Csr(A_coo->getNumRows(),
-                                                     A_coo->getNumColumns(),
-                                                     A_coo->getNnz(),
-                                                     A_coo->symmetric(),
-                                                     A_coo->expanded());
+  ReSolve::matrix::Csr* A = new ReSolve::matrix::Csr(A_coo, ReSolve::memory::HOST);
   mat1.close();
 
   // Read first rhs vector
@@ -80,7 +76,6 @@ int main(int argc, char *argv[])
   rhs1_file.close();
 
   // Convert first matrix to CSR format
-  matrix_handler->coo2csr(A_coo, A, "cpu");
   vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
   vec_rhs->setDataUpdated(ReSolve::memory::HOST);
 
@@ -94,8 +89,7 @@ int main(int argc, char *argv[])
   status = KLU->factorize();
   error_sum += status;
 
-// but DO NOT SOLVE with KLU!
-
+  // but DO NOT SOLVE the rest with KLU!
 
   matrix_type* L = KLU->getLFactor();
   matrix_type* U = KLU->getUFactor();
@@ -110,8 +104,6 @@ int main(int argc, char *argv[])
   error_sum += status;
   std::cout<<"GLU solve status: "<<status<<std::endl;      
 
-
-
   vector_type* vec_test;
   vector_type* vec_diff;
   vec_test  = new vector_type(A->getNumRows());
@@ -125,16 +117,13 @@ int main(int argc, char *argv[])
   vec_r->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
   vec_diff->update(x_data, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
 
-  // real_type normXmatrix1 = sqrt(vector_handler->dot(vec_test, vec_test, "cuda"));
   matrix_handler->setValuesChanged(true, "cuda");
   status = matrix_handler->matvec(A, vec_x, vec_r, &ONE, &MINUSONE,"csr","cuda"); 
   error_sum += status;
   
   real_type normRmatrix1 = sqrt(vector_handler->dot(vec_r, vec_r, "cuda"));
 
-
   //for testing only - control
-  
   real_type normXtrue = sqrt(vector_handler->dot(vec_x, vec_x, "cuda"));
   real_type normB1 = sqrt(vector_handler->dot(vec_rhs, vec_rhs, "cuda"));
   
@@ -148,11 +137,10 @@ int main(int argc, char *argv[])
   status = matrix_handler->matvec(A, vec_test, vec_r, &ONE, &MINUSONE, "csr", "cuda"); 
 
   vec_x->update(vec_x->getData(ReSolve::memory::DEVICE), ReSolve::memory::DEVICE, ReSolve::memory::HOST);
-  //status = matrix_handler->matvec(A, vec_test, vec_r, &ONE, &MINUSONE,"csr", "cuda"); 
   error_sum += status;
   real_type exactSol_normRmatrix1 = sqrt(vector_handler->dot(vec_r, vec_r, "cuda"));
+
   //evaluate the residual ON THE CPU using COMPUTED solution
- 
   vec_r->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
 
   status = matrix_handler->matvec(A, vec_x, vec_r, &ONE, &MINUSONE,"csr", "cpu");
@@ -189,7 +177,7 @@ int main(int argc, char *argv[])
   ReSolve::io::readAndUpdateRhs(rhs2_file, &rhs);
   rhs2_file.close();
 
-  matrix_handler->coo2csr(A_coo, A, "cuda");
+  A->updateFromCoo(A_coo, ReSolve::memory::DEVICE);
   vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
 
   status = GLU->refactorize();
diff --git a/tests/functionality/testKLU_Rf.cpp b/tests/functionality/testKLU_Rf.cpp
index 12f1927e..5aef80cf 100644
--- a/tests/functionality/testKLU_Rf.cpp
+++ b/tests/functionality/testKLU_Rf.cpp
@@ -56,11 +56,7 @@ int main(int argc, char *argv[])
     return -1;
   }
   ReSolve::matrix::Coo* A_coo = ReSolve::io::readMatrixFromFile(mat1);
-  ReSolve::matrix::Csr* A = new ReSolve::matrix::Csr(A_coo->getNumRows(),
-                                                     A_coo->getNumColumns(),
-                                                     A_coo->getNnz(),
-                                                     A_coo->symmetric(),
-                                                     A_coo->expanded());
+  ReSolve::matrix::Csr* A = new ReSolve::matrix::Csr(A_coo, ReSolve::memory::HOST);
   mat1.close();
 
   // Read first rhs vector
@@ -78,7 +74,6 @@ int main(int argc, char *argv[])
   rhs1_file.close();
 
   // Convert first matrix to CSR format
-  matrix_handler->coo2csr(A_coo, A, "cpu");
   vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
   vec_rhs->setDataUpdated(ReSolve::memory::HOST);
 
@@ -109,16 +104,13 @@ int main(int argc, char *argv[])
   vec_r->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
   vec_diff->update(x_data, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
 
-  // real_type normXmatrix1 = sqrt(vector_handler->dot(vec_test, vec_test, "cuda"));
   matrix_handler->setValuesChanged(true, "cuda");
   status = matrix_handler->matvec(A, vec_x, vec_r, &ONE, &MINUSONE,"csr","cuda"); 
   error_sum += status;
   
   real_type normRmatrix1 = sqrt(vector_handler->dot(vec_r, vec_r, "cuda"));
 
-
   //for testing only - control
-  
   real_type normXtrue = sqrt(vector_handler->dot(vec_x, vec_x, "cuda"));
   real_type normB1 = sqrt(vector_handler->dot(vec_rhs, vec_rhs, "cuda"));
   
@@ -132,8 +124,8 @@ int main(int argc, char *argv[])
   status = matrix_handler->matvec(A, vec_test, vec_r, &ONE, &MINUSONE,"csr", "cuda"); 
   error_sum += status;
   real_type exactSol_normRmatrix1 = sqrt(vector_handler->dot(vec_r, vec_r, "cuda"));
+
   //evaluate the residual ON THE CPU using COMPUTED solution
- 
   vec_r->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
 
   status = matrix_handler->matvec(A, vec_x, vec_r, &ONE, &MINUSONE,"csr", "cpu");
@@ -184,7 +176,7 @@ int main(int argc, char *argv[])
   ReSolve::io::readAndUpdateRhs(rhs2_file, &rhs);
   rhs2_file.close();
 
-  matrix_handler->coo2csr(A_coo, A, "cuda");
+  A->updateFromCoo(A_coo, ReSolve::memory::DEVICE);
   vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
 
   status = Rf->refactorize();
diff --git a/tests/functionality/testKLU_Rf_FGMRES.cpp b/tests/functionality/testKLU_Rf_FGMRES.cpp
index 6fa59a17..664ea446 100644
--- a/tests/functionality/testKLU_Rf_FGMRES.cpp
+++ b/tests/functionality/testKLU_Rf_FGMRES.cpp
@@ -61,11 +61,7 @@ int main(int argc, char *argv[])
     return -1;
   }
   ReSolve::matrix::Coo* A_coo = ReSolve::io::readMatrixFromFile(mat1);
-  ReSolve::matrix::Csr* A = new ReSolve::matrix::Csr(A_coo->getNumRows(),
-                                                     A_coo->getNumColumns(),
-                                                     A_coo->getNnz(),
-                                                     A_coo->symmetric(),
-                                                     A_coo->expanded());
+  ReSolve::matrix::Csr* A = new ReSolve::matrix::Csr(A_coo, ReSolve::memory::HOST);
   mat1.close();
 
   // Read first rhs vector
@@ -83,7 +79,6 @@ int main(int argc, char *argv[])
   rhs1_file.close();
 
   // Convert first matrix to CSR format
-  matrix_handler->coo2csr(A_coo, A, "cpu");
   vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
   vec_rhs->setDataUpdated(ReSolve::memory::HOST);
 
@@ -115,7 +110,6 @@ int main(int argc, char *argv[])
   vec_r->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
   vec_diff->update(x_data, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
 
-  // real_type normXmatrix1 = sqrt(vector_handler->dot(vec_test, vec_test, ReSolve::memory::DEVICE));
   matrix_handler->setValuesChanged(true, "cuda");
   //evaluate the residual ||b-Ax||
   status = matrix_handler->matvec(A, vec_x, vec_r, &ONE, &MINUSONE,"csr","cuda"); 
@@ -200,7 +194,7 @@ int main(int argc, char *argv[])
   ReSolve::io::readAndUpdateRhs(rhs2_file, &rhs);
   rhs2_file.close();
 
-  matrix_handler->coo2csr(A_coo, A, "cuda");
+  A->updateFromCoo(A_coo, ReSolve::memory::DEVICE);
   vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
   Rf->setNumericalProperties(1e-12, 1e-1);
 
diff --git a/tests/functionality/testKLU_RocSolver.cpp b/tests/functionality/testKLU_RocSolver.cpp
index 3e3ef4b0..9ee8ad77 100644
--- a/tests/functionality/testKLU_RocSolver.cpp
+++ b/tests/functionality/testKLU_RocSolver.cpp
@@ -56,11 +56,7 @@ int main(int argc, char *argv[])
     return -1;
   }
   ReSolve::matrix::Coo* A_coo = ReSolve::io::readMatrixFromFile(mat1);
-  ReSolve::matrix::Csr* A = new ReSolve::matrix::Csr(A_coo->getNumRows(),
-                                                     A_coo->getNumColumns(),
-                                                     A_coo->getNnz(),
-                                                     A_coo->symmetric(),
-                                                     A_coo->expanded());
+  ReSolve::matrix::Csr* A = new ReSolve::matrix::Csr(A_coo, ReSolve::memory::HOST);
   mat1.close();
 
   // Read first rhs vector
@@ -80,7 +76,6 @@ int main(int argc, char *argv[])
   rhs1_file.close();
 
   // Convert first matrix to CSR format
-  matrix_handler->coo2csr(A_coo, A, "cpu");
   vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
   vec_rhs->setDataUpdated(ReSolve::memory::HOST);
 
@@ -186,7 +181,7 @@ int main(int argc, char *argv[])
   ReSolve::io::readAndUpdateRhs(rhs2_file, &rhs);
   rhs2_file.close();
 
-  matrix_handler->coo2csr(A_coo, A, "hip");
+  A->updateFromCoo(A_coo, ReSolve::memory::DEVICE);
   vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
 
   // this hangs up
diff --git a/tests/functionality/testKLU_RocSolver_FGMRES.cpp b/tests/functionality/testKLU_RocSolver_FGMRES.cpp
index a1ecfce8..0781a83a 100644
--- a/tests/functionality/testKLU_RocSolver_FGMRES.cpp
+++ b/tests/functionality/testKLU_RocSolver_FGMRES.cpp
@@ -61,11 +61,7 @@ int main(int argc, char *argv[])
     return -1;
   }
   ReSolve::matrix::Coo* A_coo = ReSolve::io::readMatrixFromFile(mat1);
-  ReSolve::matrix::Csr* A = new ReSolve::matrix::Csr(A_coo->getNumRows(),
-                                                     A_coo->getNumColumns(),
-                                                     A_coo->getNnz(),
-                                                     A_coo->symmetric(),
-                                                     A_coo->expanded());
+  ReSolve::matrix::Csr* A = new ReSolve::matrix::Csr(A_coo, ReSolve::memory::HOST);
   mat1.close();
 
   // Read first rhs vector
@@ -83,7 +79,6 @@ int main(int argc, char *argv[])
   rhs1_file.close();
 
   // Convert first matrix to CSR format
-  matrix_handler->coo2csr(A_coo, A, "cpu");
   vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
   vec_rhs->setDataUpdated(ReSolve::memory::HOST);
 
@@ -115,7 +110,6 @@ int main(int argc, char *argv[])
   vec_r->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
   vec_diff->update(x_data, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
 
-  // real_type normXmatrix1 = sqrt(vector_handler->dot(vec_test, vec_test, ReSolve::memory::DEVICE));
   matrix_handler->setValuesChanged(true, "hip");
   //evaluate the residual ||b-Ax||
   status = matrix_handler->matvec(A, vec_x, vec_r, &ONE, &MINUSONE,"csr","hip"); 
@@ -125,7 +119,6 @@ int main(int argc, char *argv[])
 
 
   //for testing only - control
-
   real_type normXtrue = sqrt(vector_handler->dot(vec_x, vec_x, "hip"));
   real_type normB1 = sqrt(vector_handler->dot(vec_rhs, vec_rhs, "hip"));
 
@@ -139,8 +132,8 @@ int main(int argc, char *argv[])
   status = matrix_handler->matvec(A, vec_test, vec_r, &ONE, &MINUSONE,"csr", "hip"); 
   error_sum += status;
   real_type exactSol_normRmatrix1 = sqrt(vector_handler->dot(vec_r, vec_r, "hip"));
-  //evaluate the residual ON THE CPU using COMPUTED solution
 
+  //evaluate the residual ON THE CPU using COMPUTED solution
   vec_r->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
 
   status = matrix_handler->matvec(A, vec_x, vec_r, &ONE, &MINUSONE,"csr", "cpu");
@@ -197,7 +190,7 @@ int main(int argc, char *argv[])
   ReSolve::io::readAndUpdateRhs(rhs2_file, &rhs);
   rhs2_file.close();
 
-  matrix_handler->coo2csr(A_coo, A, "hip");
+  A->updateFromCoo(A_coo, ReSolve::memory::DEVICE);
   vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
 
   status = Rf->refactorize();
diff --git a/tests/functionality/testSysCuda.cpp b/tests/functionality/testSysCuda.cpp
index 76655af6..67068ac8 100644
--- a/tests/functionality/testSysCuda.cpp
+++ b/tests/functionality/testSysCuda.cpp
@@ -54,11 +54,7 @@ int main(int argc, char *argv[])
     return -1;
   }
   ReSolve::matrix::Coo* A_coo = ReSolve::io::readMatrixFromFile(mat1);
-  ReSolve::matrix::Csr* A = new ReSolve::matrix::Csr(A_coo->getNumRows(),
-                                                     A_coo->getNumColumns(),
-                                                     A_coo->getNnz(),
-                                                     A_coo->symmetric(),
-                                                     A_coo->expanded());
+  ReSolve::matrix::Csr* A = new ReSolve::matrix::Csr(A_coo, ReSolve::memory::HOST);
   mat1.close();
 
   // Read first rhs vector
@@ -78,7 +74,6 @@ int main(int argc, char *argv[])
   rhs1_file.close();
 
   // Convert first matrix to CSR format
-  matrix_handler->coo2csr(A_coo, A, "cpu");
   vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
   vec_rhs->setDataUpdated(ReSolve::memory::HOST);
 
@@ -116,7 +111,6 @@ int main(int argc, char *argv[])
   vec_r->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
   vec_diff->update(x_data, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
 
-  // real_type normXmatrix1 = sqrt(vector_handler->dot(vec_test, vec_test, "cuda"));
   matrix_handler->setValuesChanged(true, "cuda");
   status = matrix_handler->matvec(A, vec_x, vec_r, &ONE, &MINUSONE,"csr","cuda"); 
   error_sum += status;
@@ -179,7 +173,7 @@ int main(int argc, char *argv[])
   ReSolve::io::readAndUpdateRhs(rhs2_file, &rhs);
   rhs2_file.close();
 
-  matrix_handler->coo2csr(A_coo, A, "cuda");
+  A->updateFromCoo(A_coo, ReSolve::memory::DEVICE);
   vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
 
   status = solver->refactorize();
diff --git a/tests/functionality/testSysHipRefine.cpp b/tests/functionality/testSysHipRefine.cpp
index 2c92251d..91d7b918 100644
--- a/tests/functionality/testSysHipRefine.cpp
+++ b/tests/functionality/testSysHipRefine.cpp
@@ -60,11 +60,7 @@ int main(int argc, char *argv[])
     return -1;
   }
   ReSolve::matrix::Coo* A_coo = ReSolve::io::readMatrixFromFile(mat1);
-  ReSolve::matrix::Csr* A = new ReSolve::matrix::Csr(A_coo->getNumRows(),
-                                                     A_coo->getNumColumns(),
-                                                     A_coo->getNnz(),
-                                                     A_coo->symmetric(),
-                                                     A_coo->expanded());
+  ReSolve::matrix::Csr* A = new ReSolve::matrix::Csr(A_coo, ReSolve::memory::HOST);
   mat1.close();
 
   // Read first rhs vector
@@ -82,7 +78,6 @@ int main(int argc, char *argv[])
   rhs1_file.close();
 
   // Convert first matrix to CSR format
-  matrix_handler.coo2csr(A_coo, A, "cpu");
   vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::HOST);
   vec_rhs->setDataUpdated(ReSolve::memory::HOST);
 
@@ -114,7 +109,6 @@ int main(int argc, char *argv[])
   vec_r->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
   vec_diff->update(x_data, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
 
-  // real_type normXmatrix1 = sqrt(vector_handler.dot(vec_test, vec_test, ReSolve::memory::DEVICE));
   matrix_handler.setValuesChanged(true, "hip");
   //evaluate the residual ||b-Ax||
   status = matrix_handler.matvec(A, vec_x, vec_r, &ONE, &MINUSONE,"csr","hip"); 
@@ -122,9 +116,7 @@ int main(int argc, char *argv[])
 
   real_type normRmatrix1 = sqrt(vector_handler.dot(vec_r, vec_r, "hip"));
 
-
   //for testing only - control
-
   real_type normXtrue = sqrt(vector_handler.dot(vec_x, vec_x, "hip"));
   real_type normB1 = sqrt(vector_handler.dot(vec_rhs, vec_rhs, "hip"));
 
@@ -180,7 +172,7 @@ int main(int argc, char *argv[])
   ReSolve::io::readAndUpdateRhs(rhs2_file, &rhs);
   rhs2_file.close();
 
-  matrix_handler.coo2csr(A_coo, A, "hip");
+  A->updateFromCoo(A_coo, ReSolve::memory::DEVICE);
   vec_rhs->update(rhs, ReSolve::memory::HOST, ReSolve::memory::DEVICE);
 
   status = solver.refactorize();