19#ifndef OPM_COLORING_AND_REORDERING_UTILS_HPP
20#define OPM_COLORING_AND_REORDERING_UTILS_HPP
24#include <opm/common/ErrorMacros.hpp>
25#include <opm/grid/utility/SparseTable.hpp>
26#include <opm/simulators/linalg/gpuistl/detail/safe_conversion.hpp>
36inline std::vector<int>
37createReorderedToNatural(
const Opm::SparseTable<size_t>& levelSets)
41 for (
auto row : levelSets) {
42 for (
auto col : row) {
44 fmt::format(
"Internal error. globCnt = {}, res.size() = {}", globCnt, res.size()));
45 res[globCnt++] =
static_cast<int>(col);
51inline std::vector<int>
52createNaturalToReordered(
const Opm::SparseTable<size_t>& levelSets)
56 for (
auto row : levelSets) {
57 for (
auto col : row) {
59 fmt::format(
"Internal error. globCnt = {}, res.size() = {}", globCnt, res.size()));
66template <
class M,
class field_type,
class GPUM>
67inline std::unique_ptr<GPUM>
68createReorderedMatrix(
const M& naturalMatrix,
const std::vector<int>& reorderedToNatural)
70 M reorderedMatrix(naturalMatrix.N(), naturalMatrix.N(), naturalMatrix.nonzeroes(), M::row_wise);
71 for (
auto dstRowIt = reorderedMatrix.createbegin(); dstRowIt != reorderedMatrix.createend(); ++dstRowIt) {
72 auto srcRow = naturalMatrix.begin() + reorderedToNatural[dstRowIt.index()];
74 for (
auto elem = srcRow->begin(); elem != srcRow->end(); elem++) {
75 dstRowIt.insert(elem.index());
78 return std::unique_ptr<GPUM>(
new auto(GPUM::fromMatrix(reorderedMatrix,
true)));
81template <
class M,
class field_type,
class GPUM>
82inline std::tuple<std::unique_ptr<GPUM>, std::unique_ptr<GPUM>>
83extractLowerAndUpperMatrices(
const M& naturalMatrix,
const std::vector<int>& reorderedToNatural)
85 const size_t new_nnz = (naturalMatrix.nonzeroes() - naturalMatrix.N()) / 2;
87 M reorderedLower(naturalMatrix.N(), naturalMatrix.N(), new_nnz, M::row_wise);
88 M reorderedUpper(naturalMatrix.N(), naturalMatrix.N(), new_nnz, M::row_wise);
90 for (
auto lowerIt = reorderedLower.createbegin(), upperIt = reorderedUpper.createbegin();
91 lowerIt != reorderedLower.createend();
92 ++lowerIt, ++upperIt) {
94 auto srcRow = naturalMatrix.begin() + reorderedToNatural[lowerIt.index()];
96 for (
auto elem = srcRow->begin(); elem != srcRow->end(); ++elem) {
97 if (elem.index() < srcRow.index()) {
98 lowerIt.insert(elem.index());
99 }
else if (elem.index() > srcRow.index()) {
100 upperIt.insert(elem.index());
105 return {std::unique_ptr<GPUM>(
new auto(GPUM::fromMatrix(reorderedLower,
true))),
106 std::unique_ptr<GPUM>(
new auto(GPUM::fromMatrix(reorderedUpper,
true)))};
Contains wrappers to make the CuBLAS library behave as a modern C++ library with function overlading.
Definition autotuner.hpp:29
__host__ __device__ std::size_t to_size_t(int i)
to_size_t converts a (on most relevant platforms) a 32 bit signed int to a 64 bits unsigned int
Definition safe_conversion.hpp:86