From 91cf1e3bac0b2517419d2d1a66be4bdbfb38dba3 Mon Sep 17 00:00:00 2001 From: Pratik Nayak Date: Mon, 7 Oct 2019 22:28:03 +0200 Subject: [PATCH] Review update. + Update some forgotten names. + Generic improvements. --- core/matrix/dense.cpp | 16 +++++++--------- core/matrix/dense_kernels.hpp | 8 ++++---- core/matrix/sparsity_csr.cpp | 2 +- core/matrix/sparsity_csr_kernels.hpp | 2 +- cuda/matrix/dense_kernels.cu | 6 +++--- cuda/matrix/sparsity_csr_kernels.cu | 2 +- omp/matrix/dense_kernels.cpp | 6 +++--- omp/matrix/sparsity_csr_kernels.cpp | 8 ++++---- omp/test/matrix/sparsity_csr_kernels.cpp | 7 ++++--- reference/matrix/dense_kernels.cpp | 6 +++--- reference/matrix/sparsity_csr_kernels.cpp | 8 ++++---- reference/test/matrix/sparsity_csr.cpp | 4 ++-- reference/test/matrix/sparsity_csr_kernels.cpp | 4 ++-- 13 files changed, 39 insertions(+), 40 deletions(-) diff --git a/core/matrix/dense.cpp b/core/matrix/dense.cpp index 32e1c35fa02..72c23343857 100644 --- a/core/matrix/dense.cpp +++ b/core/matrix/dense.cpp @@ -76,7 +76,7 @@ GKO_REGISTER_OPERATION(convert_to_csr, dense::convert_to_csr); GKO_REGISTER_OPERATION(convert_to_ell, dense::convert_to_ell); GKO_REGISTER_OPERATION(convert_to_hybrid, dense::convert_to_hybrid); GKO_REGISTER_OPERATION(convert_to_sellp, dense::convert_to_sellp); -GKO_REGISTER_OPERATION(convert_to_sparsity, dense::convert_to_sparsity); +GKO_REGISTER_OPERATION(convert_to_sparsity_csr, dense::convert_to_sparsity_csr); } // namespace dense @@ -444,10 +444,9 @@ void Dense::move_to(Sellp *result) template void Dense::convert_to(SparsityCsr *result) const { - conversion_helper( - result, this, - dense::template make_convert_to_sparsity *&>); + conversion_helper(result, this, + dense::template make_convert_to_sparsity_csr< + decltype(result), const Dense *&>); } @@ -461,10 +460,9 @@ void Dense::move_to(SparsityCsr *result) template void Dense::convert_to(SparsityCsr *result) const { - conversion_helper( - result, this, - dense::template make_convert_to_sparsity *&>); + conversion_helper(result, this, + dense::template make_convert_to_sparsity_csr< + decltype(result), const Dense *&>); } diff --git a/core/matrix/dense_kernels.hpp b/core/matrix/dense_kernels.hpp index 5e6169d4732..4857fb81db9 100644 --- a/core/matrix/dense_kernels.hpp +++ b/core/matrix/dense_kernels.hpp @@ -98,10 +98,10 @@ namespace kernels { matrix::Sellp<_type, _prec> *other, \ const matrix::Dense<_type> *source) -#define GKO_DECLARE_DENSE_CONVERT_TO_SPARSITY_CSR_KERNEL(_type, _prec) \ - void convert_to_sparsity(std::shared_ptr exec, \ - matrix::SparsityCsr<_type, _prec> *other, \ - const matrix::Dense<_type> *source) +#define GKO_DECLARE_DENSE_CONVERT_TO_SPARSITY_CSR_KERNEL(_type, _prec) \ + void convert_to_sparsity_csr(std::shared_ptr exec, \ + matrix::SparsityCsr<_type, _prec> *other, \ + const matrix::Dense<_type> *source) #define GKO_DECLARE_DENSE_COUNT_NONZEROS_KERNEL(_type) \ void count_nonzeros(std::shared_ptr exec, \ diff --git a/core/matrix/sparsity_csr.cpp b/core/matrix/sparsity_csr.cpp index 985fd03bc1e..3c7cfa1363f 100644 --- a/core/matrix/sparsity_csr.cpp +++ b/core/matrix/sparsity_csr.cpp @@ -171,7 +171,7 @@ SparsityCsr::to_adjacency_matrix() const GKO_ASSERT_IS_SQUARE_MATRIX(this); size_type num_diagonal_elements = 0; exec->run(sparsity_csr::make_count_num_diagonal_elements( - this, num_diagonal_elements)); + this, &num_diagonal_elements)); auto adj_mat = SparsityCsr::create(exec, this->get_size(), this->get_num_nonzeros() - num_diagonal_elements); diff --git a/core/matrix/sparsity_csr_kernels.hpp b/core/matrix/sparsity_csr_kernels.hpp index 1fd40ca06a3..f9af3dcdffa 100644 --- a/core/matrix/sparsity_csr_kernels.hpp +++ b/core/matrix/sparsity_csr_kernels.hpp @@ -68,7 +68,7 @@ namespace kernels { void count_num_diagonal_elements( \ std::shared_ptr exec, \ const matrix::SparsityCsr *matrix, \ - size_type &num_diagonal_elements) + size_type *num_diagonal_elements) #define GKO_DECLARE_SPARSITY_CSR_TRANSPOSE_KERNEL(ValueType, IndexType) \ void transpose(std::shared_ptr exec, \ diff --git a/cuda/matrix/dense_kernels.cu b/cuda/matrix/dense_kernels.cu index 0b0ef1e80b3..a74c431599a 100644 --- a/cuda/matrix/dense_kernels.cu +++ b/cuda/matrix/dense_kernels.cu @@ -732,9 +732,9 @@ GKO_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE( template -void convert_to_sparsity(std::shared_ptr exec, - matrix::SparsityCsr *result, - const matrix::Dense *source) +void convert_to_sparsity_csr(std::shared_ptr exec, + matrix::SparsityCsr *result, + const matrix::Dense *source) GKO_NOT_IMPLEMENTED; GKO_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE( diff --git a/cuda/matrix/sparsity_csr_kernels.cu b/cuda/matrix/sparsity_csr_kernels.cu index 4c5ded1af43..2bdb8372030 100644 --- a/cuda/matrix/sparsity_csr_kernels.cu +++ b/cuda/matrix/sparsity_csr_kernels.cu @@ -92,7 +92,7 @@ template void count_num_diagonal_elements( std::shared_ptr exec, const matrix::SparsityCsr *matrix, - size_type &num_diagonal_elements) GKO_NOT_IMPLEMENTED; + size_type *num_diagonal_elements) GKO_NOT_IMPLEMENTED; GKO_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE( GKO_DECLARE_SPARSITY_CSR_COUNT_NUM_DIAGONAL_ELEMENTS_KERNEL); diff --git a/omp/matrix/dense_kernels.cpp b/omp/matrix/dense_kernels.cpp index 2e2047604bf..5654b40d753 100644 --- a/omp/matrix/dense_kernels.cpp +++ b/omp/matrix/dense_kernels.cpp @@ -457,9 +457,9 @@ GKO_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE( template -void convert_to_sparsity(std::shared_ptr exec, - matrix::SparsityCsr *result, - const matrix::Dense *source) +void convert_to_sparsity_csr(std::shared_ptr exec, + matrix::SparsityCsr *result, + const matrix::Dense *source) { auto num_rows = result->get_size()[0]; auto num_cols = result->get_size()[1]; diff --git a/omp/matrix/sparsity_csr_kernels.cpp b/omp/matrix/sparsity_csr_kernels.cpp index 3800faf5167..ed538793c3c 100644 --- a/omp/matrix/sparsity_csr_kernels.cpp +++ b/omp/matrix/sparsity_csr_kernels.cpp @@ -68,6 +68,7 @@ void spmv(std::shared_ptr exec, { auto row_ptrs = a->get_const_row_ptrs(); auto col_idxs = a->get_const_col_idxs(); + auto val = a->get_const_value()[0]; #pragma omp parallel for for (size_type row = 0; row < a->get_size()[0]; ++row) { @@ -76,7 +77,6 @@ void spmv(std::shared_ptr exec, } for (size_type k = row_ptrs[row]; k < static_cast(row_ptrs[row + 1]); ++k) { - auto val = one(); auto col = col_idxs[k]; for (size_type j = 0; j < c->get_size()[1]; ++j) { c->at(row, j) += val * b->at(col, j); @@ -101,6 +101,7 @@ void advanced_spmv(std::shared_ptr exec, auto col_idxs = a->get_const_col_idxs(); auto valpha = alpha->at(0, 0); auto vbeta = beta->at(0, 0); + auto val = a->get_const_value()[0]; #pragma omp parallel for for (size_type row = 0; row < a->get_size()[0]; ++row) { @@ -109,7 +110,6 @@ void advanced_spmv(std::shared_ptr exec, } for (size_type k = row_ptrs[row]; k < static_cast(row_ptrs[row + 1]); ++k) { - auto val = one(); auto col = col_idxs[k]; for (size_type j = 0; j < c->get_size()[1]; ++j) { c->at(row, j) += valpha * val * b->at(col, j); @@ -126,7 +126,7 @@ template void count_num_diagonal_elements( std::shared_ptr exec, const matrix::SparsityCsr *matrix, - size_type &num_diagonal_elements) + size_type *num_diagonal_elements) { auto num_rows = matrix->get_size()[0]; auto row_ptrs = matrix->get_const_row_ptrs(); @@ -139,7 +139,7 @@ void count_num_diagonal_elements( } } } - num_diagonal_elements = num_diag; + *num_diagonal_elements = num_diag; } GKO_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE( diff --git a/omp/test/matrix/sparsity_csr_kernels.cpp b/omp/test/matrix/sparsity_csr_kernels.cpp index 4e2a189dbf7..91852dbcb53 100644 --- a/omp/test/matrix/sparsity_csr_kernels.cpp +++ b/omp/test/matrix/sparsity_csr_kernels.cpp @@ -225,9 +225,9 @@ TEST_F(SparsityCsr, CountsNumberOfDiagElementsIsEqualToRef) gko::size_type d_num_diags = 0; gko::kernels::reference::sparsity_csr::count_num_diagonal_elements( - ref, mtx.get(), num_diags); + ref, mtx.get(), &num_diags); gko::kernels::omp::sparsity_csr::count_num_diagonal_elements( - omp, dmtx.get(), d_num_diags); + omp, dmtx.get(), &d_num_diags); ASSERT_EQ(d_num_diags, num_diags); } @@ -238,11 +238,12 @@ TEST_F(SparsityCsr, RemovesDiagElementsKernelIsEquivalentToRef) set_up_apply_data(); gko::size_type num_diags = 0; gko::kernels::reference::sparsity_csr::count_num_diagonal_elements( - ref, mtx.get(), num_diags); + ref, mtx.get(), &num_diags); auto tmp = Mtx::create(ref, mtx->get_size(), mtx->get_num_nonzeros() - num_diags); auto d_tmp = Mtx::create(omp, dmtx->get_size(), dmtx->get_num_nonzeros() - num_diags); + gko::kernels::reference::sparsity_csr::remove_diagonal_elements( ref, tmp.get(), mtx->get_const_row_ptrs(), mtx->get_const_col_idxs()); gko::kernels::omp::sparsity_csr::remove_diagonal_elements( diff --git a/reference/matrix/dense_kernels.cpp b/reference/matrix/dense_kernels.cpp index ca1c0921ad1..82a20a8b1a4 100644 --- a/reference/matrix/dense_kernels.cpp +++ b/reference/matrix/dense_kernels.cpp @@ -416,9 +416,9 @@ GKO_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE( template -void convert_to_sparsity(std::shared_ptr exec, - matrix::SparsityCsr *result, - const matrix::Dense *source) +void convert_to_sparsity_csr(std::shared_ptr exec, + matrix::SparsityCsr *result, + const matrix::Dense *source) { auto num_rows = result->get_size()[0]; auto num_cols = result->get_size()[1]; diff --git a/reference/matrix/sparsity_csr_kernels.cpp b/reference/matrix/sparsity_csr_kernels.cpp index 99698d488ff..f3c1cc46df2 100644 --- a/reference/matrix/sparsity_csr_kernels.cpp +++ b/reference/matrix/sparsity_csr_kernels.cpp @@ -65,6 +65,7 @@ void spmv(std::shared_ptr exec, { auto row_ptrs = a->get_const_row_ptrs(); auto col_idxs = a->get_const_col_idxs(); + auto val = a->get_const_value()[0]; for (size_type row = 0; row < a->get_size()[0]; ++row) { for (size_type j = 0; j < c->get_size()[1]; ++j) { @@ -72,7 +73,6 @@ void spmv(std::shared_ptr exec, } for (size_type k = row_ptrs[row]; k < static_cast(row_ptrs[row + 1]); ++k) { - auto val = one(); auto col = col_idxs[k]; for (size_type j = 0; j < c->get_size()[1]; ++j) { c->at(row, j) += val * b->at(col, j); @@ -97,6 +97,7 @@ void advanced_spmv(std::shared_ptr exec, auto col_idxs = a->get_const_col_idxs(); auto valpha = alpha->at(0, 0); auto vbeta = beta->at(0, 0); + auto val = a->get_const_value()[0]; for (size_type row = 0; row < a->get_size()[0]; ++row) { for (size_type j = 0; j < c->get_size()[1]; ++j) { @@ -104,7 +105,6 @@ void advanced_spmv(std::shared_ptr exec, } for (size_type k = row_ptrs[row]; k < static_cast(row_ptrs[row + 1]); ++k) { - auto val = one(); auto col = col_idxs[k]; for (size_type j = 0; j < c->get_size()[1]; ++j) { c->at(row, j) += valpha * val * b->at(col, j); @@ -121,7 +121,7 @@ template void count_num_diagonal_elements( std::shared_ptr exec, const matrix::SparsityCsr *matrix, - size_type &num_diagonal_elements) + size_type *num_diagonal_elements) { auto num_rows = matrix->get_size()[0]; auto row_ptrs = matrix->get_const_row_ptrs(); @@ -134,7 +134,7 @@ void count_num_diagonal_elements( } } } - num_diagonal_elements = num_diag; + *num_diagonal_elements = num_diag; } GKO_INSTANTIATE_FOR_EACH_VALUE_AND_INDEX_TYPE( diff --git a/reference/test/matrix/sparsity_csr.cpp b/reference/test/matrix/sparsity_csr.cpp index c58d15d5337..f20e8ba3134 100644 --- a/reference/test/matrix/sparsity_csr.cpp +++ b/reference/test/matrix/sparsity_csr.cpp @@ -86,7 +86,7 @@ TEST_F(SparsityCsr, CanBeCreatedFromExistingCsrMatrix) auto mtx = Mtx::create(exec, std::move(csr_mtx)); - GKO_ASSERT_MTX_NEAR(comp_mtx.get(), mtx.get(), 0.0) + GKO_ASSERT_MTX_NEAR(comp_mtx.get(), mtx.get(), 0.0); } @@ -99,7 +99,7 @@ TEST_F(SparsityCsr, CanBeCreatedFromExistingDenseMatrix) auto mtx = Mtx::create(exec, std::move(dense_mtx)); - GKO_ASSERT_MTX_NEAR(comp_mtx.get(), mtx.get(), 0.0) + GKO_ASSERT_MTX_NEAR(comp_mtx.get(), mtx.get(), 0.0); } diff --git a/reference/test/matrix/sparsity_csr_kernels.cpp b/reference/test/matrix/sparsity_csr_kernels.cpp index 4af0f932880..e63bca5453e 100644 --- a/reference/test/matrix/sparsity_csr_kernels.cpp +++ b/reference/test/matrix/sparsity_csr_kernels.cpp @@ -288,9 +288,9 @@ TEST_F(SparsityCsr, CountsCorrectNumberOfDiagonalElements) gko::size_type ms_num_diags = 0; gko::kernels::reference::sparsity_csr::count_num_diagonal_elements( - exec, mtx2.get(), m2_num_diags); + exec, mtx2.get(), &m2_num_diags); gko::kernels::reference::sparsity_csr::count_num_diagonal_elements( - exec, mtx_s.get(), ms_num_diags); + exec, mtx_s.get(), &ms_num_diags); ASSERT_EQ(m2_num_diags, 3); ASSERT_EQ(ms_num_diags, 2);