Skip to content

kokkos: fix reversed triangles in trmm, trsm and trsv #227

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
crtrott merged 2 commits into from
Jul 26, 2022
Merged

kokkos: fix reversed triangles in trmm, trsm and trsv #227

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
b601bbd
kokkos: fix reversed triangles in TRSM and TRSV
May 27, 2022
7784cee
kokkos: fix reversed triangles in TRMM #226
May 27, 2022
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
6 changes: 3 additions & 3 deletions tests/kokkos-based/triangular_matrix_left_product_kokkos.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -65,11 +65,11 @@ void updating_triangular_matrix_left_product_gold_solution(A_t A,
// because (i,j) indexing respects C updating order
// and parallelism is restricted accordingly.
for (size_type ii = 0; ii < C_ext0; ++ii) {
const size_type i = lower ? ii : C_ext0 - 1 - ii;
const size_type i = lower ? C_ext0 - 1 - ii : ii;
for (size_type j = 0; j < C_ext1; ++j) {
c_element_type t = E(i, j);
const size_type k0 = lower ? (explicit_diag ? i : i + 1) : 0;
const size_type k1 = lower ? C_ext0 :(explicit_diag ? i + 1 : i);
const size_type k0 = lower ? 0 : (explicit_diag ? i : i + 1);
const size_type k1 = lower ? (explicit_diag ? i + 1 : i) : C_ext0;
for (size_type k = k0; k < k1; ++k) {
t += A(i, k) * B(k, j);
}
Expand Down
6 changes: 3 additions & 3 deletions tests/kokkos-based/triangular_matrix_matrix_left_solve.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -66,10 +66,10 @@ void triangular_matrix_matrix_left_solve_gold_solution(

for (size_type k = 0; k < num_vectors; ++k) {
for (size_type ii = 0; ii < A_ext; ++ii) {
const size_type i = lower_triangle ? A_ext - 1 - ii : ii;
const size_type i = lower_triangle ? ii : A_ext - 1 - ii;
// A(i, j) has lower triangle in i <= j
const size_type j0 = lower_triangle ? i + 1 : 0;
const size_type j1 = lower_triangle ? A_ext : i;
const size_type j0 = lower_triangle ? 0 : i + 1;
const size_type j1 = lower_triangle ? i : A_ext;
for (size_type j = j0; j < j1; ++j) {
X(i, k) -= A(i, j) * X(j, k);
}
Expand Down
6 changes: 3 additions & 3 deletions tests/kokkos-based/triangular_matrix_matrix_right_solve.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -66,10 +66,10 @@ void triangular_matrix_matrix_right_solve_gold_solution(

for (size_type i = 0; i < num_vectors; ++i) {
for (size_type kk = 0; kk < A_ext; ++kk) {
const size_type k = lower_triangle ? kk : A_ext - 1 - kk;
const size_type k = lower_triangle ? A_ext - 1 - kk : kk;
// A(j, k) has lower triangle in j <= k
const size_type j0 = lower_triangle ? 0 : k + 1;
const size_type j1 = lower_triangle ? k : A_ext;
const size_type j0 = lower_triangle ? k + 1 : 0;
const size_type j1 = lower_triangle ? A_ext : k;
for (size_type j = j0; j < j1; ++j) {
X(i, k) -= X(i, j) * A(j, k);
}
Expand Down
6 changes: 3 additions & 3 deletions tests/kokkos-based/triangular_matrix_right_product_kokkos.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -65,12 +65,12 @@ void updating_triangular_matrix_right_product_gold_solution(A_t A,
// because (i,j) indexing respects C updating order
// and parallelism is restricted accordingly.
for (size_type jj = 0; jj < C_ext1; ++jj) {
const size_type j = lower ? C_ext1 - 1 - jj : jj;
const size_type j = lower ? jj : C_ext1 - 1 - jj;
for (size_type i = 0; i < C_ext0; ++i) {
c_element_type t = E(i, j);
// Note: lower triangle of A(k, j) means k <= j
const auto k0 = lower ? 0 : (explicit_diag ? j : j + 1);
const auto k1 = lower ? (explicit_diag ? j + 1 : j) : C.extent(1);
const auto k0 = lower ? (explicit_diag ? j : j + 1) : 0;
const auto k1 = lower ? C_ext1 : (explicit_diag ? j + 1 : j);
for (size_type k = k0; k < k1; ++k) {
t += B(i, k) * A(k, j);
}
Expand Down
6 changes: 3 additions & 3 deletions tests/kokkos-based/triangular_matrix_vector_solve.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -64,9 +64,9 @@ void triangular_matrix_vector_solve_gold_solution(A_t A, Triangle t, DiagonalSto
const size_type ext1 = A.extent(1);

for (size_type ii = 0; ii < ext0; ++ii) {
const size_type i = lower_triangle ? ext0 - 1 - ii : ii;
const size_type j0 = lower_triangle ? i + 1 : 0;
const size_type j1 = lower_triangle ? ext1 : i;
const size_type i = lower_triangle ? ii : ext0 - 1 - ii;
const size_type j0 = lower_triangle ? 0 : i + 1;
const size_type j1 = lower_triangle ? i : ext1;
for (size_type j = j0; j < j1; ++j) {
x(i) -= A(i, j) * x(j);
}
Expand Down
12 changes: 6 additions & 6 deletions ...lementations/include/experimental/__p1673_bits/kokkos-kernels/blas3_matrix_product_kk.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -261,12 +261,12 @@ void trmm_left(KokkosExecSpace &&exec, Triangle t, DiagonalStorage d,
// because (i,j) indexing respects C updating order
// and parallelism is restricted accordingly
for (size_type ii = 0; ii < C_ext0; ++ii) {
const auto i = lower ? ii : C_ext0 - 1 - ii;
const auto i = lower ? C_ext0 - 1 - ii : ii;
Kokkos::parallel_for(Kokkos::RangePolicy(exec, 0, C_ext1),
KOKKOS_LAMBDA (const size_type j) {
c_element_type t{};
const size_type k0 = lower ? (explicit_diag ? i : i + 1) : 0;
const size_type k1 = lower ? C_ext0 :(explicit_diag ? i + 1 : i);
const size_type k0 = lower ? 0 : (explicit_diag ? i : i + 1);
const size_type k1 = lower ? (explicit_diag ? i + 1 : i) : C_ext0;
for (size_type k = k0; k < k1; ++k) {
t += A_view(i, k) * B_view(k, j);
}
Expand Down Expand Up @@ -300,13 +300,13 @@ void trmm_right(KokkosExecSpace &&exec, Triangle t, DiagonalStorage d,
// because (i,j) indexing respects C updating order
// and parallelism is restricted accordingly
for (size_type jj = 0; jj < C_ext1; ++jj) {
const size_type j = lower ? C_ext1 - 1 - jj : jj;
const size_type j = lower ? jj : C_ext1 - 1 - jj;
Kokkos::parallel_for(Kokkos::RangePolicy(exec, 0, C_ext0),
KOKKOS_LAMBDA (const size_type i) {
c_element_type t{};
// Note: lower triangle of A(k, j) means k <= j
const auto k0 = lower ? 0 : (explicit_diag ? j : j + 1);
const auto k1 = lower ? (explicit_diag ? j + 1 : j) : C_ext1;
const auto k0 = lower ? (explicit_diag ? j : j + 1) : 0;
const auto k1 = lower ? C_ext1 : (explicit_diag ? j + 1 : j);
for (size_type k = k0; k < k1; ++k) {
t += B_view(i, k) * A_view(k, j);
}
Expand Down
3 changes: 1 addition & 2 deletions ...include/experimental/__p1673_bits/kokkos-kernels/blas3_triangular_matrix_matrix_solve.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -63,9 +63,8 @@ void trsm(Side /*s*/, Triangle /*t*/, DiagonalStorage /*d*/, AViewType A, XViewT
{
const auto side = std::is_same_v<Side,
std::experimental::linalg::left_side_t> ? "L" : "R";
// KK and stdBLAS use REVERSED triangle definitions
const auto uplo = std::is_same_v<Triangle,
std::experimental::linalg::lower_triangle_t> ? "U" : "L";
std::experimental::linalg::lower_triangle_t> ? "L" : "U";
const auto diag = std::is_same_v<DiagonalStorage,
std::experimental::linalg::explicit_diagonal_t> ? "N" : "U";

Expand Down