40#ifndef TPETRA_DETAILS_CRSUTILS_HPP
41#define TPETRA_DETAILS_CRSUTILS_HPP
45#include "TpetraCore_config.h"
46#include "Kokkos_Core.hpp"
48#include "Tpetra_Details_CrsPadding.hpp"
49#include "Tpetra_Details_WrappedDualView.hpp"
52#include <unordered_map>
65template<
class ViewType>
67make_uninitialized_view(
68 const std::string& name,
71 const std::string*
const prefix)
74 std::ostringstream os;
75 os << *prefix <<
"Allocate Kokkos::View " << name
76 <<
": " << size << std::endl;
77 std::cerr << os.str();
79 using Kokkos::view_alloc;
80 using Kokkos::WithoutInitializing;
81 return ViewType(view_alloc(name, WithoutInitializing), size);
84template<
class ViewType>
87 const std::string& name,
90 const std::string*
const prefix)
93 std::ostringstream os;
94 os << *prefix <<
"Allocate & initialize Kokkos::View "
95 << name <<
": " << size << std::endl;
96 std::cerr << os.str();
98 return ViewType(name, size);
101template<
class OutViewType,
class InViewType>
103assign_to_view(OutViewType& out,
104 const InViewType& in,
105 const char viewName[],
107 const std::string*
const prefix)
110 std::ostringstream os;
111 os << *prefix <<
"Assign to Kokkos::View " << viewName
112 <<
": Old size: " << out.extent(0)
113 <<
", New size: " << in.extent(0) << std::endl;
114 std::cerr << os.str();
119template<
class MemorySpace,
class ViewType>
120auto create_mirror_view(
121 const MemorySpace& memSpace,
122 const ViewType& view,
124 const std::string*
const prefix) ->
125 decltype(Kokkos::create_mirror_view(memSpace, view))
128 std::ostringstream os;
129 os << *prefix <<
"Create mirror view: "
130 <<
"view.extent(0): " << view.extent(0) << std::endl;
131 std::cerr << os.str();
133 return Kokkos::create_mirror_view(memSpace, view);
136enum class PadCrsAction {
149template<
class RowPtr,
class Indices,
class Values,
class Padding>
152 const PadCrsAction
action,
161 using execution_space =
typename Indices::t_dev::execution_space;
162 using Kokkos::view_alloc;
163 using Kokkos::WithoutInitializing;
165 std::unique_ptr<std::string>
prefix;
170 std::ostringstream
os;
171 os <<
"Proc " <<
my_rank <<
": Tpetra::...::pad_crs_arrays: ";
172 prefix = std::unique_ptr<std::string>(
new std::string(
os.str()));
174 std::cerr <<
os.str();
179 std::ostringstream
os;
195 <<
", values.extent(0): " <<
values_wdv.extent(0)
199 std::cerr <<
os.str();
204 std::ostringstream
os;
205 os << *
prefix <<
"Done; local matrix has no rows" <<
endl;
206 std::cerr <<
os.str();
216 std::ostringstream
os;
217 os << *
prefix <<
"Fill newAllocPerRow & compute increase" <<
endl;
218 std::cerr <<
os.str();
235 Kokkos::DefaultHostExecutionSpace,
size_t>;
236 Kokkos::parallel_reduce
237 (
"Tpetra::CrsGraph: Compute new allocation size per row",
264 std::ostringstream
os;
269 std::cerr <<
os.str();
280 typename Indices::t_dev::non_const_value_type;
287 "Tpetra::CrsGraph column indices",
newIndsSize, verbose,
292 if (
action == PadCrsAction::INDICES_AND_VALUES) {
301 std::ostringstream
os;
303 std::cerr <<
os.str();
306 using range_type = Kokkos::RangePolicy<execution_space, size_t>;
307 Kokkos::parallel_scan(
308 "Tpetra::CrsGraph or CrsMatrix repack",
324 const Kokkos::pair<size_t, size_t>
oldRange(
326 const Kokkos::pair<size_t, size_t>
newRange(
334 if (
action == PadCrsAction::INDICES_AND_VALUES) {
353 std::ostringstream
os;
373 std::cout <<
os.str();
383 std::ostringstream
os;
394 std::ostringstream
os;
396 std::cerr <<
os.str();
401template <
class Po
inters,
class InOutIndices,
class InIndices,
class IndexMap>
404 typename Pointers::value_type
const row,
410 std::function<
void(
size_t const,
size_t const,
size_t const)>
cb)
418 return Teuchos::OrdinalTraits<size_t>::invalid();
421 using offset_type =
typename std::decay<
decltype (
row_ptrs[0])>::type;
422 using ordinal_type =
typename std::decay<
decltype (
cur_indices[0])>::type;
424 const offset_type start =
row_ptrs[row];
450 return Teuchos::OrdinalTraits<size_t>::invalid();
478 return Teuchos::OrdinalTraits<size_t>::invalid();
500template <
class Po
inters,
class Indices1,
class Indices2,
class IndexMap,
class Callback>
503 typename Pointers::value_type
const row,
515 typename std::remove_const<typename Indices1::value_type>::type;
518 const size_t start =
static_cast<size_t> (
row_ptrs[row]);
559template<
class RowPtr,
class Indices,
class Padding>
569 using impl::pad_crs_arrays;
577template<
class RowPtr,
class Indices,
class Values,
class Padding>
588 using impl::pad_crs_arrays;
639template <
class Po
inters,
class InOutIndices,
class InIndices>
642 typename Pointers::value_type
const row,
647 std::function<
void(
const size_t,
const size_t,
const size_t)>
cb =
648 std::function<
void(
const size_t,
const size_t,
const size_t)>())
650 static_assert(std::is_same<typename std::remove_const<typename InOutIndices::value_type>::type,
651 typename std::remove_const<typename InIndices::value_type>::type>::value,
652 "Expected views to have same value type");
655 using ordinal =
typename InOutIndices::value_type;
661template <
class Po
inters,
class InOutIndices,
class InIndices>
664 typename Pointers::value_type
const row,
669 std::function<
typename InOutIndices::value_type(
const typename InIndices::value_type)>
map,
670 std::function<
void(
const size_t,
const size_t,
const size_t)>
cb =
671 std::function<
void(
const size_t,
const size_t,
const size_t)>())
708template <
class Po
inters,
class Indices1,
class Indices2,
class Callback>
711 typename Pointers::value_type
const row,
718 static_assert(std::is_same<typename std::remove_const<typename Indices1::value_type>::type,
719 typename std::remove_const<typename Indices2::value_type>::type>::value,
720 "Expected views to have same value type");
722 using ordinal =
typename Indices2::value_type;
728template <
class Po
inters,
class Indices1,
class Indices2,
class IndexMap,
class Callback>
731 typename Pointers::value_type
const row,
Declaration of Tpetra::Details::Behavior, a class that describes Tpetra's behavior.
void pad_crs_arrays(const PadCrsAction action, const RowPtr &row_ptr_beg, const RowPtr &row_ptr_end, Indices &indices_wdv, Values &values_wdv, const Padding &padding, const int my_rank, const bool verbose)
Implementation of padCrsArrays.
size_t find_crs_indices(typename Pointers::value_type const row, Pointers const &row_ptrs, const size_t curNumEntries, Indices1 const &cur_indices, Indices2 const &new_indices, IndexMap &&map, Callback &&cb)
Implementation of findCrsIndices.
size_t insert_crs_indices(typename Pointers::value_type const row, Pointers const &row_ptrs, InOutIndices &cur_indices, size_t &num_assigned, InIndices const &new_indices, IndexMap &&map, std::function< void(size_t const, size_t const, size_t const)> cb)
Implementation of insertCrsIndices.
Struct that holds views of the contents of a CrsMatrix.
static size_t verbosePrintCountThreshold()
Number of entries below which arrays, lists, etc. will be printed in debug mode.
Implementation details of Tpetra.
void padCrsArrays(const RowPtr &rowPtrBeg, const RowPtr &rowPtrEnd, Indices &indices_wdv, const Padding &padding, const int my_rank, const bool verbose)
Determine if the row pointers and indices arrays need to be resized to accommodate new entries....
void verbosePrintArray(std::ostream &out, const ArrayType &x, const char name[], const size_t maxNumToPrint)
Print min(x.size(), maxNumToPrint) entries of x.
size_t insertCrsIndices(typename Pointers::value_type const row, Pointers const &rowPtrs, InOutIndices &curIndices, size_t &numAssigned, InIndices const &newIndices, std::function< void(const size_t, const size_t, const size_t)> cb=std::function< void(const size_t, const size_t, const size_t)>())
Insert new indices in to current list of indices.
size_t findCrsIndices(typename Pointers::value_type const row, Pointers const &rowPtrs, const size_t curNumEntries, Indices1 const &curIndices, Indices2 const &newIndices, Callback &&cb)
Finds offsets in to current list of indices.
Namespace Tpetra contains the class and methods constituting the Tpetra library.