Tpetra parallel linear algebra Version of the Day
Loading...
Searching...
No Matches
Tpetra_Details_crsUtils.hpp
Go to the documentation of this file.
1// @HEADER
2// ***********************************************************************
3//
4// Tpetra: Templated Linear Algebra Services Package
5// Copyright (2008) Sandia Corporation
6//
7// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8// the U.S. Government retains certain rights in this software.
9//
10// Redistribution and use in source and binary forms, with or without
11// modification, are permitted provided that the following conditions are
12// met:
13//
14// 1. Redistributions of source code must retain the above copyright
15// notice, this list of conditions and the following disclaimer.
16//
17// 2. Redistributions in binary form must reproduce the above copyright
18// notice, this list of conditions and the following disclaimer in the
19// documentation and/or other materials provided with the distribution.
20//
21// 3. Neither the name of the Corporation nor the names of the
22// contributors may be used to endorse or promote products derived from
23// this software without specific prior written permission.
24//
25// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36//
37// ************************************************************************
38// @HEADER
39
40#ifndef TPETRA_DETAILS_CRSUTILS_HPP
41#define TPETRA_DETAILS_CRSUTILS_HPP
42#include <numeric>
43#include <type_traits>
44
45#include "TpetraCore_config.h"
46#include "Kokkos_Core.hpp"
48#include "Tpetra_Details_CrsPadding.hpp"
49#include "Tpetra_Details_WrappedDualView.hpp"
50#include <iostream>
51#include <memory>
52#include <unordered_map>
53
59
60namespace Tpetra {
61namespace Details {
62
63namespace impl {
64
65template<class ViewType>
66ViewType
67make_uninitialized_view(
68 const std::string& name,
69 const size_t size,
70 const bool verbose,
71 const std::string* const prefix)
72{
73 if (verbose) {
74 std::ostringstream os;
75 os << *prefix << "Allocate Kokkos::View " << name
76 << ": " << size << std::endl;
77 std::cerr << os.str();
78 }
79 using Kokkos::view_alloc;
80 using Kokkos::WithoutInitializing;
81 return ViewType(view_alloc(name, WithoutInitializing), size);
82}
83
84template<class ViewType>
85ViewType
86make_initialized_view(
87 const std::string& name,
88 const size_t size,
89 const bool verbose,
90 const std::string* const prefix)
91{
92 if (verbose) {
93 std::ostringstream os;
94 os << *prefix << "Allocate & initialize Kokkos::View "
95 << name << ": " << size << std::endl;
96 std::cerr << os.str();
97 }
98 return ViewType(name, size);
99}
100
101template<class OutViewType, class InViewType>
102void
103assign_to_view(OutViewType& out,
104 const InViewType& in,
105 const char viewName[],
106 const bool verbose,
107 const std::string* const prefix)
108{
109 if (verbose) {
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();
115 }
116 out = in;
117}
118
119template<class MemorySpace, class ViewType>
120auto create_mirror_view(
121 const MemorySpace& memSpace,
122 const ViewType& view,
123 const bool verbose,
124 const std::string* const prefix) ->
125 decltype(Kokkos::create_mirror_view(memSpace, view))
126{
127 if (verbose) {
128 std::ostringstream os;
129 os << *prefix << "Create mirror view: "
130 << "view.extent(0): " << view.extent(0) << std::endl;
131 std::cerr << os.str();
132 }
133 return Kokkos::create_mirror_view(memSpace, view);
134}
135
136enum class PadCrsAction {
137 INDICES_ONLY,
138 INDICES_AND_VALUES
139};
140
149template<class RowPtr, class Indices, class Values, class Padding>
150void
152 const PadCrsAction action,
153 const RowPtr& row_ptr_beg,
154 const RowPtr& row_ptr_end,
155 Indices& indices_wdv,
156 Values& values_wdv,
157 const Padding& padding,
158 const int my_rank,
159 const bool verbose)
160{
161 using execution_space = typename Indices::t_dev::execution_space;
162 using Kokkos::view_alloc;
163 using Kokkos::WithoutInitializing;
164 using std::endl;
165 std::unique_ptr<std::string> prefix;
166
167 const size_t maxNumToPrint = verbose ?
169 if (verbose) {
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()));
173 os << "Start" << endl;
174 std::cerr << os.str();
175 }
176 Kokkos::HostSpace hostSpace;
177
178 if (verbose) {
179 std::ostringstream os;
180 os << *prefix << "On input: ";
181 auto row_ptr_beg_h =
182 Kokkos::create_mirror_view(hostSpace, row_ptr_beg);
183 // DEEP_COPY REVIEW - NOT TESTED
184 Kokkos::deep_copy(row_ptr_beg_h, row_ptr_beg);
185 verbosePrintArray(os, row_ptr_beg_h, "row_ptr_beg before scan",
186 maxNumToPrint);
187 os << ", ";
188 auto row_ptr_end_h =
189 Kokkos::create_mirror_view(hostSpace, row_ptr_end);
190 // DEEP_COPY REVIEW - NOT TESTED
191 Kokkos::deep_copy(row_ptr_end_h, row_ptr_end);
192 verbosePrintArray(os, row_ptr_end_h, "row_ptr_end before scan",
193 maxNumToPrint);
194 os << ", indices.extent(0): " << indices_wdv.extent(0)
195 << ", values.extent(0): " << values_wdv.extent(0)
196 << ", padding: ";
197 padding.print(os);
198 os << endl;
199 std::cerr << os.str();
200 }
201
202 if (row_ptr_beg.size() == 0) {
203 if (verbose) {
204 std::ostringstream os;
205 os << *prefix << "Done; local matrix has no rows" << endl;
206 std::cerr << os.str();
207 }
208 return; // nothing to do
209 }
210
211 const size_t lclNumRows(row_ptr_beg.size() - 1);
212 RowPtr newAllocPerRow =
213 make_uninitialized_view<RowPtr>("newAllocPerRow", lclNumRows,
214 verbose, prefix.get());
215 if (verbose) {
216 std::ostringstream os;
217 os << *prefix << "Fill newAllocPerRow & compute increase" << endl;
218 std::cerr << os.str();
219 }
220 size_t increase = 0;
221 {
222 // Must do on host because padding uses std::map
223 auto row_ptr_end_h = create_mirror_view(
224 hostSpace, row_ptr_end, verbose, prefix.get());
225 // DEEP_COPY REVIEW - DEVICE-TO-HOSTMIRROR
226 Kokkos::deep_copy(execution_space(), row_ptr_end_h, row_ptr_end);
227 auto row_ptr_beg_h = create_mirror_view(
228 hostSpace, row_ptr_beg, verbose, prefix.get());
229 // DEEP_COPY REVIEW - DEVICE-TO-HOSTMIRROR
230 Kokkos::deep_copy(execution_space(), row_ptr_beg_h, row_ptr_beg);
231
232 auto newAllocPerRow_h = create_mirror_view(
233 hostSpace, newAllocPerRow, verbose, prefix.get());
234 using host_range_type = Kokkos::RangePolicy<
235 Kokkos::DefaultHostExecutionSpace, size_t>;
236 Kokkos::parallel_reduce
237 ("Tpetra::CrsGraph: Compute new allocation size per row",
238 host_range_type(0, lclNumRows),
239 [&] (const size_t lclRowInd, size_t& lclIncrease) {
240 const size_t start = row_ptr_beg_h[lclRowInd];
241 const size_t end = row_ptr_beg_h[lclRowInd+1];
242 TEUCHOS_ASSERT( end >= start );
243 const size_t oldAllocSize = end - start;
244 const size_t oldNumEnt = row_ptr_end_h[lclRowInd] - start;
245 TEUCHOS_ASSERT( oldNumEnt <= oldAllocSize );
246
247 // This is not a pack routine. Do not shrink! Shrinking now
248 // to fit the number of entries would ignore users' hint for
249 // the max number of entries in each row. Also, CrsPadding
250 // only counts entries and ignores any available free space.
251
252 auto result = padding.get_result(lclRowInd);
253 const size_t newNumEnt = oldNumEnt + result.numInSrcNotInTgt;
254 if (newNumEnt > oldAllocSize) {
255 lclIncrease += (newNumEnt - oldAllocSize);
256 newAllocPerRow_h[lclRowInd] = newNumEnt;
257 }
258 else {
259 newAllocPerRow_h[lclRowInd] = oldAllocSize;
260 }
261 }, increase);
262
263 if (verbose) {
264 std::ostringstream os;
265 os << *prefix << "increase: " << increase << ", ";
266 verbosePrintArray(os, newAllocPerRow_h, "newAllocPerRow",
267 maxNumToPrint);
268 os << endl;
269 std::cerr << os.str();
270 }
271
272 if (increase == 0) {
273 return;
274 }
275 // DEEP_COPY REVIEW - HOSTMIRROR-TO-DEVICE
276 Kokkos::deep_copy(execution_space(), newAllocPerRow, newAllocPerRow_h);
277 }
278
279 using inds_value_type =
280 typename Indices::t_dev::non_const_value_type;
281 using vals_value_type = typename Values::t_dev::non_const_value_type;
282
283 {
284 auto indices_old = indices_wdv.getDeviceView(Access::ReadOnly);
285 const size_t newIndsSize = size_t(indices_old.size()) + increase;
286 auto indices_new = make_uninitialized_view<typename Indices::t_dev>(
287 "Tpetra::CrsGraph column indices", newIndsSize, verbose,
288 prefix.get());
289
290 typename Values::t_dev values_new;
291 auto values_old = values_wdv.getDeviceView(Access::ReadOnly);
292 if (action == PadCrsAction::INDICES_AND_VALUES) {
293 const size_t newValsSize = newIndsSize;
294 // NOTE (mfh 10 Feb 2020) If we don't initialize values_new here,
295 // then the CrsMatrix tests fail.
296 values_new = make_initialized_view<typename Values::t_dev>(
297 "Tpetra::CrsMatrix values", newValsSize, verbose, prefix.get());
298 }
299
300 if (verbose) {
301 std::ostringstream os;
302 os << *prefix << "Repack" << endl;
303 std::cerr << os.str();
304 }
305
306 using range_type = Kokkos::RangePolicy<execution_space, size_t>;
307 Kokkos::parallel_scan(
308 "Tpetra::CrsGraph or CrsMatrix repack",
309 range_type(size_t(0), size_t(lclNumRows+1)),
310 KOKKOS_LAMBDA (const size_t lclRow, size_t& newRowBeg,
311 const bool finalPass)
312 {
313 // row_ptr_beg has lclNumRows + 1 entries.
314 // row_ptr_end has lclNumRows entries.
315 // newAllocPerRow has lclNumRows entries.
316 const size_t row_beg = row_ptr_beg[lclRow];
317 const size_t row_end =
318 lclRow < lclNumRows ? row_ptr_end[lclRow] : row_beg;
319 const size_t numEnt = row_end - row_beg;
320 const size_t newRowAllocSize =
321 lclRow < lclNumRows ? newAllocPerRow[lclRow] : size_t(0);
322 if (finalPass) {
323 if (lclRow < lclNumRows) {
324 const Kokkos::pair<size_t, size_t> oldRange(
325 row_beg, row_beg + numEnt);
326 const Kokkos::pair<size_t, size_t> newRange(
327 newRowBeg, newRowBeg + numEnt);
328 auto oldColInds = Kokkos::subview(indices_old, oldRange);
329 auto newColInds = Kokkos::subview(indices_new, newRange);
330 // memcpy works fine on device; the next step is to
331 // introduce two-level parallelism and use team copy.
332 memcpy(newColInds.data(), oldColInds.data(),
333 numEnt * sizeof(inds_value_type));
334 if (action == PadCrsAction::INDICES_AND_VALUES) {
335 auto oldVals =
336 Kokkos::subview(values_old, oldRange);
337 auto newVals = Kokkos::subview(values_new, newRange);
338 memcpy(newVals.data(), oldVals.data(),
339 numEnt * sizeof(vals_value_type));
340 }
341 }
342 // It's the final pass, so we can modify these arrays.
343 row_ptr_beg[lclRow] = newRowBeg;
344 if (lclRow < lclNumRows) {
345 row_ptr_end[lclRow] = newRowBeg + numEnt;
346 }
347 }
348 newRowBeg += newRowAllocSize;
349 });
350
351 if (verbose)
352 {
353 std::ostringstream os;
354
355 os << *prefix;
356 auto row_ptr_beg_h =
357 Kokkos::create_mirror_view(hostSpace, row_ptr_beg);
358 // DEEP_COPY REVIEW - NOT TESTED
359 Kokkos::deep_copy(row_ptr_beg_h, row_ptr_beg);
360 verbosePrintArray(os, row_ptr_beg_h, "row_ptr_beg after scan",
361 maxNumToPrint);
362 os << endl;
363
364 os << *prefix;
365 auto row_ptr_end_h =
366 Kokkos::create_mirror_view(hostSpace, row_ptr_end);
367 // DEEP_COPY REVIEW - NOT TESTED
368 Kokkos::deep_copy(row_ptr_end_h, row_ptr_end);
369 verbosePrintArray(os, row_ptr_end_h, "row_ptr_end after scan",
370 maxNumToPrint);
371 os << endl;
372
373 std::cout << os.str();
374 }
375
376 indices_wdv = Indices(indices_new);
377 values_wdv = Values(values_new);
378 }
379
380 if (verbose) {
381 auto indices_h = indices_wdv.getHostView(Access::ReadOnly);
382 auto values_h = values_wdv.getHostView(Access::ReadOnly);
383 std::ostringstream os;
384 os << "On output: ";
385 verbosePrintArray(os, indices_h, "indices", maxNumToPrint);
386 os << ", ";
387 verbosePrintArray(os, values_h, "values", maxNumToPrint);
388 os << ", padding: ";
389 padding.print(os);
390 os << endl;
391 }
392
393 if (verbose) {
394 std::ostringstream os;
395 os << *prefix << "Done" << endl;
396 std::cerr << os.str();
397 }
398}
399
401template <class Pointers, class InOutIndices, class InIndices, class IndexMap>
402size_t
404 typename Pointers::value_type const row,
405 Pointers const& row_ptrs,
406 InOutIndices& cur_indices,
407 size_t& num_assigned,
408 InIndices const& new_indices,
409 IndexMap&& map,
410 std::function<void(size_t const, size_t const, size_t const)> cb)
411{
412 if (new_indices.size() == 0) {
413 return 0;
414 }
415
416 if (cur_indices.size() == 0) {
417 // No room to insert new indices
418 return Teuchos::OrdinalTraits<size_t>::invalid();
419 }
420
421 using offset_type = typename std::decay<decltype (row_ptrs[0])>::type;
422 using ordinal_type = typename std::decay<decltype (cur_indices[0])>::type;
423
424 const offset_type start = row_ptrs[row];
425 offset_type end = start + static_cast<offset_type> (num_assigned);
426 const size_t num_avail = (row_ptrs[row + 1] < end) ? size_t (0) :
427 row_ptrs[row + 1] - end;
428 const size_t num_new_indices = static_cast<size_t> (new_indices.size ());
429 size_t num_inserted = 0;
430
431 size_t numIndicesLookup = num_assigned + num_new_indices;
432
433 // Threshold determined from test/Utils/insertCrsIndicesThreshold.cpp
434 const size_t useLookUpTableThreshold = 400;
435
436 if (numIndicesLookup <= useLookUpTableThreshold || num_new_indices == 1) {
437 // For rows with few nonzeros, can use a serial search to find duplicates
438 // Or if inserting only one index, serial search is as fast as anything else
439 for (size_t k = 0; k < num_new_indices; ++k) {
440 const ordinal_type idx = std::forward<IndexMap>(map)(new_indices[k]);
441 offset_type row_offset = start;
442 for (; row_offset < end; ++row_offset) {
443 if (idx == cur_indices[row_offset]) {
444 break;
445 }
446 }
447
448 if (row_offset == end) {
449 if (num_inserted >= num_avail) { // not enough room
450 return Teuchos::OrdinalTraits<size_t>::invalid();
451 }
452 // This index is not yet in indices
453 cur_indices[end++] = idx;
454 num_inserted++;
455 }
456 if (cb) {
457 cb(k, start, row_offset - start);
458 }
459 }
460 }
461 else {
462 // For rows with many nonzeros, use a lookup table to find duplicates
463 std::unordered_map<ordinal_type, offset_type> idxLookup(numIndicesLookup);
464
465 // Put existing indices into the lookup table
466 for (size_t k = 0; k < num_assigned; k++) {
467 idxLookup[cur_indices[start+k]] = start+k;
468 }
469
470 // Check for new indices in table; insert if not there yet
471 for (size_t k = 0; k < num_new_indices; k++) {
472 const ordinal_type idx = std::forward<IndexMap>(map)(new_indices[k]);
473 offset_type row_offset;
474
475 auto it = idxLookup.find(idx);
476 if (it == idxLookup.end()) {
477 if (num_inserted >= num_avail) { // not enough room
478 return Teuchos::OrdinalTraits<size_t>::invalid();
479 }
480 // index not found; insert it
481 row_offset = end;
482 cur_indices[end++] = idx;
483 idxLookup[idx] = row_offset;
484 num_inserted++;
485 }
486 else {
487 // index found; note its position
488 row_offset = it->second;
489 }
490 if (cb) {
491 cb(k, start, row_offset - start);
492 }
493 }
494 }
495 num_assigned += num_inserted;
496 return num_inserted;
497}
498
500template <class Pointers, class Indices1, class Indices2, class IndexMap, class Callback>
501size_t
503 typename Pointers::value_type const row,
504 Pointers const& row_ptrs,
505 const size_t curNumEntries,
506 Indices1 const& cur_indices,
507 Indices2 const& new_indices,
508 IndexMap&& map,
509 Callback&& cb)
510{
511 if (new_indices.size() == 0)
512 return 0;
513
514 using ordinal =
515 typename std::remove_const<typename Indices1::value_type>::type;
516 auto invalid_ordinal = Teuchos::OrdinalTraits<ordinal>::invalid();
517
518 const size_t start = static_cast<size_t> (row_ptrs[row]);
519 const size_t end = start + curNumEntries;
520 size_t num_found = 0;
521 for (size_t k = 0; k < new_indices.size(); k++)
522 {
523 auto row_offset = start;
524 auto idx = std::forward<IndexMap>(map)(new_indices[k]);
525 if (idx == invalid_ordinal)
526 continue;
527 for (; row_offset < end; row_offset++)
528 {
529 if (idx == cur_indices[row_offset])
530 {
531 std::forward<Callback>(cb)(k, start, row_offset - start);
532 num_found++;
533 }
534 }
535 }
536 return num_found;
537}
538
539} // namespace impl
540
541
559template<class RowPtr, class Indices, class Padding>
560void
562 const RowPtr& rowPtrBeg,
563 const RowPtr& rowPtrEnd,
564 Indices& indices_wdv,
565 const Padding& padding,
566 const int my_rank,
567 const bool verbose)
568{
569 using impl::pad_crs_arrays;
570 // send empty values array
571 Indices values_null;
572 pad_crs_arrays<RowPtr, Indices, Indices, Padding>(
573 impl::PadCrsAction::INDICES_ONLY, rowPtrBeg, rowPtrEnd,
574 indices_wdv, values_null, padding, my_rank, verbose);
575}
576
577template<class RowPtr, class Indices, class Values, class Padding>
578void
580 const RowPtr& rowPtrBeg,
581 const RowPtr& rowPtrEnd,
582 Indices& indices_wdv,
583 Values& values_wdv,
584 const Padding& padding,
585 const int my_rank,
586 const bool verbose)
587{
588 using impl::pad_crs_arrays;
589 pad_crs_arrays<RowPtr, Indices, Values, Padding>(
590 impl::PadCrsAction::INDICES_AND_VALUES, rowPtrBeg, rowPtrEnd,
591 indices_wdv, values_wdv, padding, my_rank, verbose);
592}
593
639template <class Pointers, class InOutIndices, class InIndices>
640size_t
642 typename Pointers::value_type const row,
643 Pointers const& rowPtrs,
644 InOutIndices& curIndices,
645 size_t& numAssigned,
646 InIndices const& newIndices,
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)>())
649{
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");
653
654 // Provide a unit map for the more general insert_indices
655 using ordinal = typename InOutIndices::value_type;
656 auto numInserted = impl::insert_crs_indices(row, rowPtrs, curIndices,
657 numAssigned, newIndices, [](ordinal const idx) { return idx; }, cb);
658 return numInserted;
659}
660
661template <class Pointers, class InOutIndices, class InIndices>
662size_t
664 typename Pointers::value_type const row,
665 Pointers const& rowPtrs,
666 InOutIndices& curIndices,
667 size_t& numAssigned,
668 InIndices const& newIndices,
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)>())
672{
673 auto numInserted = impl::insert_crs_indices(row, rowPtrs, curIndices,
674 numAssigned, newIndices, map, cb);
675 return numInserted;
676}
677
678
708template <class Pointers, class Indices1, class Indices2, class Callback>
709size_t
711 typename Pointers::value_type const row,
712 Pointers const& rowPtrs,
713 const size_t curNumEntries,
714 Indices1 const& curIndices,
715 Indices2 const& newIndices,
716 Callback&& cb)
717{
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");
721 // Provide a unit map for the more general find_crs_indices
722 using ordinal = typename Indices2::value_type;
723 auto numFound = impl::find_crs_indices(row, rowPtrs, curNumEntries, curIndices, newIndices,
724 [=](ordinal ind){ return ind; }, cb);
725 return numFound;
726}
727
728template <class Pointers, class Indices1, class Indices2, class IndexMap, class Callback>
729size_t
731 typename Pointers::value_type const row,
732 Pointers const& rowPtrs,
733 const size_t curNumEntries,
734 Indices1 const& curIndices,
735 Indices2 const& newIndices,
736 IndexMap&& map,
737 Callback&& cb)
738{
739 return impl::find_crs_indices(row, rowPtrs, curNumEntries, curIndices, newIndices, map, cb);
740}
741
742} // namespace Details
743} // namespace Tpetra
744
745#endif // TPETRA_DETAILS_CRSUTILS_HPP
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.
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.