Skip to content

Fix scalar jacobi generate on different precision and add DiagonalLinOpExtractable #834

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
yhmtsai merged 6 commits into from
Aug 4, 2021
Merged

Fix scalar jacobi generate on different precision and add DiagonalLinOpExtractable #834

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
6 commits
Select commit Hold shift + click to select a range
47312f1
scalar_jacobi on diff precision csr
yhmtsai Jul 14, 2021
d9392ee
scalarjacobi allows diff system matrix and add DiagonalLinOpExtractable
yhmtsai Jul 14, 2021
76b3881
move diagonalextrable impl with instantiate
yhmtsai Jul 16, 2021
824e7e4
fix executor issue
yhmtsai Jul 19, 2021
914d0f4
use template_conversion, add diagonal precision convert
yhmtsai Aug 2, 2021
28e8097
review update
yhmtsai Aug 3, 2021
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
33 changes: 33 additions & 0 deletions core/matrix/diagonal.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -146,6 +146,22 @@ std::unique_ptr<LinOp> Diagonal<ValueType>::conj_transpose() const
}


template <typename ValueType>
void Diagonal<ValueType>::convert_to(
Diagonal<next_precision<ValueType>> *result) const
{
result->values_ = this->values_;
result->set_size(this->get_size());
}


template <typename ValueType>
void Diagonal<ValueType>::move_to(Diagonal<next_precision<ValueType>> *result)
{
this->convert_to(result);
}


template <typename ValueType>
void Diagonal<ValueType>::convert_to(Csr<ValueType, int32> *result) const
{
Expand Down Expand Up @@ -302,4 +318,21 @@ GKO_INSTANTIATE_FOR_EACH_VALUE_TYPE(GKO_DECLARE_DIAGONAL_MATRIX);


} // namespace matrix


// Implement DiagonalExtractable for LinOp when Diagonal is complete class
template <typename ValueType>
std::unique_ptr<LinOp> DiagonalExtractable<ValueType>::extract_diagonal_linop()
const
{
return this->extract_diagonal();
}


#define GKO_DECLARE_DIAGONAL_EXTRACTABLE(value_type) \
std::unique_ptr<LinOp> \
DiagonalExtractable<value_type>::extract_diagonal_linop() const
GKO_INSTANTIATE_FOR_EACH_VALUE_TYPE(GKO_DECLARE_DIAGONAL_EXTRACTABLE);


} // namespace gko
18 changes: 14 additions & 4 deletions core/preconditioner/jacobi.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -40,6 +40,7 @@ OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
#include <ginkgo/core/base/executor.hpp>
#include <ginkgo/core/base/math.hpp>
#include <ginkgo/core/base/precision_dispatch.hpp>
#include <ginkgo/core/base/temporary_conversion.hpp>
#include <ginkgo/core/base/utils.hpp>
#include <ginkgo/core/matrix/csr.hpp>
#include <ginkgo/core/matrix/dense.hpp>
Expand Down Expand Up @@ -231,11 +232,20 @@ void Jacobi<ValueType, IndexType>::generate(const LinOp *system_matrix,
using csr_type = matrix::Csr<ValueType, IndexType>;
const auto exec = this->get_executor();
if (parameters_.max_block_size == 1) {
auto diag = as<DiagonalExtractable<ValueType>>(system_matrix)
->extract_diagonal();
auto temp = gko::Array<ValueType>::view(
exec, system_matrix->get_size()[0], diag->get_values());
auto diag = share(as<DiagonalLinOpExtractable>(system_matrix)
->extract_diagonal_linop());
Comment on lines +235 to +236

@pratikvn pratikvn Jul 19, 2021

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is there some reason we need the extract_diagonal_linop and extract_diagonal will not work ?

@yhmtsai yhmtsai Jul 19, 2021

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

extract_diagonal require the valuetype, but extract_diagonal_linop does not. we can get the diag first and then check which precision to work.
I am thinking do the dynamic_cast to check the matrix type and use corresponding extract_diag.

@upsj upsj Jul 19, 2021

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

An alternative to DiagonalLinOpExtractable would be using

auto csr_mtx = convert_to_with_sorting<csr_type>(exec, system_matrix, /*skip_sorting=*/true);

from the non-scalar branch and extracting from that. That would also make for more consistent behavior between scalar and block Jacobi.

@yhmtsai yhmtsai Jul 19, 2021

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I would like to avoid the converting matrix self because we need to pay more for it against converting diag

@pratikvn pratikvn Jul 19, 2021

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Could we also not convert the extracted diagonal after calling the extract_diagonal. Currently, Diagonal does not have mixed precision convert methods, but we could easily add those

@upsj upsj Jul 19, 2021

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Oh yes, I think I like that idea best of all suggestions yet :)

@yhmtsai yhmtsai Jul 19, 2021

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

but the issue is different, right?
you still need to dynamic_cast to the correct type first.
and it is the same impl in the current one if we use default behavior
if we use direct mapping, we need to add all matrix to support it (at least csr)

@pratikvn pratikvn Jul 20, 2021

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Actually, i am a bit confused as to how the current implementation works. extract_diagonal_linop should return a Diagonal of ValueType right ? Because it is just calling the extract_diagonal of Csr<ValueType, IndexType>. So, the else if branch in line 243, will never be taken ?

@yhmtsai yhmtsai Aug 2, 2021

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

extract_diagonal_linop returns LinOp (it is Diaonal<system_matrix ValueType>), the valuetype is from system matrix not from the jacobi.
I realize we miss the precision convert in Diagonal, but I guess I can not use it directly here?

@upsj upsj Aug 2, 2021

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You could use something similar to make_temporary_conversion. The framework also supports other matrix types, you only need to wrap it like you see in that function. Precision conversions in Diagonal is definitely a good idea, I forgot to suggest those when we added the format.

auto diag_vt =
::gko::detail::temporary_conversion<matrix::Diagonal<ValueType>>::

@upsj upsj Aug 2, 2021

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

same here: Is the fully qualified name necessary?

Suggested change
::gko::detail::temporary_conversion<matrix::Diagonal<ValueType>>::
detail::temporary_conversion<matrix::Diagonal<ValueType>>::

@yhmtsai yhmtsai Aug 3, 2021

Copy link
Copy Markdown
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I need to put gko::, or I will get the error
error: 'temporary_conversion' is not a member of 'gko::preconditioner::detail'
it look for the function in gko::preconditioner::detail not gko::detail

template create<matrix::Diagonal<next_precision<ValueType>>>(
diag.get());
if (!diag_vt) {
GKO_NOT_SUPPORTED(system_matrix);
}
auto temp = Array<ValueType>::view(diag_vt->get_executor(),
diag_vt->get_size()[0],
diag_vt->get_values());
this->blocks_ = temp;
this->num_blocks_ = diag_vt->get_size()[0];
} else {
auto csr_mtx = convert_to_with_sorting<csr_type>(exec, system_matrix,
skip_sorting);
Expand Down
27 changes: 26 additions & 1 deletion include/ginkgo/core/base/lin_op.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -654,6 +654,29 @@ class Preconditionable {
};


/**
* The diagonal of a LinOp can be extracted. It will be implemented by
* DiagonalExtractable<ValueType>, so the class does not need to implement it.
* extract_diagonal_linop returns a linop which extracts the elements whose col
* and row index are the same and stores the result in a min(nrows, ncols) x 1
* dense matrix.
*
* @ingroup diagonal
* @ingroup LinOp
Comment thread
yhmtsai marked this conversation as resolved.
*/
class DiagonalLinOpExtractable {
public:
virtual ~DiagonalLinOpExtractable() = default;

/**
* Extracts the diagonal entries of the matrix into a vector.
*
* @return linop the linop of diagonal format
*/
virtual std::unique_ptr<LinOp> extract_diagonal_linop() const = 0;
};


/**
* The diagonal of a LinOp implementing this interface can be extracted.
* extract_diagonal extracts the elements whose col and row index are the
Expand All @@ -662,12 +685,14 @@ class Preconditionable {
* @ingroup LinOp
*/
template <typename ValueType>
class DiagonalExtractable {
class DiagonalExtractable : public DiagonalLinOpExtractable {
public:
using value_type = ValueType;

virtual ~DiagonalExtractable() = default;

std::unique_ptr<LinOp> extract_diagonal_linop() const override;

/**
* Extracts the diagonal entries of the matrix into a vector.
*
Expand Down
7 changes: 7 additions & 0 deletions include/ginkgo/core/matrix/diagonal.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -70,6 +70,7 @@ class Diagonal
public EnableCreateMethod<Diagonal<ValueType>>,
public ConvertibleTo<Csr<ValueType, int32>>,
public ConvertibleTo<Csr<ValueType, int64>>,
public ConvertibleTo<Diagonal<next_precision<ValueType>>>,
public Transposable,
public WritableToMatrixData<ValueType, int32>,
public WritableToMatrixData<ValueType, int64>,
Expand All @@ -92,10 +93,16 @@ class Diagonal
using mat_data32 = gko::matrix_data<ValueType, int32>;
using absolute_type = remove_complex<Diagonal>;

friend class Diagonal<next_precision<ValueType>>;

std::unique_ptr<LinOp> transpose() const override;

std::unique_ptr<LinOp> conj_transpose() const override;

void convert_to(Diagonal<next_precision<ValueType>> *result) const override;

void move_to(Diagonal<next_precision<ValueType>> *result) override;

void convert_to(Csr<ValueType, int32> *result) const override;

void move_to(Csr<ValueType, int32> *result) override;
Expand Down
40 changes: 40 additions & 0 deletions reference/test/matrix/diagonal_kernels.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -108,6 +108,46 @@ class Diagonal : public ::testing::Test {
TYPED_TEST_SUITE(Diagonal, gko::test::ValueTypes);


TYPED_TEST(Diagonal, ConvertsToPrecision)
{
using ValueType = typename TestFixture::value_type;
using OtherType = typename gko::next_precision<ValueType>;
using Diagonal = typename TestFixture::Diag;
using OtherDiagonal = gko::matrix::Diagonal<OtherType>;
auto tmp = OtherDiagonal::create(this->exec);
auto res = Diagonal::create(this->exec);
// If OtherType is more precise: 0, otherwise r
auto residual = r<OtherType>::value < r<ValueType>::value
? gko::remove_complex<ValueType>{0}
: gko::remove_complex<ValueType>{r<OtherType>::value};

this->diag1->convert_to(tmp.get());
tmp->convert_to(res.get());

GKO_ASSERT_MTX_NEAR(this->diag1, res, residual);
}


TYPED_TEST(Diagonal, MovesToPrecision)
{
using ValueType = typename TestFixture::value_type;
using OtherType = typename gko::next_precision<ValueType>;
using Diagonal = typename TestFixture::Diag;
using OtherDiagonal = gko::matrix::Diagonal<OtherType>;
auto tmp = OtherDiagonal::create(this->exec);
auto res = Diagonal::create(this->exec);
// If OtherType is more precise: 0, otherwise r
auto residual = r<OtherType>::value < r<ValueType>::value
? gko::remove_complex<ValueType>{0}
: gko::remove_complex<ValueType>{r<OtherType>::value};

this->diag1->move_to(tmp.get());
tmp->move_to(res.get());

GKO_ASSERT_MTX_NEAR(this->diag1, res, residual);
}


TYPED_TEST(Diagonal, AppliesToDense)
{
using value_type = typename TestFixture::value_type;
Expand Down
23 changes: 23 additions & 0 deletions reference/test/preconditioner/jacobi.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -65,6 +65,7 @@ class Jacobi : public ::testing::Test {
Jacobi()
: exec(gko::ReferenceExecutor::create()),
bj_factory(Bj::build().with_max_block_size(3u).on(exec)),
scalar_j_factory(Bj::build().with_max_block_size(1u).on(exec)),
block_pointers(exec, 3),
block_precisions(exec, 2),
mtx(Mtx::create(exec, gko::dim<2>{5}, 13))
Expand Down Expand Up @@ -168,6 +169,7 @@ class Jacobi : public ::testing::Test {

std::shared_ptr<const gko::Executor> exec;
std::unique_ptr<typename Bj::Factory> bj_factory;
std::unique_ptr<typename Bj::Factory> scalar_j_factory;
std::unique_ptr<typename Bj::Factory> adaptive_bj_factory;
gko::Array<index_type> block_pointers;
gko::Array<gko::precision_reduction> block_precisions;
Expand Down Expand Up @@ -378,4 +380,25 @@ TYPED_TEST(Jacobi, GeneratesCorrectMatrixDataWithAdaptivePrecision)
}


TYPED_TEST(Jacobi, ScalarJacobiGeneratesOnDifferentPrecision)
{
using value_type = typename TestFixture::value_type;
using index_type = typename TestFixture::index_type;
using next_type = gko::next_precision<value_type>;
using Bj = typename TestFixture::Bj;
auto csr =
gko::share(gko::matrix::Csr<next_type, index_type>::create(this->exec));
csr->copy_from(gko::lend(this->mtx));
std::shared_ptr<Bj> bj{};

ASSERT_NO_THROW(bj = this->scalar_j_factory->generate(csr));
ASSERT_EQ(bj->get_num_blocks(), 5u);
ASSERT_EQ(bj->get_blocks()[0], value_type{4});
ASSERT_EQ(bj->get_blocks()[1], value_type{4});
ASSERT_EQ(bj->get_blocks()[2], value_type{4});
ASSERT_EQ(bj->get_blocks()[3], value_type{4});
ASSERT_EQ(bj->get_blocks()[4], value_type{4});
}


} // namespace