From 83dd0cc998e9f7d447c92aff1f0a2b3cf90a2d94 Mon Sep 17 00:00:00 2001 From: dcdemen Date: Mon, 23 Dec 2024 09:10:36 -0700 Subject: [PATCH 1/3] Split solver constructor into many functions in their own files to increase readability and testability As part of this, the types defined on Solver were simplified --- src/solver/compute_b_num_non_zero.hpp | 19 + src/solver/compute_k_col_inds.hpp | 87 +++ src/solver/compute_k_num_non_zero.hpp | 65 ++ src/solver/compute_k_row_ptrs.hpp | 89 +++ src/solver/compute_num_system_dofs.hpp | 30 + src/solver/compute_t_col_inds.hpp | 50 ++ src/solver/compute_t_num_non_zero.hpp | 31 + src/solver/compute_t_row_ptrs.hpp | 65 ++ src/solver/create_b_matrix.hpp | 47 ++ src/solver/create_b_t_matrix.hpp | 65 ++ src/solver/create_constraints_matrix_full.hpp | 28 + src/solver/create_global_matrix.hpp | 27 + src/solver/create_k_matrix.hpp | 52 ++ src/solver/create_matrix_spadd.hpp | 19 + src/solver/create_matrix_spgemm.hpp | 19 + src/solver/create_sparse_dense_solver.hpp | 27 + src/solver/create_system_matrix_full.hpp | 29 + src/solver/create_t_matrix.hpp | 47 ++ src/solver/create_transpose_matrix_full.hpp | 32 + src/solver/solver.hpp | 584 ++---------------- src/step/assemble_constraints_matrix.hpp | 2 +- src/step/assemble_constraints_residual.hpp | 11 +- tests/unit_tests/solver/CMakeLists.txt | 1 + 23 files changed, 885 insertions(+), 541 deletions(-) create mode 100644 src/solver/compute_b_num_non_zero.hpp create mode 100644 src/solver/compute_k_col_inds.hpp create mode 100644 src/solver/compute_k_num_non_zero.hpp create mode 100644 src/solver/compute_k_row_ptrs.hpp create mode 100644 src/solver/compute_num_system_dofs.hpp create mode 100644 src/solver/compute_t_col_inds.hpp create mode 100644 src/solver/compute_t_num_non_zero.hpp create mode 100644 src/solver/compute_t_row_ptrs.hpp create mode 100644 src/solver/create_b_matrix.hpp create mode 100644 src/solver/create_b_t_matrix.hpp create mode 100644 src/solver/create_constraints_matrix_full.hpp create mode 100644 src/solver/create_global_matrix.hpp create mode 100644 src/solver/create_k_matrix.hpp create mode 100644 src/solver/create_matrix_spadd.hpp create mode 100644 src/solver/create_matrix_spgemm.hpp create mode 100644 src/solver/create_sparse_dense_solver.hpp create mode 100644 src/solver/create_system_matrix_full.hpp create mode 100644 src/solver/create_t_matrix.hpp create mode 100644 src/solver/create_transpose_matrix_full.hpp diff --git a/src/solver/compute_b_num_non_zero.hpp b/src/solver/compute_b_num_non_zero.hpp new file mode 100644 index 00000000..0bd89b3c --- /dev/null +++ b/src/solver/compute_b_num_non_zero.hpp @@ -0,0 +1,19 @@ +#pragma once + +#include + +#include "compute_number_of_non_zeros_constraints.hpp" + +namespace openturbine { +[[nodiscard]] static size_t ComputeBNumNonZero( + const Kokkos::View::const_type& constraint_type, + const Kokkos::View*>::const_type& constraint_row_range +) { + auto B_num_non_zero = size_t{0U}; + Kokkos::parallel_reduce( + "ComputeNumberOfNonZeros_Constraints", constraint_type.extent(0), + ComputeNumberOfNonZeros_Constraints{constraint_type, constraint_row_range}, B_num_non_zero + ); + return B_num_non_zero; +} +} // namespace openturbine diff --git a/src/solver/compute_k_col_inds.hpp b/src/solver/compute_k_col_inds.hpp new file mode 100644 index 00000000..0248d4da --- /dev/null +++ b/src/solver/compute_k_col_inds.hpp @@ -0,0 +1,87 @@ +#pragma once + +#include + +#include "src/dof_management/freedom_signature.hpp" + +namespace openturbine { + +template +struct ComputeKColIndsFunction { + using RowPtrValueType = typename RowPtrType::value_type; + using IndicesValueType = typename IndicesType::value_type; + Kokkos::View::const_type node_freedom_allocation_table; + Kokkos::View::const_type node_freedom_map_table; + Kokkos::View::const_type num_nodes_per_element; + Kokkos::View::const_type node_state_indices; + typename RowPtrType::const_type K_row_ptrs; + IndicesType K_col_inds; + + KOKKOS_FUNCTION + bool ContainsNode(size_t element, size_t node) const { + bool contains_node = false; + for (auto n = 0U; n < num_nodes_per_element(element); ++n) { + contains_node = contains_node || (node_state_indices(element, n) == node); + } + return contains_node; + } + + KOKKOS_FUNCTION + void ComputeColInds(size_t element, size_t node, RowPtrValueType current_dof_index) const { + for (auto n = 0U; n < num_nodes_per_element(element); ++n) { + const auto node_state_index = node_state_indices(element, n); + if (node_state_index != node) { + const auto target_node_num_dof = + count_active_dofs(node_freedom_allocation_table(node_state_index)); + const auto target_node_dof_index = node_freedom_map_table(node_state_index); + for (auto k = 0U; k < target_node_num_dof; ++k, ++current_dof_index) { + const auto col_index = static_cast(target_node_dof_index + k); + K_col_inds(current_dof_index) = col_index; + } + } + } + } + + KOKKOS_FUNCTION + void operator()(size_t i) const { + const auto this_node_num_dof = count_active_dofs(node_freedom_allocation_table(i)); + const auto this_node_dof_index = node_freedom_map_table(i); + + for (auto j = 0U; j < this_node_num_dof; ++j) { + auto current_dof_index = K_row_ptrs(this_node_dof_index + j); + + for (auto k = 0U; k < this_node_num_dof; ++k, ++current_dof_index) { + K_col_inds(current_dof_index) = static_cast(this_node_dof_index + k); + } + + for (auto e = 0U; e < num_nodes_per_element.extent(0); ++e) { + if (!ContainsNode(e, i)) { + continue; + } + ComputeColInds(e, i, current_dof_index); + } + } + } +}; + +template +[[nodiscard]] inline IndicesType ComputeKColInds( + size_t K_num_non_zero, + const Kokkos::View::const_type& node_freedom_allocation_table, + const Kokkos::View::const_type& node_freedom_map_table, + const Kokkos::View::const_type& num_nodes_per_element, + const Kokkos::View::const_type& node_state_indices, + const typename RowPtrType::const_type& K_row_ptrs +) { + auto K_col_inds = IndicesType("col_inds", K_num_non_zero); + + Kokkos::parallel_for( + "ComputeKColInds", node_freedom_allocation_table.extent(0), + ComputeKColIndsFunction{ + node_freedom_allocation_table, node_freedom_map_table, num_nodes_per_element, + node_state_indices, K_row_ptrs, K_col_inds} + ); + + return K_col_inds; +} +} // namespace openturbine diff --git a/src/solver/compute_k_num_non_zero.hpp b/src/solver/compute_k_num_non_zero.hpp new file mode 100644 index 00000000..3bfdb485 --- /dev/null +++ b/src/solver/compute_k_num_non_zero.hpp @@ -0,0 +1,65 @@ +#pragma once + +#include + +#include "src/dof_management/freedom_signature.hpp" + +namespace openturbine { + +struct ComputeKNumNonZero_OffDiagonal { + Kokkos::View::const_type num_nodes_per_element; + Kokkos::View::const_type node_state_indices; + Kokkos::View::const_type node_freedom_allocation_table; + + KOKKOS_FUNCTION + void operator()(size_t i, size_t& update) const { + auto num_element_dof = 0UL; + for (auto j = 0U; j < num_nodes_per_element(i); ++j) { + const auto num_node_dof = + count_active_dofs(node_freedom_allocation_table(node_state_indices(i, j))); + num_element_dof += num_node_dof; + } + const auto num_element_non_zero = num_element_dof * num_element_dof; + auto num_diagonal_non_zero = 0UL; + for (auto j = 0U; j < num_nodes_per_element(i); ++j) { + const auto num_node_dof = + count_active_dofs(node_freedom_allocation_table(node_state_indices(i, j))); + num_diagonal_non_zero += num_node_dof * num_node_dof; + } + update += num_element_non_zero - num_diagonal_non_zero; + } +}; + +struct ComputeKNumNonZero_Diagonal { + Kokkos::View::const_type node_freedom_allocation_table; + + KOKKOS_FUNCTION + void operator()(size_t i, size_t& update) const { + const auto num_node_dof = count_active_dofs(node_freedom_allocation_table(i)); + const auto num_diagonal_non_zero = num_node_dof * num_node_dof; + update += num_diagonal_non_zero; + } +}; + +/// Computes the number of non-zero entries in the stiffness matrix K for sparse storage +[[nodiscard]] inline size_t ComputeKNumNonZero( + const Kokkos::View::const_type& num_nodes_per_element, + const Kokkos::View::const_type& node_state_indices, + const Kokkos::View::const_type& node_freedom_allocation_table +) { + auto K_num_non_zero_off_diagonal = 0UL; + Kokkos::parallel_reduce( + "ComputeKNumNonZero_OffDiagonal", num_nodes_per_element.extent(0), + ComputeKNumNonZero_OffDiagonal{ + num_nodes_per_element, node_state_indices, node_freedom_allocation_table}, + K_num_non_zero_off_diagonal + ); + auto K_num_non_zero_diagonal = 0UL; + Kokkos::parallel_reduce( + "ComputeKNumNonZero_Diagonal", node_freedom_allocation_table.extent(0), + ComputeKNumNonZero_Diagonal{node_freedom_allocation_table}, K_num_non_zero_diagonal + ); + return K_num_non_zero_off_diagonal + K_num_non_zero_diagonal; +} + +} // namespace openturbine diff --git a/src/solver/compute_k_row_ptrs.hpp b/src/solver/compute_k_row_ptrs.hpp new file mode 100644 index 00000000..3ed8ed7c --- /dev/null +++ b/src/solver/compute_k_row_ptrs.hpp @@ -0,0 +1,89 @@ +#pragma once + +#include + +#include "src/dof_management/freedom_signature.hpp" + +namespace openturbine { + +template +struct ComputeKRowEntries { + Kokkos::View::const_type node_freedom_allocation_table; + Kokkos::View::const_type node_freedom_map_table; + Kokkos::View::const_type num_nodes_per_element; + Kokkos::View::const_type node_state_indices; + RowPtrType K_row_entries; + + KOKKOS_FUNCTION + void operator()(size_t i) const { + const auto this_node_num_dof = count_active_dofs(node_freedom_allocation_table(i)); + const auto this_node_dof_index = node_freedom_map_table(i); + + auto num_entries_in_row = this_node_num_dof; + bool node_found_in_system = false; + + // contributions to non-diagonal block from coupled nodes + for (auto e = 0U; e < num_nodes_per_element.extent(0); ++e) { + bool contains_node = false; + auto num_entries_in_element = 0UL; + for (auto j = 0U; j < num_nodes_per_element(e); ++j) { + contains_node = contains_node || (node_state_indices(e, j) == i); + num_entries_in_element += + count_active_dofs(node_freedom_allocation_table(node_state_indices(e, j))); + } + if (contains_node) { + node_found_in_system = true; + num_entries_in_row += num_entries_in_element - this_node_num_dof; + } + } + if (node_found_in_system) { + for (auto j = 0U; j < this_node_num_dof; ++j) { + K_row_entries(this_node_dof_index + j) = num_entries_in_row; + } + } + } +}; + +template +struct ComputeKRowPtrsScanner { + typename RowPtrType::const_type K_row_entries; + RowPtrType K_row_ptrs; + + KOKKOS_FUNCTION + void operator()(size_t i, size_t& update, bool is_final) const { + update += K_row_entries(i); + if (is_final) { + K_row_ptrs(i + 1) = update; + } + } +}; + +/// Computes the row pointers for the sparse stiffness matrix K in CSR format +template +[[nodiscard]] inline RowPtrType ComputeKRowPtrs( + size_t K_num_rows, + const Kokkos::View::const_type& node_freedom_allocation_table, + const Kokkos::View::const_type& node_freedom_map_table, + const Kokkos::View::const_type& num_nodes_per_element, + const Kokkos::View::const_type& node_state_indices +) { + auto K_row_ptrs = RowPtrType("row_ptrs", K_num_rows + 1); + const auto K_row_entries = RowPtrType("row_ptrs", K_num_rows); + + Kokkos::parallel_for( + "ComputeKRowEntries", node_freedom_allocation_table.extent(0), + ComputeKRowEntries{ + node_freedom_allocation_table, node_freedom_map_table, num_nodes_per_element, + node_state_indices, K_row_entries} + ); + + auto result = 0UL; + Kokkos::parallel_scan( + "ComputeKRowPtrs", K_row_entries.extent(0), + ComputeKRowPtrsScanner{K_row_entries, K_row_ptrs}, result + ); + + return K_row_ptrs; +} + +} // namespace openturbine diff --git a/src/solver/compute_num_system_dofs.hpp b/src/solver/compute_num_system_dofs.hpp new file mode 100644 index 00000000..83d09071 --- /dev/null +++ b/src/solver/compute_num_system_dofs.hpp @@ -0,0 +1,30 @@ +#pragma once + +#include + +#include "src/dof_management/freedom_signature.hpp" + +namespace openturbine { + +struct ComputeNumSystemDofsReducer { + Kokkos::View::const_type node_freedom_allocation_table; + + KOKKOS_FUNCTION + void operator()(size_t i, size_t& update) const { + update += count_active_dofs(node_freedom_allocation_table(i)); + } +}; + +/// Computes the total number of active degrees of freedom in the system +[[nodiscard]] inline size_t ComputeNumSystemDofs( + const Kokkos::View::const_type& node_freedom_allocation_table +) { + auto total_system_dofs = 0UL; + Kokkos::parallel_reduce( + "ComputeNumSystemDofs", node_freedom_allocation_table.extent(0), + ComputeNumSystemDofsReducer{node_freedom_allocation_table}, total_system_dofs + ); + return total_system_dofs; +} + +} // namespace openturbine diff --git a/src/solver/compute_t_col_inds.hpp b/src/solver/compute_t_col_inds.hpp new file mode 100644 index 00000000..145e3361 --- /dev/null +++ b/src/solver/compute_t_col_inds.hpp @@ -0,0 +1,50 @@ +#pragma once + +#include + +#include "src/dof_management/freedom_signature.hpp" + +namespace openturbine { + +template +struct ComputeTColIndsFunction { + using IndicesValuesType = typename IndicesType::value_type; + Kokkos::View::const_type node_freedom_allocation_table; + Kokkos::View::const_type node_freedom_map_table; + typename RowPtrType::const_type T_row_ptrs; + IndicesType T_col_inds; + + KOKKOS_FUNCTION + void operator()(size_t i) const { + const auto this_node_num_dof = count_active_dofs(node_freedom_allocation_table(i)); + const auto this_node_dof_index = node_freedom_map_table(i); + + for (auto j = 0U; j < this_node_num_dof; ++j) { + for (auto k = 0UL, current_dof_index = T_row_ptrs(this_node_dof_index + j); + k < this_node_num_dof; ++k, ++current_dof_index) { + const auto col_index = static_cast(this_node_dof_index + k); + T_col_inds(current_dof_index) = col_index; + } + } + } +}; + +template +[[nodiscard]] inline IndicesType ComputeTColInds( + size_t T_num_non_zero, + const Kokkos::View::const_type& node_freedom_allocation_table, + const Kokkos::View::const_type& node_freedom_map_table, + const typename RowPtrType::const_type& T_row_ptrs +) { + auto T_col_inds = IndicesType("T_indices", T_num_non_zero); + + Kokkos::parallel_for( + "ComputeTColInds", node_freedom_allocation_table.extent(0), + ComputeTColIndsFunction{ + node_freedom_allocation_table, node_freedom_map_table, T_row_ptrs, T_col_inds} + ); + + return T_col_inds; +} + +} // namespace openturbine diff --git a/src/solver/compute_t_num_non_zero.hpp b/src/solver/compute_t_num_non_zero.hpp new file mode 100644 index 00000000..dec4aeb5 --- /dev/null +++ b/src/solver/compute_t_num_non_zero.hpp @@ -0,0 +1,31 @@ +#pragma once + +#include + +#include "src/dof_management/freedom_signature.hpp" + +namespace openturbine { + +struct ComputeTNumNonZeroReducer { + Kokkos::View::const_type node_freedom_allocation_table; + + KOKKOS_FUNCTION + void operator()(size_t i, size_t& update) const { + const auto num_node_dof = count_active_dofs(node_freedom_allocation_table(i)); + const auto num_diagonal_non_zero = num_node_dof * num_node_dof; + update += num_diagonal_non_zero; + } +}; + +[[nodiscard]] inline size_t ComputeTNumNonZero( + const Kokkos::View::const_type& node_freedom_allocation_table +) { + auto T_num_non_zero = 0UL; + Kokkos::parallel_reduce( + "ComputeTNumNonZero", node_freedom_allocation_table.extent(0), + ComputeTNumNonZeroReducer{node_freedom_allocation_table}, T_num_non_zero + ); + + return T_num_non_zero; +} +} // namespace openturbine diff --git a/src/solver/compute_t_row_ptrs.hpp b/src/solver/compute_t_row_ptrs.hpp new file mode 100644 index 00000000..9fa3e085 --- /dev/null +++ b/src/solver/compute_t_row_ptrs.hpp @@ -0,0 +1,65 @@ +#pragma once + +#include + +#include "src/dof_management/freedom_signature.hpp" + +namespace openturbine { + +template +struct ComputeTRowEntries { + Kokkos::View::const_type node_freedom_allocation_table; + Kokkos::View::const_type node_freedom_map_table; + RowPtrType T_row_entries; + + KOKKOS_FUNCTION + void operator()(size_t i) const { + const auto this_node_num_dof = count_active_dofs(node_freedom_allocation_table(i)); + const auto this_node_dof_index = node_freedom_map_table(i); + + auto num_entries_in_row = this_node_num_dof; + + for (auto j = 0U; j < this_node_num_dof; ++j) { + T_row_entries(this_node_dof_index + j) = num_entries_in_row; + } + } +}; + +template +struct ComputeTRowPtrsScanner { + typename RowPtrType::const_type T_row_entries; + RowPtrType T_row_ptrs; + + KOKKOS_FUNCTION + void operator()(size_t i, size_t& update, bool is_final) const { + update += T_row_entries(i); + if (is_final) { + T_row_ptrs(i + 1) = update; + } + } +}; + +template +[[nodiscard]] inline RowPtrType ComputeTRowPtrs( + size_t T_num_rows, + const Kokkos::View::const_type& node_freedom_allocation_table, + const Kokkos::View::const_type& node_freedom_map_table +) { + auto T_row_ptrs = RowPtrType("T_row_ptrs", T_num_rows + 1); + const auto T_row_entries = RowPtrType("row_entries", T_num_rows); + + Kokkos::parallel_for( + "ComputeTRowEntries", node_freedom_allocation_table.extent(0), + ComputeTRowEntries{ + node_freedom_allocation_table, node_freedom_map_table, T_row_entries} + ); + + auto result = 0UL; + Kokkos::parallel_scan( + "ComputeTRowPtrs", T_row_entries.extent(0), + ComputeTRowPtrsScanner{T_row_entries, T_row_ptrs}, result + ); + + return T_row_ptrs; +} +} // namespace openturbine diff --git a/src/solver/create_b_matrix.hpp b/src/solver/create_b_matrix.hpp new file mode 100644 index 00000000..ebf883fc --- /dev/null +++ b/src/solver/create_b_matrix.hpp @@ -0,0 +1,47 @@ +#pragma once + +#include +#include + +#include "compute_b_num_non_zero.hpp" +#include "populate_sparse_row_ptrs_col_inds_constraints.hpp" + +namespace openturbine { + +/// Creates the constraint matrix B in sparse CRS format +template +[[nodiscard]] inline CrsMatrixType CreateBMatrix( + size_t system_dofs, size_t constraint_dofs, + const Kokkos::View::const_type& constraint_type, + const Kokkos::View::const_type& constraint_base_node_freedom_table, + const Kokkos::View::const_type& constraint_target_node_freedom_table, + const Kokkos::View*>::const_type& constraint_row_range +) { + using ValuesType = typename CrsMatrixType::values_type::non_const_type; + using RowPtrType = typename CrsMatrixType::staticcrsgraph_type::row_map_type::non_const_type; + using IndicesType = typename CrsMatrixType::staticcrsgraph_type::entries_type::non_const_type; + + const auto B_num_rows = constraint_dofs; + const auto B_num_columns = system_dofs; + const auto B_num_non_zero = ComputeBNumNonZero(constraint_type, constraint_row_range); + + const auto B_row_ptrs = RowPtrType("b_row_ptrs", B_num_rows + 1); + const auto B_col_ind = IndicesType("b_indices", B_num_non_zero); + Kokkos::parallel_for( + "PopulateSparseRowPtrsColInds_Constraints", 1, + PopulateSparseRowPtrsColInds_Constraints{ + constraint_type, constraint_base_node_freedom_table, + constraint_target_node_freedom_table, constraint_row_range, B_row_ptrs, B_col_ind} + ); + const auto B_values = ValuesType("B values", B_num_non_zero); + KokkosSparse::sort_crs_matrix(B_row_ptrs, B_col_ind, B_values); + return { + "B", + static_cast(B_num_rows), + static_cast(B_num_columns), + B_num_non_zero, + B_values, + B_row_ptrs, + B_col_ind}; +} +} // namespace openturbine diff --git a/src/solver/create_b_t_matrix.hpp b/src/solver/create_b_t_matrix.hpp new file mode 100644 index 00000000..373d3c28 --- /dev/null +++ b/src/solver/create_b_t_matrix.hpp @@ -0,0 +1,65 @@ +#pragma once + +#include +#include + +#include "compute_b_num_non_zero.hpp" +#include "populate_sparse_row_ptrs_col_inds_constraints.hpp" +#include "populate_sparse_row_ptrs_col_inds_transpose.hpp" + +namespace openturbine { + +/// Creates the constraint matrix B_t in sparse CRS format +template +[[nodiscard]] inline CrsMatrixType CreateBtMatrix( + size_t system_dofs, size_t constraint_dofs, + const Kokkos::View::const_type& constraint_type, + const Kokkos::View::const_type& constraint_base_node_freedom_table, + const Kokkos::View::const_type& constraint_target_node_freedom_table, + const Kokkos::View*>::const_type& constraint_row_range +) { + using ValuesType = typename CrsMatrixType::values_type::non_const_type; + using RowPtrType = typename CrsMatrixType::staticcrsgraph_type::row_map_type::non_const_type; + using IndicesType = typename CrsMatrixType::staticcrsgraph_type::entries_type::non_const_type; + + const auto B_num_rows = constraint_dofs; + const auto B_num_columns = system_dofs; + const auto B_num_non_zero = ComputeBNumNonZero(constraint_type, constraint_row_range); + + const auto B_row_ptrs = RowPtrType("b_row_ptrs", B_num_rows + 1); + const auto B_col_ind = IndicesType("b_indices", B_num_non_zero); + Kokkos::parallel_for( + "PopulateSparseRowPtrsColInds_Constraints", 1, + PopulateSparseRowPtrsColInds_Constraints{ + constraint_type, constraint_base_node_freedom_table, + constraint_target_node_freedom_table, constraint_row_range, B_row_ptrs, B_col_ind} + ); + const auto B_values = ValuesType("B values", B_num_non_zero); + KokkosSparse::sort_crs_matrix(B_row_ptrs, B_col_ind, B_values); + + const auto B_t_num_rows = system_dofs; + const auto B_t_num_columns = constraint_dofs; + const auto B_t_num_non_zero = ComputeBNumNonZero(constraint_type, constraint_row_range); + + auto col_count = IndicesType("col_count", B_num_columns); + auto tmp_row_ptrs = RowPtrType("tmp_row_ptrs", B_t_num_rows + 1); + auto B_t_row_ptrs = RowPtrType("b_t_row_ptrs", B_t_num_rows + 1); + auto B_t_col_inds = IndicesType("B_t_indices", B_t_num_non_zero); + auto B_t_values = ValuesType("B_t values", B_t_num_non_zero); + Kokkos::parallel_for( + "PopulateSparseRowPtrsColInds_Transpose", 1, + PopulateSparseRowPtrsColInds_Transpose{ + B_num_rows, B_num_columns, B_row_ptrs, B_col_ind, col_count, tmp_row_ptrs, B_t_row_ptrs, + B_t_col_inds} + ); + KokkosSparse::sort_crs_matrix(B_t_row_ptrs, B_t_col_inds, B_t_values); + return { + "B_t", + static_cast(B_t_num_rows), + static_cast(B_t_num_columns), + B_t_num_non_zero, + B_t_values, + B_t_row_ptrs, + B_t_col_inds}; +} +} // namespace openturbine diff --git a/src/solver/create_constraints_matrix_full.hpp b/src/solver/create_constraints_matrix_full.hpp new file mode 100644 index 00000000..f80a4b21 --- /dev/null +++ b/src/solver/create_constraints_matrix_full.hpp @@ -0,0 +1,28 @@ +#pragma once + +#include +#include + +namespace openturbine { + +template +[[nodiscard]] inline CrsMatrixType CreateConstraintsMatrixFull( + size_t num_system_dofs, size_t num_dofs, const CrsMatrixType& constraints_matrix +) { + using RowPtrType = typename CrsMatrixType::staticcrsgraph_type::row_map_type::non_const_type; + auto constraints_matrix_full_row_ptrs = + RowPtrType("constraints_matrix_full_row_ptrs", num_dofs + 1); + Kokkos::deep_copy( + Kokkos::subview( + constraints_matrix_full_row_ptrs, Kokkos::pair(num_system_dofs, num_dofs + 1) + ), + constraints_matrix.graph.row_map + ); + return CrsMatrixType( + "constraints_matrix_full", static_cast(num_dofs), static_cast(num_dofs), + constraints_matrix.nnz(), constraints_matrix.values, constraints_matrix_full_row_ptrs, + constraints_matrix.graph.entries + ); +} + +} // namespace openturbine diff --git a/src/solver/create_global_matrix.hpp b/src/solver/create_global_matrix.hpp new file mode 100644 index 00000000..189770ad --- /dev/null +++ b/src/solver/create_global_matrix.hpp @@ -0,0 +1,27 @@ +#pragma once + +#include +#include +#include + +namespace openturbine { + +template +[[nodiscard]] inline Teuchos::RCP CreateGlobalMatrix( + const typename GlobalCrsMatrixType::local_matrix_device_type& full_matrix +) { + using CrsMatrixType = typename GlobalCrsMatrixType::local_matrix_device_type; + using LocalOrdinalType = typename GlobalCrsMatrixType::local_ordinal_type; + using GlobalOrdinalType = typename GlobalCrsMatrixType::global_ordinal_type; + auto comm = Teuchos::createSerialComm(); + auto row_map = Tpetra::createLocalMap( + static_cast(full_matrix.numRows()), comm + ); + auto col_map = Tpetra::createLocalMap( + static_cast(full_matrix.numCols()), comm + ); + + return Teuchos::make_rcp(row_map, col_map, CrsMatrixType("A", full_matrix)); +} + +} // namespace openturbine diff --git a/src/solver/create_k_matrix.hpp b/src/solver/create_k_matrix.hpp new file mode 100644 index 00000000..ad97bdb0 --- /dev/null +++ b/src/solver/create_k_matrix.hpp @@ -0,0 +1,52 @@ +#pragma once + +#include +#include + +#include "compute_k_col_inds.hpp" +#include "compute_k_num_non_zero.hpp" +#include "compute_k_row_ptrs.hpp" + +#include "src/dof_management/freedom_signature.hpp" + +namespace openturbine { + +/// Creates the system stiffness matrix K in sparse CRS format +template +[[nodiscard]] inline CrsMatrixType CreateKMatrix( + size_t system_dofs, + const Kokkos::View::const_type& node_freedom_allocation_table, + const Kokkos::View::const_type& node_freedom_map_table, + const Kokkos::View::const_type& num_nodes_per_element, + const Kokkos::View::const_type& node_state_indices +) { + using ValuesType = typename CrsMatrixType::values_type::non_const_type; + using RowPtrType = typename CrsMatrixType::staticcrsgraph_type::row_map_type::non_const_type; + using IndicesType = typename CrsMatrixType::staticcrsgraph_type::entries_type::non_const_type; + + const auto K_num_rows = system_dofs; + const auto K_num_columns = K_num_rows; + const auto K_num_non_zero = + ComputeKNumNonZero(num_nodes_per_element, node_state_indices, node_freedom_allocation_table); + + const auto K_row_ptrs = ComputeKRowPtrs( + K_num_rows, node_freedom_allocation_table, node_freedom_map_table, num_nodes_per_element, + node_state_indices + ); + const auto K_col_inds = ComputeKColInds( + K_num_non_zero, node_freedom_allocation_table, node_freedom_map_table, num_nodes_per_element, + node_state_indices, K_row_ptrs + ); + const auto K_values = ValuesType("K values", K_num_non_zero); + + KokkosSparse::sort_crs_matrix(K_row_ptrs, K_col_inds, K_values); + return { + "K", + static_cast(K_num_rows), + static_cast(K_num_columns), + K_num_non_zero, + K_values, + K_row_ptrs, + K_col_inds}; +} +} // namespace openturbine diff --git a/src/solver/create_matrix_spadd.hpp b/src/solver/create_matrix_spadd.hpp new file mode 100644 index 00000000..dbcbac99 --- /dev/null +++ b/src/solver/create_matrix_spadd.hpp @@ -0,0 +1,19 @@ +#pragma once + +#include +#include + +namespace openturbine { + +template +[[nodiscard]] inline CrsMatrixType CreateMatrixSpadd( + const CrsMatrixType& A, const CrsMatrixType& B, KernelHandle& handle +) { + auto C = CrsMatrixType{}; + handle.create_spadd_handle(true, true); + KokkosSparse::spadd_symbolic(&handle, A, B, C); + KokkosSparse::spadd_numeric(&handle, 1., A, 1., B, C); + return C; +} + +} // namespace openturbine diff --git a/src/solver/create_matrix_spgemm.hpp b/src/solver/create_matrix_spgemm.hpp new file mode 100644 index 00000000..2704cad1 --- /dev/null +++ b/src/solver/create_matrix_spgemm.hpp @@ -0,0 +1,19 @@ +#pragma once + +#include +#include + +namespace openturbine { + +template +[[nodiscard]] inline CrsMatrixType CreateMatrixSpgemm( + const CrsMatrixType& A, const CrsMatrixType& B, KernelHandle& handle +) { + auto C = CrsMatrixType{}; + handle.create_spgemm_handle(); + KokkosSparse::spgemm_symbolic(handle, A, false, B, false, C); + KokkosSparse::spgemm_numeric(handle, A, false, B, false, C); + return C; +} + +} // namespace openturbine diff --git a/src/solver/create_sparse_dense_solver.hpp b/src/solver/create_sparse_dense_solver.hpp new file mode 100644 index 00000000..06dee656 --- /dev/null +++ b/src/solver/create_sparse_dense_solver.hpp @@ -0,0 +1,27 @@ +#pragma once + +#include + +#include + +namespace openturbine { + +template +[[nodiscard]] inline Teuchos::RCP> +CreateSparseDenseSolver( + Teuchos::RCP& A, Teuchos::RCP& x_mv, + Teuchos::RCP& b +) { + using ExecutionSpace = typename GlobalCrsMatrixType::execution_space; + + const auto solver_name = (std::is_same_v) + ? std::string{"klu2"} + : std::string{"basker"}; + + auto amesos_solver = + Amesos2::create(solver_name, A, x_mv, b); + amesos_solver->symbolicFactorization(); + return amesos_solver; +} + +} // namespace openturbine diff --git a/src/solver/create_system_matrix_full.hpp b/src/solver/create_system_matrix_full.hpp new file mode 100644 index 00000000..321d6a6b --- /dev/null +++ b/src/solver/create_system_matrix_full.hpp @@ -0,0 +1,29 @@ +#pragma once + +#include +#include + +#include "fill_unshifted_row_ptrs.hpp" + +namespace openturbine { + +template +[[nodiscard]] inline CrsMatrixType CreateSystemMatrixFull( + size_t num_system_dofs, size_t num_dofs, const CrsMatrixType& system_matrix +) { + using RowPtrType = typename CrsMatrixType::staticcrsgraph_type::row_map_type::non_const_type; + + auto system_matrix_full_row_ptrs = RowPtrType("system_matrix_full_row_ptrs", num_dofs + 1); + Kokkos::parallel_for( + "FillUnshiftedRowPtrs", num_dofs + 1, + FillUnshiftedRowPtrs{ + num_system_dofs, system_matrix.graph.row_map, system_matrix_full_row_ptrs} + ); + return CrsMatrixType( + "system_matrix_full", static_cast(num_dofs), static_cast(num_dofs), + system_matrix.nnz(), system_matrix.values, system_matrix_full_row_ptrs, + system_matrix.graph.entries + ); +} + +} // namespace openturbine diff --git a/src/solver/create_t_matrix.hpp b/src/solver/create_t_matrix.hpp new file mode 100644 index 00000000..acfda8b8 --- /dev/null +++ b/src/solver/create_t_matrix.hpp @@ -0,0 +1,47 @@ +#pragma once + +#include +#include + +#include "compute_t_col_inds.hpp" +#include "compute_t_num_non_zero.hpp" +#include "compute_t_row_ptrs.hpp" + +#include "src/dof_management/freedom_signature.hpp" + +namespace openturbine { + +/// Creates the tangent operator matrix T in sparse CRS format +template +[[nodiscard]] inline CrsMatrixType CreateTMatrix( + size_t system_dofs, + const Kokkos::View::const_type& node_freedom_allocation_table, + const Kokkos::View::const_type& node_freedom_map_table +) { + using ValuesType = typename CrsMatrixType::values_type::non_const_type; + using RowPtrType = typename CrsMatrixType::staticcrsgraph_type::row_map_type::non_const_type; + using IndicesType = typename CrsMatrixType::staticcrsgraph_type::entries_type::non_const_type; + + const auto T_num_rows = system_dofs; + const auto T_num_columns = T_num_rows; + const auto T_num_non_zero = ComputeTNumNonZero(node_freedom_allocation_table); + + const auto T_row_ptrs = ComputeTRowPtrs( + T_num_rows, node_freedom_allocation_table, node_freedom_map_table + ); + const auto T_col_inds = ComputeTColInds( + T_num_non_zero, node_freedom_allocation_table, node_freedom_map_table, T_row_ptrs + ); + const auto T_values = ValuesType("T values", T_num_non_zero); + + KokkosSparse::sort_crs_matrix(T_row_ptrs, T_col_inds, T_values); + return { + "T", + static_cast(T_num_rows), + static_cast(T_num_columns), + T_num_non_zero, + T_values, + T_row_ptrs, + T_col_inds}; +} +} // namespace openturbine diff --git a/src/solver/create_transpose_matrix_full.hpp b/src/solver/create_transpose_matrix_full.hpp new file mode 100644 index 00000000..2006cd34 --- /dev/null +++ b/src/solver/create_transpose_matrix_full.hpp @@ -0,0 +1,32 @@ +#pragma once + +#include +#include + +#include "fill_unshifted_row_ptrs.hpp" + +namespace openturbine { + +template +[[nodiscard]] inline CrsMatrixType CreateTransposeMatrixFull( + size_t num_system_dofs, size_t num_dofs, const CrsMatrixType& B_t +) { + using RowPtrType = typename CrsMatrixType::staticcrsgraph_type::row_map_type::non_const_type; + using IndicesType = typename CrsMatrixType::staticcrsgraph_type::entries_type::non_const_type; + + auto transpose_matrix_full_row_ptrs = RowPtrType("transpose_matrix_full_row_ptrs", num_dofs + 1); + Kokkos::parallel_for( + "FillUnshiftedRowPtrs", num_dofs + 1, + FillUnshiftedRowPtrs{ + num_system_dofs, B_t.graph.row_map, transpose_matrix_full_row_ptrs} + ); + auto transpose_matrix_full_indices = IndicesType("transpose_matrix_full_indices", B_t.nnz()); + Kokkos::deep_copy(transpose_matrix_full_indices, static_cast(num_system_dofs)); + KokkosBlas::axpy(1, B_t.graph.entries, transpose_matrix_full_indices); + return CrsMatrixType( + "transpose_matrix_full", static_cast(num_dofs), static_cast(num_dofs), B_t.nnz(), + B_t.values, transpose_matrix_full_row_ptrs, transpose_matrix_full_indices + ); +} + +} // namespace openturbine diff --git a/src/solver/solver.hpp b/src/solver/solver.hpp index 20d716a8..77d29ab4 100644 --- a/src/solver/solver.hpp +++ b/src/solver/solver.hpp @@ -8,10 +8,19 @@ #include #include -#include "compute_number_of_non_zeros_constraints.hpp" +#include "compute_num_system_dofs.hpp" +#include "create_b_matrix.hpp" +#include "create_b_t_matrix.hpp" +#include "create_constraints_matrix_full.hpp" +#include "create_global_matrix.hpp" +#include "create_k_matrix.hpp" +#include "create_matrix_spadd.hpp" +#include "create_matrix_spgemm.hpp" +#include "create_sparse_dense_solver.hpp" +#include "create_system_matrix_full.hpp" +#include "create_t_matrix.hpp" +#include "create_transpose_matrix_full.hpp" #include "fill_unshifted_row_ptrs.hpp" -#include "populate_sparse_row_ptrs_col_inds_constraints.hpp" -#include "populate_sparse_row_ptrs_col_inds_transpose.hpp" #include "src/constraints/constraint_type.hpp" @@ -34,26 +43,22 @@ struct Solver { using GlobalCrsMatrixType = Tpetra::CrsMatrix<>; using ExecutionSpace = GlobalCrsMatrixType::execution_space; using MemorySpace = GlobalCrsMatrixType::memory_space; - using GlobalMapType = GlobalCrsMatrixType::map_type; using GlobalMultiVectorType = Tpetra::MultiVector<>; - using DualViewType = GlobalMultiVectorType::dual_view_type; using CrsMatrixType = GlobalCrsMatrixType::local_matrix_device_type; - using ValuesType = CrsMatrixType::values_type::non_const_type; - using RowPtrType = CrsMatrixType::staticcrsgraph_type::row_map_type::non_const_type; - using IndicesType = CrsMatrixType::staticcrsgraph_type::entries_type::non_const_type; - using ScalarType = GlobalCrsMatrixType::scalar_type; - using LocalOrdinalType = GlobalCrsMatrixType::local_ordinal_type; - using GlobalOrdinalType = GlobalCrsMatrixType::global_ordinal_type; using KernelHandle = typename KokkosKernels::Experimental::KokkosKernelsHandle< - RowPtrType::value_type, IndicesType::value_type, ValuesType::value_type, ExecutionSpace, - MemorySpace, MemorySpace>; - using SpmvHandle = - KokkosSparse::SPMVHandle; + CrsMatrixType::const_size_type, CrsMatrixType::const_ordinal_type, + CrsMatrixType::const_value_type, ExecutionSpace, MemorySpace, MemorySpace>; size_t num_system_nodes; //< Number of system nodes size_t num_system_dofs; //< Number of system degrees of freedom size_t num_dofs; //< Number of degrees of freedom + KernelHandle system_spgemm_handle; + KernelHandle constraints_spgemm_handle; + KernelHandle system_spadd_handle; + KernelHandle spc_spadd_handle; + KernelHandle full_system_spadd_handle; + CrsMatrixType K; //< Stiffness matrix CrsMatrixType T; //< Tangent operator CrsMatrixType B; //< Constraints matrix @@ -74,422 +79,8 @@ struct Solver { std::vector convergence_err; - KernelHandle system_spgemm_handle; - KernelHandle constraints_spgemm_handle; - KernelHandle system_spadd_handle; - KernelHandle spc_spadd_handle; - KernelHandle full_system_spadd_handle; - Teuchos::RCP> amesos_solver; - /// Computes the total number of active degrees of freedom in the system - [[nodiscard]] static size_t ComputeNumSystemDofs( - const Kokkos::View::const_type& node_freedom_allocation_table - ) { - auto total_system_dofs = 0UL; - Kokkos::parallel_reduce( - "ComputeNumSystemDofs", node_freedom_allocation_table.extent(0), - KOKKOS_LAMBDA(size_t i, size_t & update) { - update += count_active_dofs(node_freedom_allocation_table(i)); - }, - total_system_dofs - ); - return total_system_dofs; - } - - /// Computes the number of non-zero entries in the stiffness matrix K for sparse storage - [[nodiscard]] static size_t ComputeKNumNonZero( - const Kokkos::View::const_type& num_nodes_per_element, - const Kokkos::View::const_type& node_state_indices, - const Kokkos::View::const_type& node_freedom_allocation_table - ) { - auto K_num_non_zero_off_diagonal = 0UL; - Kokkos::parallel_reduce( - "ComputeKNumNonZero_OffDiagonal", num_nodes_per_element.extent(0), - KOKKOS_LAMBDA(size_t i, size_t & update) { - auto num_element_dof = 0UL; - for (auto j = 0U; j < num_nodes_per_element(i); ++j) { - const auto num_node_dof = - count_active_dofs(node_freedom_allocation_table(node_state_indices(i, j))); - num_element_dof += num_node_dof; - } - const auto num_element_non_zero = num_element_dof * num_element_dof; - auto num_diagonal_non_zero = 0UL; - for (auto j = 0U; j < num_nodes_per_element(i); ++j) { - const auto num_node_dof = - count_active_dofs(node_freedom_allocation_table(node_state_indices(i, j))); - num_diagonal_non_zero += num_node_dof * num_node_dof; - } - update += num_element_non_zero - num_diagonal_non_zero; - }, - K_num_non_zero_off_diagonal - ); - auto K_num_non_zero_diagonal = 0UL; - Kokkos::parallel_reduce( - "ComputeKNumNonZero_Diagonal", node_freedom_allocation_table.extent(0), - KOKKOS_LAMBDA(size_t i, size_t & update) { - const auto num_node_dof = count_active_dofs(node_freedom_allocation_table(i)); - const auto num_diagonal_non_zero = num_node_dof * num_node_dof; - update += num_diagonal_non_zero; - }, - K_num_non_zero_diagonal - ); - return K_num_non_zero_off_diagonal + K_num_non_zero_diagonal; - } - - /// Computes the row pointers for the sparse stiffness matrix K in CSR format - [[nodiscard]] static RowPtrType ComputeKRowPtrs( - size_t K_num_rows, - const Kokkos::View::const_type& node_freedom_allocation_table, - const Kokkos::View::const_type& node_freedom_map_table, - const Kokkos::View::const_type& num_nodes_per_element, - const Kokkos::View::const_type& node_state_indices - ) { - auto K_row_ptrs = RowPtrType("row_ptrs", K_num_rows + 1); - const auto K_row_entries = RowPtrType("row_ptrs", K_num_rows); - - Kokkos::parallel_for( - "ComputeKRowEntries", node_freedom_allocation_table.extent(0), - KOKKOS_LAMBDA(size_t i) { - const auto this_node_num_dof = count_active_dofs(node_freedom_allocation_table(i)); - const auto this_node_dof_index = node_freedom_map_table(i); - - auto num_entries_in_row = this_node_num_dof; - bool node_found_in_system = false; - - // contributions to non-diagonal block from coupled nodes - for (auto e = 0U; e < num_nodes_per_element.extent(0); ++e) { - bool contains_node = false; - auto num_entries_in_element = 0UL; - for (auto j = 0U; j < num_nodes_per_element(e); ++j) { - contains_node = contains_node || (node_state_indices(e, j) == i); - num_entries_in_element += - count_active_dofs(node_freedom_allocation_table(node_state_indices(e, j)) - ); - } - if (contains_node) { - node_found_in_system = true; - num_entries_in_row += num_entries_in_element - this_node_num_dof; - } - } - if (node_found_in_system) { - for (auto j = 0U; j < this_node_num_dof; ++j) { - K_row_entries(this_node_dof_index + j) = num_entries_in_row; - } - } - } - ); - - auto result = 0UL; - Kokkos::parallel_scan( - "ComputeKRowPtrs", K_row_entries.extent(0), - KOKKOS_LAMBDA(size_t i, size_t & update, bool is_final) { - update += K_row_entries(i); - if (is_final) { - K_row_ptrs(i + 1) = update; - } - }, - result - ); - - return K_row_ptrs; - } - - [[nodiscard]] static IndicesType ComputeKColInds( - size_t K_num_non_zero, - const Kokkos::View::const_type& node_freedom_allocation_table, - const Kokkos::View::const_type& node_freedom_map_table, - const Kokkos::View::const_type& num_nodes_per_element, - const Kokkos::View::const_type& node_state_indices, const RowPtrType& K_row_ptrs - ) { - auto K_col_inds = IndicesType("col_inds", K_num_non_zero); - - Kokkos::parallel_for( - "ComputeKColInds", node_freedom_allocation_table.extent(0), - KOKKOS_LAMBDA(size_t i) { - const auto this_node_num_dof = count_active_dofs(node_freedom_allocation_table(i)); - const auto this_node_dof_index = node_freedom_map_table(i); - - for (auto j = 0U; j < this_node_num_dof; ++j) { - auto current_dof_index = K_row_ptrs(this_node_dof_index + j); - - for (auto k = 0U; k < this_node_num_dof; ++k, ++current_dof_index) { - K_col_inds(current_dof_index) = static_cast(this_node_dof_index + k); - } - - for (auto e = 0U; e < num_nodes_per_element.extent(0); ++e) { - bool contains_node = false; - for (auto n = 0U; n < num_nodes_per_element(e); ++n) { - contains_node = contains_node || (node_state_indices(e, n) == i); - } - if (!contains_node) { - continue; - } - for (auto n = 0U; n < num_nodes_per_element(e); ++n) { - if (node_state_indices(e, n) != i) { - const auto target_node_num_dof = count_active_dofs( - node_freedom_allocation_table(node_state_indices(e, n)) - ); - const auto target_node_dof_index = - node_freedom_map_table(node_state_indices(e, n)); - for (auto k = 0U; k < target_node_num_dof; - ++k, ++current_dof_index) { - K_col_inds(current_dof_index) = - static_cast(target_node_dof_index + k); - } - } - } - } - } - } - ); - - return K_col_inds; - } - - [[nodiscard]] static size_t ComputeTNumNonZero( - const Kokkos::View::const_type& node_freedom_allocation_table - ) { - auto T_num_non_zero = 0UL; - Kokkos::parallel_reduce( - "ComputeTNumNonZero", node_freedom_allocation_table.extent(0), - KOKKOS_LAMBDA(size_t i, size_t & update) { - const auto num_node_dof = count_active_dofs(node_freedom_allocation_table(i)); - const auto num_diagonal_non_zero = num_node_dof * num_node_dof; - update += num_diagonal_non_zero; - }, - T_num_non_zero - ); - - return T_num_non_zero; - } - - [[nodiscard]] static RowPtrType ComputeTRowPtrs( - size_t T_num_rows, - const Kokkos::View::const_type& node_freedom_allocation_table, - const Kokkos::View::const_type& node_freedom_map_table - ) { - auto T_row_ptrs = RowPtrType("T_row_ptrs", T_num_rows + 1); - const auto T_row_entries = RowPtrType("row_entries", T_num_rows); - - Kokkos::parallel_for( - "ComputeTRowEntries", node_freedom_allocation_table.extent(0), - KOKKOS_LAMBDA(size_t i) { - const auto this_node_num_dof = count_active_dofs(node_freedom_allocation_table(i)); - const auto this_node_dof_index = node_freedom_map_table(i); - - auto num_entries_in_row = this_node_num_dof; - - for (auto j = 0U; j < this_node_num_dof; ++j) { - T_row_entries(this_node_dof_index + j) = num_entries_in_row; - } - } - ); - - auto result = 0UL; - Kokkos::parallel_scan( - "ComputeTRowPtrs", T_row_entries.extent(0), - KOKKOS_LAMBDA(size_t i, size_t & update, bool is_final) { - update += T_row_entries(i); - if (is_final) { - T_row_ptrs(i + 1) = update; - } - }, - result - ); - - return T_row_ptrs; - } - - [[nodiscard]] static IndicesType ComputeTColInds( - size_t T_num_non_zero, - const Kokkos::View::const_type& node_freedom_allocation_table, - const Kokkos::View::const_type& node_freedom_map_table, const RowPtrType& T_row_ptrs - ) { - auto T_col_inds = IndicesType("T_indices", T_num_non_zero); - - Kokkos::parallel_for( - "ComputeTColInds", node_freedom_allocation_table.extent(0), - KOKKOS_LAMBDA(size_t i) { - const auto this_node_num_dof = count_active_dofs(node_freedom_allocation_table(i)); - const auto this_node_dof_index = node_freedom_map_table(i); - - for (auto j = 0U; j < this_node_num_dof; ++j) { - for (auto k = 0UL, current_dof_index = T_row_ptrs(this_node_dof_index + j); - k < this_node_num_dof; ++k, ++current_dof_index) { - T_col_inds(current_dof_index) = static_cast(this_node_dof_index + k); - } - } - } - ); - - return T_col_inds; - } - - [[nodiscard]] static size_t ComputeBNumNonZero( - const Kokkos::View::const_type& constraint_type, - const Kokkos::View*>::const_type& constraint_row_range - ) { - auto B_num_non_zero = size_t{0U}; - Kokkos::parallel_reduce( - "ComputeNumberOfNonZeros_Constraints", constraint_type.extent(0), - ComputeNumberOfNonZeros_Constraints{constraint_type, constraint_row_range}, - B_num_non_zero - ); - return B_num_non_zero; - } - - /// Creates the system stiffness matrix K in sparse CRS format - [[nodiscard]] static CrsMatrixType CreateKMatrix( - size_t system_dofs, - const Kokkos::View::const_type& node_freedom_allocation_table, - const Kokkos::View::const_type& node_freedom_map_table, - const Kokkos::View::const_type& num_nodes_per_element, - const Kokkos::View::const_type& node_state_indices - ) { - const auto K_num_rows = system_dofs; - const auto K_num_columns = K_num_rows; - const auto K_num_non_zero = ComputeKNumNonZero( - num_nodes_per_element, node_state_indices, node_freedom_allocation_table - ); - - const auto K_row_ptrs = ComputeKRowPtrs( - K_num_rows, node_freedom_allocation_table, node_freedom_map_table, num_nodes_per_element, - node_state_indices - ); - const auto K_col_inds = ComputeKColInds( - K_num_non_zero, node_freedom_allocation_table, node_freedom_map_table, - num_nodes_per_element, node_state_indices, K_row_ptrs - ); - const auto K_values = ValuesType("K values", K_num_non_zero); - - KokkosSparse::sort_crs_matrix(K_row_ptrs, K_col_inds, K_values); - return { - "K", - static_cast(K_num_rows), - static_cast(K_num_columns), - K_num_non_zero, - K_values, - K_row_ptrs, - K_col_inds - }; - } - - /// Creates the tangent operator matrix T in sparse CRS format - [[nodiscard]] static CrsMatrixType CreateTMatrix( - size_t system_dofs, - const Kokkos::View::const_type& node_freedom_allocation_table, - const Kokkos::View::const_type& node_freedom_map_table - ) { - const auto T_num_rows = system_dofs; - const auto T_num_columns = T_num_rows; - const auto T_num_non_zero = ComputeTNumNonZero(node_freedom_allocation_table); - - const auto T_row_ptrs = - ComputeTRowPtrs(T_num_rows, node_freedom_allocation_table, node_freedom_map_table); - const auto T_col_inds = ComputeTColInds( - T_num_non_zero, node_freedom_allocation_table, node_freedom_map_table, T_row_ptrs - ); - const auto T_values = ValuesType("T values", T_num_non_zero); - - KokkosSparse::sort_crs_matrix(T_row_ptrs, T_col_inds, T_values); - return { - "T", - static_cast(T_num_rows), - static_cast(T_num_columns), - T_num_non_zero, - T_values, - T_row_ptrs, - T_col_inds - }; - } - - /// Creates the constraint matrix B in sparse CRS format - [[nodiscard]] static CrsMatrixType CreateBMatrix( - size_t system_dofs, size_t constraint_dofs, - const Kokkos::View::const_type& constraint_type, - const Kokkos::View::const_type& constraint_base_node_freedom_table, - const Kokkos::View::const_type& constraint_target_node_freedom_table, - const Kokkos::View*>::const_type& constraint_row_range - ) { - const auto B_num_rows = constraint_dofs; - const auto B_num_columns = system_dofs; - const auto B_num_non_zero = ComputeBNumNonZero(constraint_type, constraint_row_range); - - const auto B_row_ptrs = RowPtrType("b_row_ptrs", B_num_rows + 1); - const auto B_col_ind = IndicesType("b_indices", B_num_non_zero); - Kokkos::parallel_for( - "PopulateSparseRowPtrsColInds_Constraints", 1, - PopulateSparseRowPtrsColInds_Constraints{ - constraint_type, constraint_base_node_freedom_table, - constraint_target_node_freedom_table, constraint_row_range, B_row_ptrs, B_col_ind - } - ); - const auto B_values = ValuesType("B values", B_num_non_zero); - KokkosSparse::sort_crs_matrix(B_row_ptrs, B_col_ind, B_values); - return { - "B", - static_cast(B_num_rows), - static_cast(B_num_columns), - B_num_non_zero, - B_values, - B_row_ptrs, - B_col_ind - }; - } - - /// Creates the constraint matrix B_t in sparse CRS format - [[nodiscard]] static CrsMatrixType CreateBtMatrix( - size_t system_dofs, size_t constraint_dofs, - const Kokkos::View::const_type& constraint_type, - const Kokkos::View::const_type& constraint_base_node_freedom_table, - const Kokkos::View::const_type& constraint_target_node_freedom_table, - const Kokkos::View*>::const_type& constraint_row_range - ) { - const auto B_num_rows = constraint_dofs; - const auto B_num_columns = system_dofs; - const auto B_num_non_zero = ComputeBNumNonZero(constraint_type, constraint_row_range); - - const auto B_row_ptrs = RowPtrType("b_row_ptrs", B_num_rows + 1); - const auto B_col_ind = IndicesType("b_indices", B_num_non_zero); - Kokkos::parallel_for( - "PopulateSparseRowPtrsColInds_Constraints", 1, - PopulateSparseRowPtrsColInds_Constraints{ - constraint_type, constraint_base_node_freedom_table, - constraint_target_node_freedom_table, constraint_row_range, B_row_ptrs, B_col_ind - } - ); - const auto B_values = ValuesType("B values", B_num_non_zero); - KokkosSparse::sort_crs_matrix(B_row_ptrs, B_col_ind, B_values); - - const auto B_t_num_rows = system_dofs; - const auto B_t_num_columns = constraint_dofs; - const auto B_t_num_non_zero = ComputeBNumNonZero(constraint_type, constraint_row_range); - - auto col_count = IndicesType("col_count", B_num_columns); - auto tmp_row_ptrs = RowPtrType("tmp_row_ptrs", B_t_num_rows + 1); - auto B_t_row_ptrs = RowPtrType("b_t_row_ptrs", B_t_num_rows + 1); - auto B_t_col_inds = IndicesType("B_t_indices", B_t_num_non_zero); - auto B_t_values = ValuesType("B_t values", B_t_num_non_zero); - Kokkos::parallel_for( - "PopulateSparseRowPtrsColInds_Transpose", 1, - PopulateSparseRowPtrsColInds_Transpose{ - B_num_rows, B_num_columns, B_row_ptrs, B_col_ind, col_count, tmp_row_ptrs, - B_t_row_ptrs, B_t_col_inds - } - ); - KokkosSparse::sort_crs_matrix(B_t_row_ptrs, B_t_col_inds, B_t_values); - return { - "B_t", - static_cast(B_t_num_rows), - static_cast(B_t_num_columns), - B_t_num_non_zero, - B_t_values, - B_t_row_ptrs, - B_t_col_inds - }; - } - /** @brief Constructs a sparse linear systems solver for OpenTurbine * * @param node_IDs View containing the global IDs for each node @@ -517,126 +108,45 @@ struct Solver { : num_system_nodes(node_IDs.extent(0)), num_system_dofs(ComputeNumSystemDofs(node_freedom_allocation_table)), num_dofs(num_system_dofs + num_constraint_dofs), - K(CreateKMatrix( + K(CreateKMatrix( num_system_dofs, node_freedom_allocation_table, node_freedom_map_table, num_nodes_per_element, node_state_indices )), - T(CreateTMatrix(num_system_dofs, node_freedom_allocation_table, node_freedom_map_table)), - B(CreateBMatrix( + T(CreateTMatrix( + num_system_dofs, node_freedom_allocation_table, node_freedom_map_table + )), + B(CreateBMatrix( num_system_dofs, num_constraint_dofs, constraint_type, constraint_base_node_freedom_table, constraint_target_node_freedom_table, constraint_row_range )), - B_t(CreateBtMatrix( + B_t(CreateBtMatrix( num_system_dofs, num_constraint_dofs, constraint_type, constraint_base_node_freedom_table, constraint_target_node_freedom_table, constraint_row_range )), + static_system_matrix(CreateMatrixSpgemm(K, T, system_spgemm_handle)), + system_matrix(CreateMatrixSpadd(K, static_system_matrix, system_spadd_handle)), + constraints_matrix(CreateMatrixSpgemm(B, T, constraints_spgemm_handle)), + system_matrix_full(CreateSystemMatrixFull(num_system_dofs, num_dofs, system_matrix)), + constraints_matrix_full( + CreateConstraintsMatrixFull(num_system_dofs, num_dofs, constraints_matrix) + ), + transpose_matrix_full(CreateTransposeMatrixFull(num_system_dofs, num_dofs, B_t)), + system_plus_constraints( + CreateMatrixSpadd(system_matrix_full, constraints_matrix_full, spc_spadd_handle) + ), + full_matrix(CreateMatrixSpadd( + system_plus_constraints, transpose_matrix_full, full_system_spadd_handle + )), + A(CreateGlobalMatrix(full_matrix)), + b(Tpetra::createMultiVector(A->getRangeMap(), 1)), + x_mv(Tpetra::createMultiVector(A->getDomainMap(), 1)), R("R", num_dofs), - x("x", num_dofs) { - system_spgemm_handle.create_spgemm_handle(); - KokkosSparse::spgemm_symbolic( - system_spgemm_handle, K, false, T, false, static_system_matrix - ); - KokkosSparse::spgemm_numeric(system_spgemm_handle, K, false, T, false, static_system_matrix); - - constraints_spgemm_handle.create_spgemm_handle(); - KokkosSparse::spgemm_symbolic( - constraints_spgemm_handle, B, false, T, false, constraints_matrix - ); - KokkosSparse::spgemm_numeric( - constraints_spgemm_handle, B, false, T, false, constraints_matrix - ); - - system_spadd_handle.create_spadd_handle(true, true); - KokkosSparse::spadd_symbolic(&system_spadd_handle, K, static_system_matrix, system_matrix); - KokkosSparse::spadd_numeric( - &system_spadd_handle, 1., K, 1., static_system_matrix, system_matrix - ); - - auto system_matrix_full_row_ptrs = RowPtrType("system_matrix_full_row_ptrs", num_dofs + 1); - Kokkos::parallel_for( - "FillUnshiftedRowPtrs", num_dofs + 1, - FillUnshiftedRowPtrs{ - num_system_dofs, system_matrix.graph.row_map, system_matrix_full_row_ptrs - } - ); - system_matrix_full = CrsMatrixType( - "system_matrix_full", static_cast(num_dofs), static_cast(num_dofs), - system_matrix.nnz(), system_matrix.values, system_matrix_full_row_ptrs, - system_matrix.graph.entries - ); - - auto constraints_matrix_full_row_ptrs = - RowPtrType("constraints_matrix_full_row_ptrs", num_dofs + 1); - Kokkos::deep_copy( - Kokkos::subview( - constraints_matrix_full_row_ptrs, Kokkos::pair(num_system_dofs, num_dofs + 1) - ), - constraints_matrix.graph.row_map - ); - constraints_matrix_full = CrsMatrixType( - "constraints_matrix_full", static_cast(num_dofs), static_cast(num_dofs), - constraints_matrix.nnz(), constraints_matrix.values, constraints_matrix_full_row_ptrs, - constraints_matrix.graph.entries - ); - - auto transpose_matrix_full_row_ptrs = - RowPtrType("transpose_matrix_full_row_ptrs", num_dofs + 1); - Kokkos::parallel_for( - "FillUnshiftedRowPtrs", num_dofs + 1, - FillUnshiftedRowPtrs{ - num_system_dofs, B_t.graph.row_map, transpose_matrix_full_row_ptrs - } - ); - auto transpose_matrix_full_indices = IndicesType("transpose_matrix_full_indices", B_t.nnz()); - Kokkos::deep_copy(transpose_matrix_full_indices, static_cast(num_system_dofs)); - KokkosBlas::axpy(1, B_t.graph.entries, transpose_matrix_full_indices); - transpose_matrix_full = CrsMatrixType( - "transpose_matrix_full", static_cast(num_dofs), static_cast(num_dofs), - B_t.nnz(), B_t.values, transpose_matrix_full_row_ptrs, transpose_matrix_full_indices - ); - - spc_spadd_handle.create_spadd_handle(true, true); - KokkosSparse::spadd_symbolic( - &spc_spadd_handle, system_matrix_full, constraints_matrix_full, system_plus_constraints - ); - KokkosSparse::spadd_numeric( - &spc_spadd_handle, 1., system_matrix_full, 1., constraints_matrix_full, - system_plus_constraints - ); - - full_system_spadd_handle.create_spadd_handle(true, true); - KokkosSparse::spadd_symbolic( - &full_system_spadd_handle, system_plus_constraints, transpose_matrix_full, full_matrix - ); - KokkosSparse::spadd_numeric( - &full_system_spadd_handle, 1., system_plus_constraints, 1., transpose_matrix_full, - full_matrix - ); - - auto comm = Teuchos::createSerialComm(); - auto row_map = Tpetra::createLocalMap( - static_cast(full_matrix.numRows()), comm - ); - auto col_map = Tpetra::createLocalMap( - static_cast(full_matrix.numCols()), comm - ); - - A = Teuchos::make_rcp( - row_map, col_map, CrsMatrixType("A", full_matrix) - ); - b = Tpetra::createMultiVector(A->getRangeMap(), 1); - x_mv = Tpetra::createMultiVector(A->getDomainMap(), 1); - - const auto solver_name = (std::is_same_v) - ? std::string{"klu2"} - : std::string{"basker"}; - - amesos_solver = - Amesos2::create(solver_name, A, x_mv, b); - amesos_solver->symbolicFactorization(); - } + x("x", num_dofs), + amesos_solver( + CreateSparseDenseSolver(A, x_mv, b) + ) {} }; } // namespace openturbine diff --git a/src/step/assemble_constraints_matrix.hpp b/src/step/assemble_constraints_matrix.hpp index d0f71619..a59c33a7 100644 --- a/src/step/assemble_constraints_matrix.hpp +++ b/src/step/assemble_constraints_matrix.hpp @@ -35,7 +35,7 @@ inline void AssembleConstraintsMatrix(Solver& solver, Constraints& constraints) auto trans_region = Kokkos::Profiling::ScopedRegion("Transpose Constraints Matrix"); auto B_num_rows = solver.B.numRows(); auto constraint_transpose_policy = Kokkos::TeamPolicy<>(B_num_rows, Kokkos::AUTO()); - auto tmp_row_map = Solver::RowPtrType( + auto tmp_row_map = Solver::CrsMatrixType::row_map_type::non_const_type( Kokkos::view_alloc(Kokkos::WithoutInitializing, "tmp_row_map"), solver.B_t.graph.row_map.extent(0) ); diff --git a/src/step/assemble_constraints_residual.hpp b/src/step/assemble_constraints_residual.hpp index 9870575e..e3449bbe 100644 --- a/src/step/assemble_constraints_residual.hpp +++ b/src/step/assemble_constraints_residual.hpp @@ -24,17 +24,22 @@ inline void AssembleConstraintsResidual(Solver& solver, Constraints& constraints } ); - auto R = Solver::ValuesType( + using CrsMatrixType = Solver::CrsMatrixType; + using VectorType = CrsMatrixType::values_type::non_const_type; + + auto R = VectorType( Kokkos::view_alloc(Kokkos::WithoutInitializing, "R_local"), solver.num_system_dofs ); Kokkos::deep_copy( R, Kokkos::subview(solver.R, Kokkos::make_pair(size_t{0U}, solver.num_system_dofs)) ); - auto lambda = Solver::ValuesType( + auto lambda = VectorType( Kokkos::view_alloc(Kokkos::WithoutInitializing, "lambda"), constraints.lambda.extent(0) ); Kokkos::deep_copy(lambda, constraints.lambda); - auto spmv_handle = Solver::SpmvHandle(); + + auto spmv_handle = + KokkosSparse::SPMVHandle(); KokkosSparse::spmv(&spmv_handle, "T", 1., solver.B, lambda, 1., R); Kokkos::deep_copy( Kokkos::subview(solver.R, Kokkos::make_pair(size_t{0U}, solver.num_system_dofs)), R diff --git a/tests/unit_tests/solver/CMakeLists.txt b/tests/unit_tests/solver/CMakeLists.txt index dc26546b..78dba690 100644 --- a/tests/unit_tests/solver/CMakeLists.txt +++ b/tests/unit_tests/solver/CMakeLists.txt @@ -2,6 +2,7 @@ target_sources( openturbine_unit_tests PRIVATE +# test_compute_k_num_non_zero.cpp test_compute_number_of_non_zeros_constraints.cpp test_copy_into_sparse_matrix.cpp test_fill_unshifted_row_ptrs.cpp From 82d4ba854de197824d0a77aedefb3586e76e1dd4 Mon Sep 17 00:00:00 2001 From: dcdemen Date: Mon, 23 Dec 2024 16:33:08 -0700 Subject: [PATCH 2/3] Added unit tests for creating the K and T sparse matrix structures --- src/solver/compute_k_col_inds.hpp | 6 +- src/solver/create_k_matrix.hpp | 4 +- tests/unit_tests/solver/CMakeLists.txt | 7 +- .../solver/test_compute_k_col_inds.cpp | 234 ++++++++++++++++++ .../solver/test_compute_k_num_non_zero.cpp | 94 +++++-- .../solver/test_compute_k_row_ptrs.cpp | 195 +++++++++++++++ .../solver/test_compute_t_col_inds.cpp | 118 +++++++++ .../solver/test_compute_t_num_non_zero.cpp | 27 ++ .../solver/test_compute_t_row_ptrs.cpp | 83 +++++++ .../solver/test_create_k_col_inds.cpp | 83 ------- .../solver/test_create_k_row_ptrs.cpp | 53 ---- .../solver/test_create_t_col_inds.cpp | 42 ---- .../solver/test_create_t_row_ptrs.cpp | 28 --- 13 files changed, 737 insertions(+), 237 deletions(-) create mode 100644 tests/unit_tests/solver/test_compute_k_col_inds.cpp create mode 100644 tests/unit_tests/solver/test_compute_k_row_ptrs.cpp create mode 100644 tests/unit_tests/solver/test_compute_t_col_inds.cpp create mode 100644 tests/unit_tests/solver/test_compute_t_num_non_zero.cpp create mode 100644 tests/unit_tests/solver/test_compute_t_row_ptrs.cpp delete mode 100644 tests/unit_tests/solver/test_create_k_col_inds.cpp delete mode 100644 tests/unit_tests/solver/test_create_k_row_ptrs.cpp delete mode 100644 tests/unit_tests/solver/test_create_t_col_inds.cpp delete mode 100644 tests/unit_tests/solver/test_create_t_row_ptrs.cpp diff --git a/src/solver/compute_k_col_inds.hpp b/src/solver/compute_k_col_inds.hpp index 0248d4da..e27f42be 100644 --- a/src/solver/compute_k_col_inds.hpp +++ b/src/solver/compute_k_col_inds.hpp @@ -27,7 +27,8 @@ struct ComputeKColIndsFunction { } KOKKOS_FUNCTION - void ComputeColInds(size_t element, size_t node, RowPtrValueType current_dof_index) const { + RowPtrValueType ComputeColInds(size_t element, size_t node, RowPtrValueType current_dof_index) + const { for (auto n = 0U; n < num_nodes_per_element(element); ++n) { const auto node_state_index = node_state_indices(element, n); if (node_state_index != node) { @@ -40,6 +41,7 @@ struct ComputeKColIndsFunction { } } } + return current_dof_index; } KOKKOS_FUNCTION @@ -58,7 +60,7 @@ struct ComputeKColIndsFunction { if (!ContainsNode(e, i)) { continue; } - ComputeColInds(e, i, current_dof_index); + current_dof_index = ComputeColInds(e, i, current_dof_index); } } } diff --git a/src/solver/create_k_matrix.hpp b/src/solver/create_k_matrix.hpp index ad97bdb0..e9c284de 100644 --- a/src/solver/create_k_matrix.hpp +++ b/src/solver/create_k_matrix.hpp @@ -21,8 +21,8 @@ template const Kokkos::View::const_type& node_state_indices ) { using ValuesType = typename CrsMatrixType::values_type::non_const_type; - using RowPtrType = typename CrsMatrixType::staticcrsgraph_type::row_map_type::non_const_type; - using IndicesType = typename CrsMatrixType::staticcrsgraph_type::entries_type::non_const_type; + using RowPtrType = typename CrsMatrixType::row_map_type::non_const_type; + using IndicesType = typename CrsMatrixType::index_type::non_const_type; const auto K_num_rows = system_dofs; const auto K_num_columns = K_num_rows; diff --git a/tests/unit_tests/solver/CMakeLists.txt b/tests/unit_tests/solver/CMakeLists.txt index 78dba690..c0d166e9 100644 --- a/tests/unit_tests/solver/CMakeLists.txt +++ b/tests/unit_tests/solver/CMakeLists.txt @@ -2,8 +2,13 @@ target_sources( openturbine_unit_tests PRIVATE -# test_compute_k_num_non_zero.cpp + test_compute_k_col_inds.cpp + test_compute_k_row_ptrs.cpp + test_compute_k_num_non_zero.cpp test_compute_number_of_non_zeros_constraints.cpp + test_compute_t_num_non_zero.cpp + test_compute_t_col_inds.cpp + test_compute_t_row_ptrs.cpp test_copy_into_sparse_matrix.cpp test_fill_unshifted_row_ptrs.cpp test_populate_sparse_row_ptrs_col_inds_constraints.cpp diff --git a/tests/unit_tests/solver/test_compute_k_col_inds.cpp b/tests/unit_tests/solver/test_compute_k_col_inds.cpp new file mode 100644 index 00000000..2c655a65 --- /dev/null +++ b/tests/unit_tests/solver/test_compute_k_col_inds.cpp @@ -0,0 +1,234 @@ +#include +#include +#include + +#include "src/solver/compute_k_col_inds.hpp" + +namespace openturbine::tests { + +TEST(ComputeKColInds, OneElementOneNode) { + constexpr auto K_num_non_zero = 36UL; + + const auto node_freedom_allocation_table = + Kokkos::View("node_freedom_allocation_table"); + Kokkos::deep_copy(node_freedom_allocation_table, FreedomSignature::AllComponents); + + const auto node_freedom_map_table = Kokkos::View("node_freedom_map_table"); + Kokkos::deep_copy(node_freedom_map_table, 0UL); + + const auto num_nodes_per_element = Kokkos::View("num_nodes_per_element"); + Kokkos::deep_copy(num_nodes_per_element, 1UL); + + const auto node_state_indices = Kokkos::View("node_state_indices"); + Kokkos::deep_copy(node_state_indices, 0UL); + + const auto K_row_ptrs = Kokkos::View("K_row_ptrs"); + constexpr auto K_row_ptrs_host_data = std::array{0, 6, 12, 18, 24, 30, 36}; + const auto K_row_ptrs_host = + Kokkos::View::const_type(K_row_ptrs_host_data.data()); + const auto K_row_ptrs_mirror = Kokkos::create_mirror(K_row_ptrs); + Kokkos::deep_copy(K_row_ptrs_mirror, K_row_ptrs_host); + Kokkos::deep_copy(K_row_ptrs, K_row_ptrs_mirror); + + const auto col_inds = ComputeKColInds, Kokkos::View>( + K_num_non_zero, node_freedom_allocation_table, node_freedom_map_table, num_nodes_per_element, + node_state_indices, K_row_ptrs + ); + + const auto col_inds_mirror = Kokkos::create_mirror(col_inds); + Kokkos::deep_copy(col_inds_mirror, col_inds); + + for (auto row = 0U; row < 6U; ++row) { + for (auto col = 0U; col < 6U; ++col) { + EXPECT_EQ(col_inds_mirror(row * 6U + col), col); + } + } +} + +TEST(ComputeKColInds, OneElementTwoNodes) { + constexpr auto K_num_non_zero = 144UL; + + const auto node_freedom_allocation_table = + Kokkos::View("node_freedom_allocation_table"); + Kokkos::deep_copy(node_freedom_allocation_table, FreedomSignature::AllComponents); + + const auto node_freedom_map_table = Kokkos::View("node_freedom_map_table"); + constexpr auto node_freedom_map_table_host_data = std::array{0UL, 6UL}; + const auto node_freedom_map_table_host = Kokkos::View::const_type( + node_freedom_map_table_host_data.data() + ); + const auto node_freedom_map_table_mirror = Kokkos::create_mirror(node_freedom_map_table); + Kokkos::deep_copy(node_freedom_map_table_mirror, node_freedom_map_table_host); + Kokkos::deep_copy(node_freedom_map_table, node_freedom_map_table_mirror); + + const auto num_nodes_per_element = Kokkos::View("num_nodes_per_element"); + Kokkos::deep_copy(num_nodes_per_element, 2UL); + + const auto node_state_indices = Kokkos::View("node_state_indices"); + constexpr auto node_state_indices_host_data = std::array{0UL, 1UL}; + const auto node_state_indices_host = + Kokkos::View::const_type(node_state_indices_host_data.data() + ); + const auto node_state_indices_mirror = Kokkos::create_mirror(node_state_indices); + Kokkos::deep_copy(node_state_indices_mirror, node_state_indices_host); + Kokkos::deep_copy(node_state_indices, node_state_indices_mirror); + + const auto K_row_ptrs = Kokkos::View("K_row_ptrs"); + constexpr auto K_row_ptrs_host_data = + std::array{0, 12, 24, 36, 48, 60, 72, 84, 96, 108, 120, 132, 144}; + const auto K_row_ptrs_host = + Kokkos::View::const_type(K_row_ptrs_host_data.data()); + const auto K_row_ptrs_mirror = Kokkos::create_mirror(K_row_ptrs); + Kokkos::deep_copy(K_row_ptrs_mirror, K_row_ptrs_host); + Kokkos::deep_copy(K_row_ptrs, K_row_ptrs_mirror); + + const auto K_values = Kokkos::View("K_values", K_num_non_zero); + + const auto col_inds = ComputeKColInds, Kokkos::View>( + K_num_non_zero, node_freedom_allocation_table, node_freedom_map_table, num_nodes_per_element, + node_state_indices, K_row_ptrs + ); + + KokkosSparse::sort_crs_matrix(K_row_ptrs, col_inds, K_values); + + const auto col_inds_mirror = Kokkos::create_mirror(col_inds); + Kokkos::deep_copy(col_inds_mirror, col_inds); + + for (auto row = 0U; row < 12U; ++row) { + for (auto col = 0U; col < 12U; ++col) { + EXPECT_EQ(col_inds_mirror(row * 12U + col), col); + } + } +} + +TEST(ComputeKColInds, TwoElementsTwoNodesNoOverlap) { + constexpr auto K_num_non_zero = 288UL; + + const auto node_freedom_allocation_table = + Kokkos::View("node_freedom_allocation_table"); + Kokkos::deep_copy(node_freedom_allocation_table, FreedomSignature::AllComponents); + + const auto node_freedom_map_table = Kokkos::View("node_freedom_map_table"); + constexpr auto node_freedom_map_table_host_data = std::array{0UL, 6UL, 12UL, 18UL}; + const auto node_freedom_map_table_host = Kokkos::View::const_type( + node_freedom_map_table_host_data.data() + ); + const auto node_freedom_map_table_mirror = Kokkos::create_mirror(node_freedom_map_table); + Kokkos::deep_copy(node_freedom_map_table_mirror, node_freedom_map_table_host); + Kokkos::deep_copy(node_freedom_map_table, node_freedom_map_table_mirror); + + const auto num_nodes_per_element = Kokkos::View("num_nodes_per_element"); + Kokkos::deep_copy(num_nodes_per_element, 2UL); + + const auto node_state_indices = Kokkos::View("node_state_indices"); + constexpr auto node_state_indices_host_data = std::array{0UL, 1UL, 2UL, 3UL}; + const auto node_state_indices_host = + Kokkos::View::const_type(node_state_indices_host_data.data() + ); + const auto node_state_indices_mirror = Kokkos::create_mirror(node_state_indices); + Kokkos::deep_copy(node_state_indices_mirror, node_state_indices_host); + Kokkos::deep_copy(node_state_indices, node_state_indices_mirror); + + const auto K_row_ptrs = Kokkos::View("K_row_ptrs"); + constexpr auto K_row_ptrs_host_data = + std::array{0, 12, 24, 36, 48, 60, 72, 84, 96, 108, 120, 132, 144, + 156, 168, 180, 192, 204, 216, 228, 240, 252, 264, 276, 288}; + const auto K_row_ptrs_host = + Kokkos::View::const_type(K_row_ptrs_host_data.data()); + const auto K_row_ptrs_mirror = Kokkos::create_mirror(K_row_ptrs); + Kokkos::deep_copy(K_row_ptrs_mirror, K_row_ptrs_host); + Kokkos::deep_copy(K_row_ptrs, K_row_ptrs_mirror); + + const auto K_values = Kokkos::View("K_values", K_num_non_zero); + + const auto col_inds = ComputeKColInds, Kokkos::View>( + K_num_non_zero, node_freedom_allocation_table, node_freedom_map_table, num_nodes_per_element, + node_state_indices, K_row_ptrs + ); + + KokkosSparse::sort_crs_matrix(K_row_ptrs, col_inds, K_values); + + const auto col_inds_mirror = Kokkos::create_mirror(col_inds); + Kokkos::deep_copy(col_inds_mirror, col_inds); + + for (auto row = 0U; row < 12U; ++row) { + for (auto col = 0U; col < 12U; ++col) { + EXPECT_EQ(col_inds_mirror(row * 12U + col), col); + } + } + + for (auto row = 0U; row < 12U; ++row) { + for (auto col = 0U; col < 12U; ++col) { + EXPECT_EQ(col_inds_mirror(144U + row * 12U + col), 12U + col); + } + } +} + +TEST(ComputeKColInds, TwoElementsTwoNodesOverlap) { + constexpr auto K_num_non_zero = 252UL; + + const auto node_freedom_allocation_table = + Kokkos::View("node_freedom_allocation_table"); + Kokkos::deep_copy(node_freedom_allocation_table, FreedomSignature::AllComponents); + + const auto node_freedom_map_table = Kokkos::View("node_freedom_map_table"); + constexpr auto node_freedom_map_table_host_data = std::array{0UL, 6UL, 12UL}; + const auto node_freedom_map_table_host = Kokkos::View::const_type( + node_freedom_map_table_host_data.data() + ); + const auto node_freedom_map_table_mirror = Kokkos::create_mirror(node_freedom_map_table); + Kokkos::deep_copy(node_freedom_map_table_mirror, node_freedom_map_table_host); + Kokkos::deep_copy(node_freedom_map_table, node_freedom_map_table_mirror); + + const auto num_nodes_per_element = Kokkos::View("num_nodes_per_element"); + Kokkos::deep_copy(num_nodes_per_element, 2UL); + + const auto node_state_indices = Kokkos::View("node_state_indices"); + constexpr auto node_state_indices_host_data = std::array{0UL, 1UL, 1UL, 2UL}; + const auto node_state_indices_host = + Kokkos::View::const_type(node_state_indices_host_data.data() + ); + const auto node_state_indices_mirror = Kokkos::create_mirror(node_state_indices); + Kokkos::deep_copy(node_state_indices_mirror, node_state_indices_host); + Kokkos::deep_copy(node_state_indices, node_state_indices_mirror); + + const auto K_row_ptrs = Kokkos::View("K_row_ptrs"); + constexpr auto K_row_ptrs_host_data = std::array{ + 0, 12, 24, 36, 48, 60, 72, 90, 108, 126, 144, 162, 180, 192, 204, 216, 228, 240, 252}; + const auto K_row_ptrs_host = + Kokkos::View::const_type(K_row_ptrs_host_data.data()); + const auto K_row_ptrs_mirror = Kokkos::create_mirror(K_row_ptrs); + Kokkos::deep_copy(K_row_ptrs_mirror, K_row_ptrs_host); + Kokkos::deep_copy(K_row_ptrs, K_row_ptrs_mirror); + + const auto K_values = Kokkos::View("K_values", K_num_non_zero); + + const auto col_inds = ComputeKColInds, Kokkos::View>( + K_num_non_zero, node_freedom_allocation_table, node_freedom_map_table, num_nodes_per_element, + node_state_indices, K_row_ptrs + ); + + KokkosSparse::sort_crs_matrix(K_row_ptrs, col_inds, K_values); + + const auto col_inds_mirror = Kokkos::create_mirror(col_inds); + Kokkos::deep_copy(col_inds_mirror, col_inds); + + for (auto row = 0U; row < 6U; ++row) { + for (auto col = 0U; col < 12U; ++col) { + EXPECT_EQ(col_inds_mirror(row * 12U + col), col); + } + } + + for (auto row = 0U; row < 6U; ++row) { + for (auto col = 0U; col < 18U; ++col) { + EXPECT_EQ(col_inds_mirror(72U + row * 18U + col), col); + } + } + + for (auto row = 0U; row < 6U; ++row) { + for (auto col = 0U; col < 12U; ++col) { + EXPECT_EQ(col_inds_mirror(180U + row * 12U + col), 6U + col); + } + } +} +} // namespace openturbine::tests diff --git a/tests/unit_tests/solver/test_compute_k_num_non_zero.cpp b/tests/unit_tests/solver/test_compute_k_num_non_zero.cpp index c9dd8298..68ac307c 100644 --- a/tests/unit_tests/solver/test_compute_k_num_non_zero.cpp +++ b/tests/unit_tests/solver/test_compute_k_num_non_zero.cpp @@ -1,36 +1,78 @@ #include #include -#include "src/solver/compute_number_of_non_zeros.hpp" +#include "src/solver/compute_k_num_non_zero.hpp" namespace openturbine::tests { -TEST(ComputeNumberOfNonZeros, SingleElement) { - auto num_nodes = Kokkos::View("num_nodes"); - auto num_nodes_host = Kokkos::create_mirror(num_nodes); - num_nodes_host(0) = size_t{5U}; - Kokkos::deep_copy(num_nodes, num_nodes_host); - - auto num_non_zero = size_t{0U}; - Kokkos::parallel_reduce(1, ComputeNumberOfNonZeros{num_nodes}, num_non_zero); - constexpr auto num_dof = 5 * 6; - constexpr auto expected_num_non_zero = num_dof * num_dof; - EXPECT_EQ(num_non_zero, expected_num_non_zero); +TEST(ComputeKNumNonZero, OneElementFiveNodes) { + auto num_nodes_per_element = Kokkos::View("num_nodes_per_element"); + Kokkos::deep_copy(num_nodes_per_element, 5UL); + + auto node_state_indices = Kokkos::View("node_state_indices"); + constexpr auto node_state_indices_host_data = std::array{0UL, 1UL, 2UL, 3UL, 4UL}; + const auto node_state_indices_host = + Kokkos::View::const_type(node_state_indices_host_data.data() + ); + const auto node_state_indices_mirror = Kokkos::create_mirror(node_state_indices); + Kokkos::deep_copy(node_state_indices_mirror, node_state_indices_host); + Kokkos::deep_copy(node_state_indices, node_state_indices_mirror); + + auto node_freedom_allocation_table = + Kokkos::View("node_freedom_allocation_table"); + Kokkos::deep_copy(node_freedom_allocation_table, FreedomSignature::AllComponents); + + const auto num_non_zero = + ComputeKNumNonZero(num_nodes_per_element, node_state_indices, node_freedom_allocation_table); + + EXPECT_EQ(num_non_zero, 900UL); +} + +TEST(ComputeKNumNonZero, TwoElementsFiveNodesNoOverlap) { + auto num_nodes_per_element = Kokkos::View("num_nodes_per_element"); + Kokkos::deep_copy(num_nodes_per_element, 5UL); + + auto node_state_indices = Kokkos::View("node_state_indices"); + constexpr auto node_state_indices_host_data = + std::array{0UL, 1UL, 2UL, 3UL, 4UL, 5UL, 6UL, 7UL, 8UL, 9UL}; + const auto node_state_indices_host = + Kokkos::View::const_type(node_state_indices_host_data.data() + ); + const auto node_state_indices_mirror = Kokkos::create_mirror(node_state_indices); + Kokkos::deep_copy(node_state_indices_mirror, node_state_indices_host); + Kokkos::deep_copy(node_state_indices, node_state_indices_mirror); + + auto node_freedom_allocation_table = + Kokkos::View("node_freedom_allocation_table"); + Kokkos::deep_copy(node_freedom_allocation_table, FreedomSignature::AllComponents); + + const auto num_non_zero = + ComputeKNumNonZero(num_nodes_per_element, node_state_indices, node_freedom_allocation_table); + + EXPECT_EQ(num_non_zero, 1800UL); } -TEST(ComputeNumberOfNonZeros, TwoElements) { - auto num_nodes = Kokkos::View("num_nodes"); - auto num_nodes_host = Kokkos::create_mirror(num_nodes); - num_nodes_host(0) = size_t{5U}; - num_nodes_host(1) = size_t{3U}; - Kokkos::deep_copy(num_nodes, num_nodes_host); - - auto num_non_zero = size_t{0U}; - Kokkos::parallel_reduce(2, ComputeNumberOfNonZeros{num_nodes}, num_non_zero); - constexpr auto num_dof_elem1 = 5U * 6U; - constexpr auto num_dof_elem2 = 3U * 6U; - constexpr auto expected_num_non_zero = - num_dof_elem1 * num_dof_elem1 + num_dof_elem2 * num_dof_elem2; - EXPECT_EQ(num_non_zero, expected_num_non_zero); +TEST(ComputeKNumNonZero, TwoElementsFiveNodesOverlap) { + auto num_nodes_per_element = Kokkos::View("num_nodes_per_element"); + Kokkos::deep_copy(num_nodes_per_element, 5UL); + + auto node_state_indices = Kokkos::View("node_state_indices"); + constexpr auto node_state_indices_host_data = + std::array{0UL, 1UL, 2UL, 3UL, 4UL, 4UL, 5UL, 6UL, 7UL, 8UL}; + const auto node_state_indices_host = + Kokkos::View::const_type(node_state_indices_host_data.data() + ); + const auto node_state_indices_mirror = Kokkos::create_mirror(node_state_indices); + Kokkos::deep_copy(node_state_indices_mirror, node_state_indices_host); + Kokkos::deep_copy(node_state_indices, node_state_indices_mirror); + + auto node_freedom_allocation_table = + Kokkos::View("node_freedom_allocation_table"); + Kokkos::deep_copy(node_freedom_allocation_table, FreedomSignature::AllComponents); + + const auto num_non_zero = + ComputeKNumNonZero(num_nodes_per_element, node_state_indices, node_freedom_allocation_table); + + EXPECT_EQ(num_non_zero, 1764UL); } } // namespace openturbine::tests diff --git a/tests/unit_tests/solver/test_compute_k_row_ptrs.cpp b/tests/unit_tests/solver/test_compute_k_row_ptrs.cpp new file mode 100644 index 00000000..b3f0748f --- /dev/null +++ b/tests/unit_tests/solver/test_compute_k_row_ptrs.cpp @@ -0,0 +1,195 @@ +#include +#include + +#include "src/solver/compute_k_row_ptrs.hpp" + +namespace openturbine::tests { + +TEST(ComputeKRowPtrs, OneElementOneNode) { + constexpr auto K_num_rows = 6U; + + const auto node_freedom_allocation_table = + Kokkos::View("node_freedom_allocation_table"); + Kokkos::deep_copy(node_freedom_allocation_table, FreedomSignature::AllComponents); + + const auto node_freedom_map_table = Kokkos::View("node_freedom_map_table"); + Kokkos::deep_copy(node_freedom_map_table, 0UL); + + const auto num_nodes_per_element = Kokkos::View("num_nodes_per_element"); + Kokkos::deep_copy(num_nodes_per_element, 1UL); + + const auto node_state_indices = Kokkos::View("node_state_indices"); + Kokkos::deep_copy(node_state_indices, 0UL); + + const auto row_ptrs = ComputeKRowPtrs>( + K_num_rows, node_freedom_allocation_table, node_freedom_map_table, num_nodes_per_element, + node_state_indices + ); + + const auto row_ptrs_mirror = Kokkos::create_mirror(row_ptrs); + Kokkos::deep_copy(row_ptrs_mirror, row_ptrs); + + for (auto row = 0U; row < 7U; ++row) { + EXPECT_EQ(row_ptrs_mirror(row), row * 6U); + } +} + +TEST(ComputeKRowPtrs, OneElementTwoNodes) { + constexpr auto K_num_rows = 12U; + + const auto node_freedom_allocation_table = + Kokkos::View("node_freedom_allocation_table"); + Kokkos::deep_copy(node_freedom_allocation_table, FreedomSignature::AllComponents); + + const auto node_freedom_map_table = Kokkos::View("node_freedom_map_table"); + constexpr auto node_freedom_map_table_host_data = std::array{0UL, 6UL}; + const auto node_freedom_map_table_host = Kokkos::View::const_type( + node_freedom_map_table_host_data.data() + ); + const auto node_freedom_map_table_mirror = Kokkos::create_mirror(node_freedom_map_table); + Kokkos::deep_copy(node_freedom_map_table_mirror, node_freedom_map_table_host); + Kokkos::deep_copy(node_freedom_map_table, node_freedom_map_table_mirror); + + const auto num_nodes_per_element = Kokkos::View("num_nodes_per_element"); + Kokkos::deep_copy(num_nodes_per_element, 2UL); + + const auto node_state_indices = Kokkos::View("node_state_indices"); + constexpr auto node_state_indices_host_data = std::array{0UL, 1UL}; + const auto node_state_indices_host = + Kokkos::View::const_type(node_state_indices_host_data.data() + ); + const auto node_state_indices_mirror = Kokkos::create_mirror(node_state_indices); + Kokkos::deep_copy(node_state_indices_mirror, node_state_indices_host); + Kokkos::deep_copy(node_state_indices, node_state_indices_mirror); + + const auto row_ptrs = ComputeKRowPtrs>( + K_num_rows, node_freedom_allocation_table, node_freedom_map_table, num_nodes_per_element, + node_state_indices + ); + + const auto row_ptrs_mirror = Kokkos::create_mirror(row_ptrs); + Kokkos::deep_copy(row_ptrs_mirror, row_ptrs); + + for (auto row = 0U; row < 13U; ++row) { + EXPECT_EQ(row_ptrs_mirror(row), row * 12U); + } +} + +TEST(ComputeKRowPtrs, TwoElementTwoNodesNoOverlap) { + constexpr auto K_num_rows = 24U; + + const auto node_freedom_allocation_table = + Kokkos::View("node_freedom_allocation_table"); + Kokkos::deep_copy(node_freedom_allocation_table, FreedomSignature::AllComponents); + + const auto node_freedom_map_table = Kokkos::View("node_freedom_map_table"); + constexpr auto node_freedom_map_table_host_data = std::array{0UL, 6UL, 12UL, 18UL}; + const auto node_freedom_map_table_host = Kokkos::View::const_type( + node_freedom_map_table_host_data.data() + ); + const auto node_freedom_map_table_mirror = Kokkos::create_mirror(node_freedom_map_table); + Kokkos::deep_copy(node_freedom_map_table_mirror, node_freedom_map_table_host); + Kokkos::deep_copy(node_freedom_map_table, node_freedom_map_table_mirror); + + const auto num_nodes_per_element = Kokkos::View("num_nodes_per_element"); + Kokkos::deep_copy(num_nodes_per_element, 2UL); + + const auto node_state_indices = Kokkos::View("node_state_indices"); + constexpr auto node_state_indices_host_data = std::array{0UL, 1UL, 2UL, 3UL}; + const auto node_state_indices_host = + Kokkos::View::const_type(node_state_indices_host_data.data() + ); + const auto node_state_indices_mirror = Kokkos::create_mirror(node_state_indices); + Kokkos::deep_copy(node_state_indices_mirror, node_state_indices_host); + Kokkos::deep_copy(node_state_indices, node_state_indices_mirror); + + const auto row_ptrs = ComputeKRowPtrs>( + K_num_rows, node_freedom_allocation_table, node_freedom_map_table, num_nodes_per_element, + node_state_indices + ); + + const auto row_ptrs_mirror = Kokkos::create_mirror(row_ptrs); + Kokkos::deep_copy(row_ptrs_mirror, row_ptrs); + + for (auto row = 0U; row < 25U; ++row) { + EXPECT_EQ(row_ptrs_mirror(row), row * 12U); + } +} + +TEST(ComputeKRowPtrs, TwoElementTwoNodesOverlap) { + constexpr auto K_num_rows = 18U; + + const auto node_freedom_allocation_table = + Kokkos::View("node_freedom_allocation_table"); + Kokkos::deep_copy(node_freedom_allocation_table, FreedomSignature::AllComponents); + + const auto node_freedom_map_table = Kokkos::View("node_freedom_map_table"); + constexpr auto node_freedom_map_table_host_data = std::array{0UL, 6UL, 12UL}; + const auto node_freedom_map_table_host = Kokkos::View::const_type( + node_freedom_map_table_host_data.data() + ); + const auto node_freedom_map_table_mirror = Kokkos::create_mirror(node_freedom_map_table); + Kokkos::deep_copy(node_freedom_map_table_mirror, node_freedom_map_table_host); + Kokkos::deep_copy(node_freedom_map_table, node_freedom_map_table_mirror); + + const auto num_nodes_per_element = Kokkos::View("num_nodes_per_element"); + Kokkos::deep_copy(num_nodes_per_element, 2UL); + + const auto node_state_indices = Kokkos::View("node_state_indices"); + constexpr auto node_state_indices_host_data = std::array{0UL, 1UL, 1UL, 2UL}; + const auto node_state_indices_host = + Kokkos::View::const_type(node_state_indices_host_data.data() + ); + const auto node_state_indices_mirror = Kokkos::create_mirror(node_state_indices); + Kokkos::deep_copy(node_state_indices_mirror, node_state_indices_host); + Kokkos::deep_copy(node_state_indices, node_state_indices_mirror); + + const auto row_ptrs = ComputeKRowPtrs>( + K_num_rows, node_freedom_allocation_table, node_freedom_map_table, num_nodes_per_element, + node_state_indices + ); + + const auto row_ptrs_mirror = Kokkos::create_mirror(row_ptrs); + Kokkos::deep_copy(row_ptrs_mirror, row_ptrs); + + for (auto row = 0U; row < 7U; ++row) { + EXPECT_EQ(row_ptrs_mirror(row), row * 12U); + } + + for (auto row = 0U; row < 6U; ++row) { + EXPECT_EQ(row_ptrs_mirror(row + 7U), 72U + (row + 1) * 18U); + } + + for (auto row = 0U; row < 6U; ++row) { + EXPECT_EQ(row_ptrs_mirror(row + 13U), 180U + (row + 1) * 12U); + } +} +/* +TEST(PopulateSparseRowPtrs, TwoElements) { + auto num_nodes = Kokkos::View("num_nodes"); + auto num_nodes_host = Kokkos::create_mirror(num_nodes); + num_nodes_host(0) = size_t{5U}; + num_nodes_host(1) = size_t{3U}; + Kokkos::deep_copy(num_nodes, num_nodes_host); + + constexpr auto elem1_rows = 5 * 6; + constexpr auto elem2_rows = 3 * 6; + constexpr auto num_rows = elem1_rows + elem2_rows; + auto row_ptrs = Kokkos::View("row_ptrs"); + + Kokkos::parallel_for(1, PopulateSparseRowPtrs{num_nodes, row_ptrs}); + + auto row_ptrs_host = Kokkos::create_mirror(row_ptrs); + Kokkos::deep_copy(row_ptrs_host, row_ptrs); + for (int row = 0; row < elem1_rows + 1; ++row) { + EXPECT_EQ(row_ptrs_host(row), elem1_rows * row); + } + for (int row = elem1_rows * (elem1_rows + 1); row < elem1_rows + 1 + elem2_rows; ++row) { + EXPECT_EQ( + row_ptrs_host(row), + elem1_rows * (elem1_rows + 1) + elem2_rows * (row - elem1_rows * (elem1_rows + 1)) + ); + } +} +*/ +} // namespace openturbine::tests diff --git a/tests/unit_tests/solver/test_compute_t_col_inds.cpp b/tests/unit_tests/solver/test_compute_t_col_inds.cpp new file mode 100644 index 00000000..d802f6fc --- /dev/null +++ b/tests/unit_tests/solver/test_compute_t_col_inds.cpp @@ -0,0 +1,118 @@ +#include +#include + +#include "src/solver/compute_t_col_inds.hpp" + +namespace openturbine::tests { + +TEST(ComputeTColInds, OneNode) { + constexpr auto T_num_non_zero = 36UL; + + const auto node_freedom_allocation_table = + Kokkos::View("node_freedom_allocation_table"); + Kokkos::deep_copy(node_freedom_allocation_table, FreedomSignature::AllComponents); + + const auto node_freedom_map_table = Kokkos::View("node_freedom_map_table"); + Kokkos::deep_copy(node_freedom_map_table, 0UL); + + const auto T_row_ptrs = Kokkos::View("T_row_ptrs"); + constexpr auto T_row_ptrs_host_data = std::array{0UL, 6UL, 12UL, 18UL, 24UL, 30UL, 36UL}; + const auto T_row_ptrs_host = + Kokkos::View::const_type(T_row_ptrs_host_data.data()); + const auto T_row_ptrs_mirror = Kokkos::create_mirror(T_row_ptrs); + Kokkos::deep_copy(T_row_ptrs_mirror, T_row_ptrs_host); + Kokkos::deep_copy(T_row_ptrs, T_row_ptrs_mirror); + + const auto col_inds = ComputeTColInds, Kokkos::View>( + T_num_non_zero, node_freedom_allocation_table, node_freedom_map_table, T_row_ptrs + ); + + const auto col_inds_mirror = Kokkos::create_mirror(col_inds); + Kokkos::deep_copy(col_inds_mirror, col_inds); + + for (auto row = 0U; row < 6U; ++row) { + for (auto col = 0U; col < 6U; ++col) { + EXPECT_EQ(col_inds_mirror(row * 6U + col), col); + } + } +} + +TEST(ComputeTColInds, TwoNodes) { + constexpr auto T_num_non_zero = 72UL; + + const auto node_freedom_allocation_table = + Kokkos::View("node_freedom_allocation_table"); + Kokkos::deep_copy(node_freedom_allocation_table, FreedomSignature::AllComponents); + + const auto node_freedom_map_table = Kokkos::View("node_freedom_map_table"); + constexpr auto node_freedom_map_table_host_data = std::array{0UL, 6UL}; + const auto node_freedom_map_table_host = + Kokkos::View::const_type(node_freedom_map_table_host_data.data()); + const auto node_freedom_map_table_mirror = Kokkos::create_mirror(node_freedom_map_table); + Kokkos::deep_copy(node_freedom_map_table_mirror, node_freedom_map_table_host); + Kokkos::deep_copy(node_freedom_map_table, node_freedom_map_table_mirror); + + const auto T_row_ptrs = Kokkos::View("T_row_ptrs"); + constexpr auto T_row_ptrs_host_data = + std::array{0UL, 6UL, 12UL, 18UL, 24UL, 30UL, 36UL, 42UL, 48UL, 54UL, 60UL, 66UL, 72UL}; + const auto T_row_ptrs_host = + Kokkos::View::const_type(T_row_ptrs_host_data.data()); + const auto T_row_ptrs_mirror = Kokkos::create_mirror(T_row_ptrs); + Kokkos::deep_copy(T_row_ptrs_mirror, T_row_ptrs_host); + Kokkos::deep_copy(T_row_ptrs, T_row_ptrs_mirror); + + const auto col_inds = ComputeTColInds, Kokkos::View>( + T_num_non_zero, node_freedom_allocation_table, node_freedom_map_table, T_row_ptrs + ); + + const auto col_inds_mirror = Kokkos::create_mirror(col_inds); + Kokkos::deep_copy(col_inds_mirror, col_inds); + + for (auto row = 0U; row < 6U; ++row) { + for (auto col = 0U; col < 6U; ++col) { + EXPECT_EQ(col_inds_mirror(row * 6U + col), col); + } + } + + for (auto row = 0U; row < 6U; ++row) { + for (auto col = 0U; col < 6U; ++col) { + EXPECT_EQ(col_inds_mirror(36U + row * 6U + col), col + 6U); + } + } +} +/* +TEST(PopulateTangentIndices, TwoNodes) { + constexpr auto num_system_nodes = 2U; + constexpr auto num_system_dofs = num_system_nodes * 6U * 6U; + constexpr auto node_state_indices_host_data = std::array{size_t{0U}, size_t{1U}}; + const auto node_state_indices_host = + Kokkos::View( + node_state_indices_host_data.data() + ); + const auto node_state_indices = Kokkos::View("node_state_indices"); + Kokkos::deep_copy(node_state_indices, node_state_indices_host); + + const auto indices = Kokkos::View("indices"); + + Kokkos::parallel_for( + "PopulateTangentIndices", 1, + PopulateTangentIndices{num_system_nodes, node_state_indices, indices} + ); + + constexpr auto indices_exact_data = std::array{ + 0, 1, 2, 3, 4, 5, 0, 1, 2, 3, 4, 5, 0, 1, 2, 3, 4, 5, 0, 1, 2, 3, 4, 5, + 0, 1, 2, 3, 4, 5, 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 6, 7, 8, 9, 10, 11, + 6, 7, 8, 9, 10, 11, 6, 7, 8, 9, 10, 11, 6, 7, 8, 9, 10, 11, 6, 7, 8, 9, 10, 11, + }; + const auto indices_exact = + Kokkos::View(indices_exact_data.data()); + + const auto indices_host = Kokkos::create_mirror(indices); + Kokkos::deep_copy(indices_host, indices); + + for (auto i = 0U; i < num_system_dofs; ++i) { + EXPECT_EQ(indices_host(i), indices_exact(i)); + } +} +*/ +} // namespace openturbine::tests diff --git a/tests/unit_tests/solver/test_compute_t_num_non_zero.cpp b/tests/unit_tests/solver/test_compute_t_num_non_zero.cpp new file mode 100644 index 00000000..b8527156 --- /dev/null +++ b/tests/unit_tests/solver/test_compute_t_num_non_zero.cpp @@ -0,0 +1,27 @@ +#include +#include + +#include "src/solver/compute_t_num_non_zero.hpp" + +namespace openturbine::tests { +TEST(ComputeTNumNonZero, OneNode) { + const auto node_freedom_allocation_table = + Kokkos::View("node_freedom_allocation_table"); + Kokkos::deep_copy(node_freedom_allocation_table, FreedomSignature::AllComponents); + + const auto num_non_zero = ComputeTNumNonZero(node_freedom_allocation_table); + + EXPECT_EQ(num_non_zero, 36UL); +} + +TEST(ComputeTNumNonZero, TwoNodes) { + const auto node_freedom_allocation_table = + Kokkos::View("node_freedom_allocation_table"); + Kokkos::deep_copy(node_freedom_allocation_table, FreedomSignature::AllComponents); + + const auto num_non_zero = ComputeTNumNonZero(node_freedom_allocation_table); + + EXPECT_EQ(num_non_zero, 72UL); +} + +} // namespace openturbine::tests diff --git a/tests/unit_tests/solver/test_compute_t_row_ptrs.cpp b/tests/unit_tests/solver/test_compute_t_row_ptrs.cpp new file mode 100644 index 00000000..1b14a1a9 --- /dev/null +++ b/tests/unit_tests/solver/test_compute_t_row_ptrs.cpp @@ -0,0 +1,83 @@ +#include +#include + +#include "src/solver/compute_t_row_ptrs.hpp" + +namespace openturbine::tests { + +TEST(ComputeTRowPtrs, OneNode) { + constexpr auto T_num_rows = 6U; + + const auto node_freedom_allocation_table = + Kokkos::View("node_freedom_allocation_table"); + Kokkos::deep_copy(node_freedom_allocation_table, FreedomSignature::AllComponents); + + const auto node_freedom_map_table = Kokkos::View("node_freedom_map_table"); + Kokkos::deep_copy(node_freedom_map_table, 0UL); + + const auto row_ptrs = ComputeTRowPtrs>( + T_num_rows, node_freedom_allocation_table, node_freedom_map_table + ); + + const auto row_ptrs_mirror = Kokkos::create_mirror(row_ptrs); + Kokkos::deep_copy(row_ptrs_mirror, row_ptrs); + + for (auto i = 0U; i < 7U; ++i) { + EXPECT_EQ(row_ptrs_mirror(i), i * 6U); + } +} + +TEST(ComputeTRowPtrs, TwoNodes) { + constexpr auto T_num_rows = 12U; + + const auto node_freedom_allocation_table = + Kokkos::View("node_freedom_allocation_table"); + Kokkos::deep_copy(node_freedom_allocation_table, FreedomSignature::AllComponents); + + const auto node_freedom_map_table = Kokkos::View("node_freedom_map_table"); + constexpr auto node_freedom_map_table_host_data = std::array{0UL, 6UL}; + const auto node_freedom_map_table_host = Kokkos::View::const_type( + node_freedom_map_table_host_data.data() + ); + const auto node_freedom_map_table_mirror = Kokkos::create_mirror(node_freedom_map_table); + Kokkos::deep_copy(node_freedom_map_table_mirror, node_freedom_map_table_host); + Kokkos::deep_copy(node_freedom_map_table, node_freedom_map_table_mirror); + + const auto row_ptrs = ComputeTRowPtrs>( + T_num_rows, node_freedom_allocation_table, node_freedom_map_table + ); + + const auto row_ptrs_mirror = Kokkos::create_mirror(row_ptrs); + Kokkos::deep_copy(row_ptrs_mirror, row_ptrs); + + for (auto i = 0U; i < 7U; ++i) { + EXPECT_EQ(row_ptrs_mirror(i), i * 6U); + } + + for (auto i = 0U; i < 6U; ++i) { + EXPECT_EQ(row_ptrs_mirror(i + 7U), 36U + (i + 1) * 6U); + } +} + +/* +TEST(PopulateTangentRowPtrs, TwoNodeSystem) { + constexpr auto num_system_nodes = 2U; + constexpr auto num_system_dofs = 6U * num_system_nodes; + constexpr auto num_entries = num_system_dofs + 1U; + const auto row_ptrs = Kokkos::View("row_ptrs"); + Kokkos::parallel_for( + "PopulateTangentRowPtrs", 1, PopulateTangentRowPtrs{num_system_nodes, row_ptrs} + ); + + constexpr auto exact_row_ptrs_data = + std::array{0U, 6U, 12U, 18U, 24U, 30U, 36U, 42U, 48U, 54U, 60U, 66U, 72U, 78U}; + const auto exact_row_ptrs = + Kokkos::View(exact_row_ptrs_data.data()); + const auto row_ptrs_host = Kokkos::create_mirror(row_ptrs); + Kokkos::deep_copy(row_ptrs_host, row_ptrs); + for (auto i = 0U; i < num_entries; ++i) { + EXPECT_EQ(row_ptrs_host(i), exact_row_ptrs(i)); + } +} +*/ +} // namespace openturbine::tests diff --git a/tests/unit_tests/solver/test_create_k_col_inds.cpp b/tests/unit_tests/solver/test_create_k_col_inds.cpp deleted file mode 100644 index d583031f..00000000 --- a/tests/unit_tests/solver/test_create_k_col_inds.cpp +++ /dev/null @@ -1,83 +0,0 @@ -#include -#include - -#include "src/solver/populate_sparse_indices.hpp" - -namespace openturbine::tests { - -TEST(PopulateSparseIndices, SingleElement) { - constexpr auto num_nodes = size_t{5U}; - constexpr auto num_dof = num_nodes * 6U; - - auto num_nodes_per_element = Kokkos::View("num_nodes_per_element"); - auto num_nodes_per_element_host = Kokkos::create_mirror(num_nodes_per_element); - num_nodes_per_element_host(0) = num_nodes; - Kokkos::deep_copy(num_nodes_per_element, num_nodes_per_element_host); - - auto node_state_indices = Kokkos::View("node_state_indices"); - auto node_state_indices_host_data = std::array{0U, 1U, 2U, 3U, 4U}; - auto node_state_indices_host = - Kokkos::View(node_state_indices_host_data.data()); - auto node_state_indices_mirror = Kokkos::create_mirror(node_state_indices); - Kokkos::deep_copy(node_state_indices_mirror, node_state_indices_host); - Kokkos::deep_copy(node_state_indices, node_state_indices_mirror); - auto indices = Kokkos::View("indices"); - - Kokkos::parallel_for( - 1, PopulateSparseIndices{num_nodes_per_element, node_state_indices, indices} - ); - - auto indices_host = Kokkos::create_mirror(indices); - Kokkos::deep_copy(indices_host, indices); - for (auto row = 0U; row < num_dof; ++row) { - for (auto column = 0U; column < num_dof; ++column) { - ASSERT_EQ(indices_host(row * num_dof + column), column); - } - } -} - -TEST(PopulateSparseIndices, TwoElements) { - constexpr auto elem1_num_nodes = size_t{5U}; - constexpr auto elem1_num_dof = elem1_num_nodes * 6U; - constexpr auto elem2_num_nodes = size_t{3U}; - constexpr auto elem2_num_dof = elem2_num_nodes * 6U; - constexpr auto max_num_nodes = elem1_num_nodes; - - auto num_nodes_per_element = Kokkos::View("num_nodes_per_element"); - auto num_nodes_per_element_host = Kokkos::create_mirror(num_nodes_per_element); - num_nodes_per_element_host(0) = elem1_num_nodes; - num_nodes_per_element_host(1) = elem2_num_nodes; - Kokkos::deep_copy(num_nodes_per_element, num_nodes_per_element_host); - - auto node_state_indices = Kokkos::View("node_state_indices"); - auto node_state_indices_host_data = - std::array{0U, 1U, 2U, 3U, 4U, 5U, 6U, 7U, 0U, 0U}; - auto node_state_indices_host = - Kokkos::View(node_state_indices_host_data.data() - ); - auto node_state_indices_mirror = Kokkos::create_mirror(node_state_indices); - Kokkos::deep_copy(node_state_indices_mirror, node_state_indices_host); - Kokkos::deep_copy(node_state_indices, node_state_indices_mirror); - auto indices = - Kokkos::View("indices"); - - Kokkos::parallel_for( - 1, PopulateSparseIndices{num_nodes_per_element, node_state_indices, indices} - ); - - auto indices_host = Kokkos::create_mirror(indices); - Kokkos::deep_copy(indices_host, indices); - for (auto row = 0U; row < elem1_num_dof; ++row) { - for (auto column = 0U; column < elem1_num_dof; ++column) { - EXPECT_EQ(indices_host(row * elem1_num_dof + column), column); - } - } - - for (auto row = 0U; row < elem2_num_dof; ++row) { - for (auto column = 0U; column < elem2_num_dof; ++column) { - const auto index = elem1_num_dof * elem1_num_dof + row * elem2_num_dof + column; - EXPECT_EQ(indices_host(index), column + elem1_num_dof); - } - } -} -} // namespace openturbine::tests diff --git a/tests/unit_tests/solver/test_create_k_row_ptrs.cpp b/tests/unit_tests/solver/test_create_k_row_ptrs.cpp deleted file mode 100644 index 7eba7bd9..00000000 --- a/tests/unit_tests/solver/test_create_k_row_ptrs.cpp +++ /dev/null @@ -1,53 +0,0 @@ -#include -#include - -#include "src/solver/populate_sparse_row_ptrs.hpp" - -namespace openturbine::tests { - -TEST(PopulateSparseRowPtrs, SingleElement) { - auto num_nodes = Kokkos::View("num_nodes"); - auto num_nodes_host = Kokkos::create_mirror(num_nodes); - num_nodes_host(0) = size_t{5U}; - Kokkos::deep_copy(num_nodes, num_nodes_host); - - constexpr auto num_rows = 5 * 6; - auto row_ptrs = Kokkos::View("row_ptrs"); - - Kokkos::parallel_for(1, PopulateSparseRowPtrs{num_nodes, row_ptrs}); - - auto row_ptrs_host = Kokkos::create_mirror(row_ptrs); - Kokkos::deep_copy(row_ptrs_host, row_ptrs); - for (int row = 0; row < num_rows + 1; ++row) { - EXPECT_EQ(row_ptrs_host(row), num_rows * row); - } -} - -TEST(PopulateSparseRowPtrs, TwoElements) { - auto num_nodes = Kokkos::View("num_nodes"); - auto num_nodes_host = Kokkos::create_mirror(num_nodes); - num_nodes_host(0) = size_t{5U}; - num_nodes_host(1) = size_t{3U}; - Kokkos::deep_copy(num_nodes, num_nodes_host); - - constexpr auto elem1_rows = 5 * 6; - constexpr auto elem2_rows = 3 * 6; - constexpr auto num_rows = elem1_rows + elem2_rows; - auto row_ptrs = Kokkos::View("row_ptrs"); - - Kokkos::parallel_for(1, PopulateSparseRowPtrs{num_nodes, row_ptrs}); - - auto row_ptrs_host = Kokkos::create_mirror(row_ptrs); - Kokkos::deep_copy(row_ptrs_host, row_ptrs); - for (int row = 0; row < elem1_rows + 1; ++row) { - EXPECT_EQ(row_ptrs_host(row), elem1_rows * row); - } - for (int row = elem1_rows * (elem1_rows + 1); row < elem1_rows + 1 + elem2_rows; ++row) { - EXPECT_EQ( - row_ptrs_host(row), - elem1_rows * (elem1_rows + 1) + elem2_rows * (row - elem1_rows * (elem1_rows + 1)) - ); - } -} - -} // namespace openturbine::tests diff --git a/tests/unit_tests/solver/test_create_t_col_inds.cpp b/tests/unit_tests/solver/test_create_t_col_inds.cpp deleted file mode 100644 index c79bf1cf..00000000 --- a/tests/unit_tests/solver/test_create_t_col_inds.cpp +++ /dev/null @@ -1,42 +0,0 @@ -#include -#include - -#include "src/solver/populate_tangent_indices.hpp" - -namespace openturbine::tests { - -TEST(PopulateTangentIndices, TwoNodes) { - constexpr auto num_system_nodes = 2U; - constexpr auto num_system_dofs = num_system_nodes * 6U * 6U; - constexpr auto node_state_indices_host_data = std::array{size_t{0U}, size_t{1U}}; - const auto node_state_indices_host = - Kokkos::View( - node_state_indices_host_data.data() - ); - const auto node_state_indices = Kokkos::View("node_state_indices"); - Kokkos::deep_copy(node_state_indices, node_state_indices_host); - - const auto indices = Kokkos::View("indices"); - - Kokkos::parallel_for( - "PopulateTangentIndices", 1, - PopulateTangentIndices{num_system_nodes, node_state_indices, indices} - ); - - constexpr auto indices_exact_data = std::array{ - 0, 1, 2, 3, 4, 5, 0, 1, 2, 3, 4, 5, 0, 1, 2, 3, 4, 5, 0, 1, 2, 3, 4, 5, - 0, 1, 2, 3, 4, 5, 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 6, 7, 8, 9, 10, 11, - 6, 7, 8, 9, 10, 11, 6, 7, 8, 9, 10, 11, 6, 7, 8, 9, 10, 11, 6, 7, 8, 9, 10, 11, - }; - const auto indices_exact = - Kokkos::View(indices_exact_data.data()); - - const auto indices_host = Kokkos::create_mirror(indices); - Kokkos::deep_copy(indices_host, indices); - - for (auto i = 0U; i < num_system_dofs; ++i) { - EXPECT_EQ(indices_host(i), indices_exact(i)); - } -} - -} // namespace openturbine::tests diff --git a/tests/unit_tests/solver/test_create_t_row_ptrs.cpp b/tests/unit_tests/solver/test_create_t_row_ptrs.cpp deleted file mode 100644 index a70241d1..00000000 --- a/tests/unit_tests/solver/test_create_t_row_ptrs.cpp +++ /dev/null @@ -1,28 +0,0 @@ -#include -#include - -#include "src/solver/populate_tangent_row_ptrs.hpp" - -namespace openturbine::tests { - -TEST(PopulateTangentRowPtrs, TwoNodeSystem) { - constexpr auto num_system_nodes = 2U; - constexpr auto num_system_dofs = 6U * num_system_nodes; - constexpr auto num_entries = num_system_dofs + 1U; - const auto row_ptrs = Kokkos::View("row_ptrs"); - Kokkos::parallel_for( - "PopulateTangentRowPtrs", 1, PopulateTangentRowPtrs{num_system_nodes, row_ptrs} - ); - - constexpr auto exact_row_ptrs_data = - std::array{0U, 6U, 12U, 18U, 24U, 30U, 36U, 42U, 48U, 54U, 60U, 66U, 72U, 78U}; - const auto exact_row_ptrs = - Kokkos::View(exact_row_ptrs_data.data()); - const auto row_ptrs_host = Kokkos::create_mirror(row_ptrs); - Kokkos::deep_copy(row_ptrs_host, row_ptrs); - for (auto i = 0U; i < num_entries; ++i) { - EXPECT_EQ(row_ptrs_host(i), exact_row_ptrs(i)); - } -} - -} // namespace openturbine::tests From 91e1d1c25129403d73ff318b4ec0734afb525b18 Mon Sep 17 00:00:00 2001 From: dcdemen Date: Mon, 23 Dec 2024 16:42:00 -0700 Subject: [PATCH 3/3] clang-format --- src/solver/compute_k_col_inds.hpp | 3 ++- src/solver/compute_k_num_non_zero.hpp | 3 ++- src/solver/compute_k_row_ptrs.hpp | 3 ++- src/solver/compute_t_col_inds.hpp | 3 ++- src/solver/compute_t_row_ptrs.hpp | 3 ++- src/solver/create_b_matrix.hpp | 6 ++++-- src/solver/create_b_t_matrix.hpp | 9 ++++++--- src/solver/create_k_matrix.hpp | 3 ++- src/solver/create_system_matrix_full.hpp | 3 ++- src/solver/create_t_matrix.hpp | 3 ++- src/solver/create_transpose_matrix_full.hpp | 3 ++- tests/unit_tests/solver/test_compute_k_col_inds.cpp | 5 +++-- 12 files changed, 31 insertions(+), 16 deletions(-) diff --git a/src/solver/compute_k_col_inds.hpp b/src/solver/compute_k_col_inds.hpp index e27f42be..0051a10f 100644 --- a/src/solver/compute_k_col_inds.hpp +++ b/src/solver/compute_k_col_inds.hpp @@ -81,7 +81,8 @@ template "ComputeKColInds", node_freedom_allocation_table.extent(0), ComputeKColIndsFunction{ node_freedom_allocation_table, node_freedom_map_table, num_nodes_per_element, - node_state_indices, K_row_ptrs, K_col_inds} + node_state_indices, K_row_ptrs, K_col_inds + } ); return K_col_inds; diff --git a/src/solver/compute_k_num_non_zero.hpp b/src/solver/compute_k_num_non_zero.hpp index 3bfdb485..b147d6a2 100644 --- a/src/solver/compute_k_num_non_zero.hpp +++ b/src/solver/compute_k_num_non_zero.hpp @@ -51,7 +51,8 @@ struct ComputeKNumNonZero_Diagonal { Kokkos::parallel_reduce( "ComputeKNumNonZero_OffDiagonal", num_nodes_per_element.extent(0), ComputeKNumNonZero_OffDiagonal{ - num_nodes_per_element, node_state_indices, node_freedom_allocation_table}, + num_nodes_per_element, node_state_indices, node_freedom_allocation_table + }, K_num_non_zero_off_diagonal ); auto K_num_non_zero_diagonal = 0UL; diff --git a/src/solver/compute_k_row_ptrs.hpp b/src/solver/compute_k_row_ptrs.hpp index 3ed8ed7c..443895c4 100644 --- a/src/solver/compute_k_row_ptrs.hpp +++ b/src/solver/compute_k_row_ptrs.hpp @@ -74,7 +74,8 @@ template "ComputeKRowEntries", node_freedom_allocation_table.extent(0), ComputeKRowEntries{ node_freedom_allocation_table, node_freedom_map_table, num_nodes_per_element, - node_state_indices, K_row_entries} + node_state_indices, K_row_entries + } ); auto result = 0UL; diff --git a/src/solver/compute_t_col_inds.hpp b/src/solver/compute_t_col_inds.hpp index 145e3361..6e32c7e2 100644 --- a/src/solver/compute_t_col_inds.hpp +++ b/src/solver/compute_t_col_inds.hpp @@ -41,7 +41,8 @@ template Kokkos::parallel_for( "ComputeTColInds", node_freedom_allocation_table.extent(0), ComputeTColIndsFunction{ - node_freedom_allocation_table, node_freedom_map_table, T_row_ptrs, T_col_inds} + node_freedom_allocation_table, node_freedom_map_table, T_row_ptrs, T_col_inds + } ); return T_col_inds; diff --git a/src/solver/compute_t_row_ptrs.hpp b/src/solver/compute_t_row_ptrs.hpp index 9fa3e085..a109102d 100644 --- a/src/solver/compute_t_row_ptrs.hpp +++ b/src/solver/compute_t_row_ptrs.hpp @@ -51,7 +51,8 @@ template Kokkos::parallel_for( "ComputeTRowEntries", node_freedom_allocation_table.extent(0), ComputeTRowEntries{ - node_freedom_allocation_table, node_freedom_map_table, T_row_entries} + node_freedom_allocation_table, node_freedom_map_table, T_row_entries + } ); auto result = 0UL; diff --git a/src/solver/create_b_matrix.hpp b/src/solver/create_b_matrix.hpp index ebf883fc..0e996167 100644 --- a/src/solver/create_b_matrix.hpp +++ b/src/solver/create_b_matrix.hpp @@ -31,7 +31,8 @@ template "PopulateSparseRowPtrsColInds_Constraints", 1, PopulateSparseRowPtrsColInds_Constraints{ constraint_type, constraint_base_node_freedom_table, - constraint_target_node_freedom_table, constraint_row_range, B_row_ptrs, B_col_ind} + constraint_target_node_freedom_table, constraint_row_range, B_row_ptrs, B_col_ind + } ); const auto B_values = ValuesType("B values", B_num_non_zero); KokkosSparse::sort_crs_matrix(B_row_ptrs, B_col_ind, B_values); @@ -42,6 +43,7 @@ template B_num_non_zero, B_values, B_row_ptrs, - B_col_ind}; + B_col_ind + }; } } // namespace openturbine diff --git a/src/solver/create_b_t_matrix.hpp b/src/solver/create_b_t_matrix.hpp index 373d3c28..b402c379 100644 --- a/src/solver/create_b_t_matrix.hpp +++ b/src/solver/create_b_t_matrix.hpp @@ -32,7 +32,8 @@ template "PopulateSparseRowPtrsColInds_Constraints", 1, PopulateSparseRowPtrsColInds_Constraints{ constraint_type, constraint_base_node_freedom_table, - constraint_target_node_freedom_table, constraint_row_range, B_row_ptrs, B_col_ind} + constraint_target_node_freedom_table, constraint_row_range, B_row_ptrs, B_col_ind + } ); const auto B_values = ValuesType("B values", B_num_non_zero); KokkosSparse::sort_crs_matrix(B_row_ptrs, B_col_ind, B_values); @@ -50,7 +51,8 @@ template "PopulateSparseRowPtrsColInds_Transpose", 1, PopulateSparseRowPtrsColInds_Transpose{ B_num_rows, B_num_columns, B_row_ptrs, B_col_ind, col_count, tmp_row_ptrs, B_t_row_ptrs, - B_t_col_inds} + B_t_col_inds + } ); KokkosSparse::sort_crs_matrix(B_t_row_ptrs, B_t_col_inds, B_t_values); return { @@ -60,6 +62,7 @@ template B_t_num_non_zero, B_t_values, B_t_row_ptrs, - B_t_col_inds}; + B_t_col_inds + }; } } // namespace openturbine diff --git a/src/solver/create_k_matrix.hpp b/src/solver/create_k_matrix.hpp index e9c284de..04991558 100644 --- a/src/solver/create_k_matrix.hpp +++ b/src/solver/create_k_matrix.hpp @@ -47,6 +47,7 @@ template K_num_non_zero, K_values, K_row_ptrs, - K_col_inds}; + K_col_inds + }; } } // namespace openturbine diff --git a/src/solver/create_system_matrix_full.hpp b/src/solver/create_system_matrix_full.hpp index 321d6a6b..74744110 100644 --- a/src/solver/create_system_matrix_full.hpp +++ b/src/solver/create_system_matrix_full.hpp @@ -17,7 +17,8 @@ template Kokkos::parallel_for( "FillUnshiftedRowPtrs", num_dofs + 1, FillUnshiftedRowPtrs{ - num_system_dofs, system_matrix.graph.row_map, system_matrix_full_row_ptrs} + num_system_dofs, system_matrix.graph.row_map, system_matrix_full_row_ptrs + } ); return CrsMatrixType( "system_matrix_full", static_cast(num_dofs), static_cast(num_dofs), diff --git a/src/solver/create_t_matrix.hpp b/src/solver/create_t_matrix.hpp index acfda8b8..25910b3c 100644 --- a/src/solver/create_t_matrix.hpp +++ b/src/solver/create_t_matrix.hpp @@ -42,6 +42,7 @@ template T_num_non_zero, T_values, T_row_ptrs, - T_col_inds}; + T_col_inds + }; } } // namespace openturbine diff --git a/src/solver/create_transpose_matrix_full.hpp b/src/solver/create_transpose_matrix_full.hpp index 2006cd34..50a34e9c 100644 --- a/src/solver/create_transpose_matrix_full.hpp +++ b/src/solver/create_transpose_matrix_full.hpp @@ -18,7 +18,8 @@ template Kokkos::parallel_for( "FillUnshiftedRowPtrs", num_dofs + 1, FillUnshiftedRowPtrs{ - num_system_dofs, B_t.graph.row_map, transpose_matrix_full_row_ptrs} + num_system_dofs, B_t.graph.row_map, transpose_matrix_full_row_ptrs + } ); auto transpose_matrix_full_indices = IndicesType("transpose_matrix_full_indices", B_t.nnz()); Kokkos::deep_copy(transpose_matrix_full_indices, static_cast(num_system_dofs)); diff --git a/tests/unit_tests/solver/test_compute_k_col_inds.cpp b/tests/unit_tests/solver/test_compute_k_col_inds.cpp index 2c655a65..94713894 100644 --- a/tests/unit_tests/solver/test_compute_k_col_inds.cpp +++ b/tests/unit_tests/solver/test_compute_k_col_inds.cpp @@ -193,8 +193,9 @@ TEST(ComputeKColInds, TwoElementsTwoNodesOverlap) { Kokkos::deep_copy(node_state_indices, node_state_indices_mirror); const auto K_row_ptrs = Kokkos::View("K_row_ptrs"); - constexpr auto K_row_ptrs_host_data = std::array{ - 0, 12, 24, 36, 48, 60, 72, 90, 108, 126, 144, 162, 180, 192, 204, 216, 228, 240, 252}; + constexpr auto K_row_ptrs_host_data = + std::array{0, 12, 24, 36, 48, 60, 72, 90, 108, 126, + 144, 162, 180, 192, 204, 216, 228, 240, 252}; const auto K_row_ptrs_host = Kokkos::View::const_type(K_row_ptrs_host_data.data()); const auto K_row_ptrs_mirror = Kokkos::create_mirror(K_row_ptrs);