17#ifdef DEAL_II_TRILINOS_WITH_TPETRA
34 namespace TpetraWrappers
38 template <
typename MemorySpace>
44 if (
static_cast<std::size_t
>(this->a_row) == sparsity_pattern->n_rows())
51 if (!sparsity_pattern->is_compressed())
52 sparsity_pattern->compress();
55 std::make_shared<std::vector<::types::signed_global_dof_index>>(
56 sparsity_pattern->row_length(this->a_row));
58 if (colnum_cache->size() > 0)
62# if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
65 column_indices_view(colnum_cache->data(), colnum_cache->size());
67 Teuchos::ArrayView<::types::signed_global_dof_index>
68 column_indices_view(colnum_cache->data(), colnum_cache->size());
70 sparsity_pattern->graph->getGlobalRowCopy(this->a_row,
84 template <
typename MemorySpace>
96 graph->fillComplete();
101 template <
typename MemorySpace>
107 reinit(m, n, n_entries_per_row);
112 template <
typename MemorySpace>
116 const std::vector<size_type> &n_entries_per_row)
118 reinit(m, n, n_entries_per_row);
123 template <
typename MemorySpace>
127 , column_space_map(std::move(other.column_space_map))
128 , graph(std::move(other.graph))
129 , nonlocal_graph(std::move(other.nonlocal_graph))
135 template <
typename MemorySpace>
143 Utilities::Trilinos::tpetra_comm_self()))
145 TpetraTypes::GraphType<
MemorySpace>>(column_space_map,
149 (void)input_sparsity;
152 "Copy constructor only works for empty sparsity patterns."));
157 template <
typename MemorySpace>
159 const IndexSet ¶llel_partitioning,
163 reinit(parallel_partitioning,
164 parallel_partitioning,
171 template <
typename MemorySpace>
173 const IndexSet ¶llel_partitioning,
175 const std::vector<size_type> &n_entries_per_row)
177 reinit(parallel_partitioning,
178 parallel_partitioning,
185 template <
typename MemorySpace>
187 const IndexSet &row_parallel_partitioning,
188 const IndexSet &col_parallel_partitioning,
192 reinit(row_parallel_partitioning,
193 col_parallel_partitioning,
200 template <
typename MemorySpace>
202 const IndexSet &row_parallel_partitioning,
203 const IndexSet &col_parallel_partitioning,
205 const std::vector<size_type> &n_entries_per_row)
207 reinit(row_parallel_partitioning,
208 col_parallel_partitioning,
215 template <
typename MemorySpace>
217 const IndexSet &row_parallel_partitioning,
218 const IndexSet &col_parallel_partitioning,
223 reinit(row_parallel_partitioning,
224 col_parallel_partitioning,
227 n_max_entries_per_row);
232 template <
typename MemorySpace>
246 template <
typename MemorySpace>
251 const std::vector<size_type> &n_entries_per_row)
261 namespace SparsityPatternImpl
263 template <
typename MemorySpace>
266 template <
typename MemorySpace>
276 Assert(row_map->isOneToOne(),
277 ExcMessage(
"Row map must be 1-to-1, i.e., no overlap between "
278 "the maps of different processors."));
279 Assert(col_map->isOneToOne(),
280 ExcMessage(
"Column map must be 1-to-1, i.e., no overlap between "
281 "the maps of different processors."));
283 nonlocal_graph.reset();
285 column_space_map = col_map;
303 template <
typename MemorySpace>
313 Assert(row_map->isOneToOne(),
314 ExcMessage(
"Row map must be 1-to-1, i.e., no overlap between "
315 "the maps of different processors."));
316 Assert(col_map->isOneToOne(),
317 ExcMessage(
"Column map must be 1-to-1, i.e., no overlap between "
318 "the maps of different processors."));
321 nonlocal_graph.reset();
324 row_map->getGlobalNumElements());
326 column_space_map = col_map;
330 Kokkos::DualView<size_t *, typename MemorySpace::kokkos_space>
331 local_entries_per_row(
"local_entries_per_row",
332 row_map->getMaxGlobalIndex() -
333 row_map->getMinGlobalIndex());
335 auto local_entries_per_row_host =
336 local_entries_per_row
337 .template view<Kokkos::DefaultHostExecutionSpace>();
339 std::uint64_t total_size = 0;
340 for (
unsigned int i = 0; i < local_entries_per_row.extent(0); ++i)
342 local_entries_per_row_host(i) =
343 n_entries_per_row[row_map->getMinGlobalIndex() + i];
344 total_size += local_entries_per_row_host[i];
346 local_entries_per_row
347 .template modify<Kokkos::DefaultHostExecutionSpace>();
348 local_entries_per_row
349 .template sync<typename MemorySpace::kokkos_space>();
352 total_size <
static_cast<std::uint64_t
>(
356 "You are requesting to store more elements than global ordinal type allows."));
361 local_entries_per_row);
366 template <
typename SparsityPatternType,
typename MemorySpace>
371 const SparsityPatternType &sp,
372 [[maybe_unused]]
const bool exchange_data,
377 nonlocal_graph.reset();
386 Assert(row_map->isContiguous() ==
true,
388 "This function only works if the row map is contiguous."));
392 row_map->getMaxGlobalIndex() + 1;
393 Teuchos::Array<size_t> n_entries_per_row(last_row - first_row);
396 n_entries_per_row[row - first_row] = sp.row_length(row);
399 std::accumulate(n_entries_per_row.begin(),
400 n_entries_per_row.end(),
402 static_cast<std::uint64_t
>(std::numeric_limits<int>::max()),
404 "The TrilinosWrappers use Tpetra internally, and "
405 "Trilinos/Tpetra was compiled with 'local ordinate = int'. "
406 "Therefore, 'signed int' is used to represent local indices, "
407 "and only 2,147,483,647 nonzero matrix entries can be stored "
408 "on a single process, but you are requesting more than "
409 "that. Either use more MPI processes or recompile Trilinos "
410 "with 'local ordinate = long long' "));
412# if DEAL_II_TRILINOS_VERSION_GTE(12, 16, 0)
413 if (row_map->getComm()->getSize() > 1)
422 if (row_map->getComm()->getSize() > 1)
425 row_map, Teuchos::arcpFromArray(n_entries_per_row));
429 row_map, col_map, Teuchos::arcpFromArray(n_entries_per_row));
436 std::vector<TrilinosWrappers::types::int_type> row_indices;
445 row_indices.resize(row_length, -1);
447 typename SparsityPatternType::iterator p = sp.begin(row);
450 for (
int col = 0; col < row_length;)
452 row_indices[col++] = p->column();
453 if (col < row_length)
457 graph->insertGlobalIndices(row, row_length, row_indices.data());
460 graph->globalAssemble();
465 template <
typename MemorySpace>
472 parallel_partitioning.
size());
473 Teuchos::RCP<TpetraTypes::MapType<MemorySpace>> map =
474 parallel_partitioning
475 .template make_tpetra_map_rcp<TpetraTypes::NodeType<MemorySpace>>(
476 communicator,
false);
477 SparsityPatternImpl::reinit_sp<MemorySpace>(
478 map, map, n_entries_per_row, column_space_map, graph, nonlocal_graph);
483 template <
typename MemorySpace>
486 const IndexSet ¶llel_partitioning,
488 const std::vector<size_type> &n_entries_per_row)
491 parallel_partitioning.
size());
492 Teuchos::RCP<TpetraTypes::MapType<MemorySpace>> map =
493 parallel_partitioning
494 .template make_tpetra_map_rcp<TpetraTypes::NodeType<MemorySpace>>(
495 communicator,
false);
496 SparsityPatternImpl::reinit_sp<MemorySpace>(
497 map, map, n_entries_per_row, column_space_map, graph, nonlocal_graph);
502 template <
typename MemorySpace>
505 const IndexSet &row_parallel_partitioning,
506 const IndexSet &col_parallel_partitioning,
511 col_parallel_partitioning.
size());
512 Teuchos::RCP<TpetraTypes::MapType<MemorySpace>> row_map =
513 row_parallel_partitioning
514 .template make_tpetra_map_rcp<TpetraTypes::NodeType<MemorySpace>>(
515 communicator,
false);
516 Teuchos::RCP<TpetraTypes::MapType<MemorySpace>> col_map =
517 col_parallel_partitioning
518 .template make_tpetra_map_rcp<TpetraTypes::NodeType<MemorySpace>>(
519 communicator,
false);
520 SparsityPatternImpl::reinit_sp<MemorySpace>(row_map,
530 template <
typename MemorySpace>
533 const IndexSet &row_parallel_partitioning,
534 const IndexSet &col_parallel_partitioning,
536 const std::vector<size_type> &n_entries_per_row)
539 col_parallel_partitioning.
size());
540 Teuchos::RCP<TpetraTypes::MapType<MemorySpace>> row_map =
541 row_parallel_partitioning
542 .template make_tpetra_map_rcp<TpetraTypes::NodeType<MemorySpace>>(
543 communicator,
false);
544 Teuchos::RCP<TpetraTypes::MapType<MemorySpace>> col_map =
545 col_parallel_partitioning
546 .template make_tpetra_map_rcp<TpetraTypes::NodeType<MemorySpace>>(
547 communicator,
false);
548 SparsityPatternImpl::reinit_sp<MemorySpace>(row_map,
558 template <
typename MemorySpace>
561 const IndexSet &row_parallel_partitioning,
562 const IndexSet &col_parallel_partitioning,
568 col_parallel_partitioning.
size());
569 Teuchos::RCP<TpetraTypes::MapType<MemorySpace>> row_map =
570 row_parallel_partitioning
571 .template make_tpetra_map_rcp<TpetraTypes::NodeType<MemorySpace>>(
572 communicator,
false);
573 Teuchos::RCP<TpetraTypes::MapType<MemorySpace>> col_map =
574 col_parallel_partitioning
575 .template make_tpetra_map_rcp<TpetraTypes::NodeType<MemorySpace>>(
576 communicator,
false);
577 SparsityPatternImpl::reinit_sp<MemorySpace>(row_map,
584 IndexSet nonlocal_partitioner = writable_rows;
586 row_parallel_partitioning.
size());
590 IndexSet tmp = writable_rows & row_parallel_partitioning;
591 Assert(tmp == row_parallel_partitioning,
593 "The set of writable rows passed to this method does not "
594 "contain the locally owned rows, which is not allowed."));
597 nonlocal_partitioner.subtract_set(row_parallel_partitioning);
600 Teuchos::RCP<TpetraTypes::MapType<MemorySpace>> nonlocal_map =
602 .template make_tpetra_map_rcp<TpetraTypes::NodeType<MemorySpace>>(
605 TpetraTypes::GraphType<MemorySpace>>(nonlocal_map, col_map, 0);
613 template <
typename MemorySpace>
614 template <
typename SparsityPatternType>
618 const IndexSet &col_parallel_partitioning,
619 const SparsityPatternType &nontrilinos_sparsity_pattern,
621 const bool exchange_data)
624 col_parallel_partitioning.
size());
625 Teuchos::RCP<TpetraTypes::MapType<MemorySpace>> row_map =
626 row_parallel_partitioning
627 .template make_tpetra_map_rcp<TpetraTypes::NodeType<MemorySpace>>(
628 communicator,
false);
629 Teuchos::RCP<TpetraTypes::MapType<MemorySpace>> col_map =
630 col_parallel_partitioning
631 .template make_tpetra_map_rcp<TpetraTypes::NodeType<MemorySpace>>(
632 communicator,
false);
633 SparsityPatternImpl::reinit_sp<SparsityPatternType, MemorySpace>(
636 nontrilinos_sparsity_pattern,
645 template <
typename MemorySpace>
646 template <
typename SparsityPatternType>
649 const IndexSet ¶llel_partitioning,
650 const SparsityPatternType &nontrilinos_sparsity_pattern,
652 const bool exchange_data)
655 parallel_partitioning.
size());
657 parallel_partitioning.
size());
659 parallel_partitioning.
size());
660 Teuchos::RCP<TpetraTypes::MapType<MemorySpace>> map =
661 parallel_partitioning
662 .template make_tpetra_map_rcp<TpetraTypes::NodeType<MemorySpace>>(
663 communicator,
false);
664 SparsityPatternImpl::reinit_sp<SparsityPatternType, MemorySpace>(
667 nontrilinos_sparsity_pattern,
660 Teuchos::RCP<TpetraTypes::MapType<MemorySpace>> map = {
…}
676 template <
typename MemorySpace>
687 template <
typename MemorySpace>
702 nonlocal_graph.reset();
707 template <
typename MemorySpace>
708 template <
typename SparsityPatternType>
713 Teuchos::RCP<TpetraTypes::MapType<MemorySpace>> rows =
719 Teuchos::RCP<TpetraTypes::MapType<MemorySpace>> columns =
726 SparsityPatternImpl::reinit_sp<SparsityPatternType, MemorySpace>(
727 rows, columns, sp,
false, column_space_map, graph, nonlocal_graph);
732 template <
typename MemorySpace>
749 graph->fillComplete();
751 nonlocal_graph.reset();
756 template <
typename MemorySpace>
761 if (nonlocal_graph.get() !=
nullptr)
763# if DEAL_II_TRILINOS_VERSION_GTE(14, 0, 0)
764 if (nonlocal_graph->getRowMap()->getLocalNumElements() > 0 &&
765 column_space_map->getGlobalNumElements() > 0)
767 if (nonlocal_graph->getRowMap()->getNodeNumElements() > 0 &&
768 column_space_map->getGlobalNumElements() > 0)
774 nonlocal_graph->getRowMap()->getGlobalElement(0);
779 if (column_space_map->getGlobalNumElements() ==
780 graph->getRangeMap()->getGlobalNumElements())
786# if DEAL_II_TRILINOS_VERSION_GTE(14, 0, 0)
787 else if (column_space_map->getLocalNumElements() > 0)
789 else if (column_space_map->getNodeNumElements() > 0)
791 column = column_space_map->getGlobalElement(0);
792 nonlocal_graph->insertGlobalIndices(row, 1, &column);
794# if DEAL_II_TRILINOS_VERSION_GTE(14, 0, 0)
795 Assert(nonlocal_graph->getRowMap()->getLocalNumElements() == 0 ||
796 column_space_map->getGlobalNumElements() == 0,
799 Assert(nonlocal_graph->getRowMap()->getNodeNumElements() == 0 ||
800 column_space_map->getGlobalNumElements() == 0,
804 nonlocal_graph->fillComplete(column_space_map, graph->getRowMap());
806 graph->fillComplete(column_space_map, graph->getRowMap());
817 template <
typename MemorySpace>
821 return graph->getRowMap()->getLocalElement(i) !=
822 Teuchos::OrdinalTraits<int>::invalid();
827 template <
typename MemorySpace>
832# if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
833 if (!row_is_stored_locally(i))
837 const auto trilinos_i = graph->getRowMap()->getLocalElement(i);
838 const auto trilinos_j = graph->getColMap()->getLocalElement(j);
844 graph->getLocalRowView(trilinos_i, col_indices);
848 std::find(col_indices.data(),
849 col_indices.data() + col_indices.size(),
853 return static_cast<std::size_t
>(local_col_index) != col_indices.size();
864 template <
typename MemorySpace>
869 for (
int i = 0; i < static_cast<int>(local_size()); ++i)
871# if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
875 Teuchos::ArrayView<const int> indices;
878 graph->getLocalRowView(i, indices);
879 const auto num_entries = indices.size();
880 for (
unsigned int j = 0; j < static_cast<unsigned int>(num_entries);
897 template <
typename MemorySpace>
901# if DEAL_II_TRILINOS_VERSION_GTE(14, 0, 0)
902 return graph->getLocalNumRows();
904 return graph->getNodeNumRows();
910 template <
typename MemorySpace>
911 std::pair<typename SparsityPattern<MemorySpace>::size_type,
915 const size_type begin = graph->getRowMap()->getMinGlobalIndex();
916 const size_type end = graph->getRowMap()->getMaxGlobalIndex() + 1;
923 template <
typename MemorySpace>
927 return graph->getGlobalNumEntries();
932 template <
typename MemorySpace>
936# if DEAL_II_TRILINOS_VERSION_GTE(14, 0, 0)
937 return graph->getLocalMaxNumRowEntries();
939 return graph->getNodeMaxNumRowEntries();
945 template <
typename MemorySpace>
954 graph->getRowMap()->getLocalElement(row);
959 return graph->getNumEntriesInLocalRow(local_row);
966 template <
typename MemorySpace>
969 const ::types::global_dof_index &row,
971 const bool indices_are_sorted)
973 add_entries(row, columns.
begin(), columns.
end(), indices_are_sorted);
978 template <
typename MemorySpace>
979 Teuchos::RCP<const typename TpetraTypes::MapType<MemorySpace>>
982 return graph->getDomainMap();
987 template <
typename MemorySpace>
988 Teuchos::RCP<const typename TpetraTypes::MapType<MemorySpace>>
991 return graph->getRangeMap();
996 template <
typename MemorySpace>
1001 graph->getRangeMap()->getComm());
1006 template <
typename MemorySpace>
1007 Teuchos::RCP<const Teuchos::Comm<int>>
1010 return graph->getRangeMap()->getComm();
1018 template <
typename MemorySpace>
1022 const bool write_extended_trilinos_info)
const
1024 if (write_extended_trilinos_info)
1028# if DEAL_II_TRILINOS_VERSION_GTE(14, 0, 0)
1029 for (
unsigned int i = 0; i < graph->getLocalNumRows(); ++i)
1031 for (
unsigned int i = 0; i < graph->getNodeNumRows(); ++i)
1034# if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
1038 Teuchos::ArrayView<const int> indices;
1040 graph->getLocalRowView(i, indices);
1041 int num_entries = indices.size();
1042 for (
int j = 0; j < num_entries; ++j)
1043 out <<
"(" << graph->getRowMap()->getGlobalElement(i) <<
","
1044 << graph->getColMap()->getGlobalElement(indices[j]) <<
") "
1054 template <
typename MemorySpace>
1060 for (
unsigned int row = 0; row < local_size(); ++row)
1062# if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0)
1066 Teuchos::ArrayView<const int> indices;
1068 graph->getLocalRowView(row, indices);
1069 int num_entries = indices.size();
1073 const ::types::signed_global_dof_index num_entries_ =
1082 out <<
static_cast<int>(
1083 graph->getColMap()->getGlobalElement(indices[j]))
1085 << -
static_cast<int>(graph->getRowMap()->getGlobalElement(row))
936# if DEAL_II_TRILINOS_VERSION_GTE(14, 0, 0) {
…}
1093 template <
typename MemorySpace>
1108 const ::SparsityPattern &);
1111 const ::DynamicSparsityPattern &);
1116 const ::SparsityPattern &,
1122 const ::DynamicSparsityPattern &,
1131 const ::SparsityPattern &,
1138 const ::DynamicSparsityPattern &,
1147 const ::SparsityPattern &);
1150 const ::DynamicSparsityPattern &);
1155 const ::SparsityPattern &,
1161 const ::DynamicSparsityPattern &,
1170 const ::SparsityPattern &,
1177 const ::DynamicSparsityPattern &,
918 return {begin, end}; {
…}
838 const auto trilinos_j = graph->getColMap()->getLocalElement(j); {
…}
832# if DEAL_II_TRILINOS_VERSION_GTE(13, 2, 0) {
…}
817 template <
typename MemorySpace> {
…}
765 column_space_map->getGlobalNumElements() > 0) {
…}
751 nonlocal_graph.reset(); {
…}
628 communicator,
false); {
…}
586 row_parallel_partitioning.
size()); {
…}
550 n_entries_per_row, {
…}
497 map, map, n_entries_per_row, column_space_map, graph, nonlocal_graph); {
…}
408 "on a single process, but you are requesting more than " {
…}
399 std::accumulate(n_entries_per_row.begin(), {
…}
348 local_entries_per_row {
…}
342 local_entries_per_row_host(i) = {
…}
261 namespace SparsityPatternImpl {
…}
Teuchos::RCP< TpetraTypes::GraphType< MemorySpace > > graph
Teuchos::RCP< const Teuchos::Comm< int > > get_teuchos_mpi_communicator() const
std::pair< size_type, size_type > local_range() const
Teuchos::RCP< const TpetraTypes::MapType< MemorySpace > > range_partitioner() const
void copy_from(const SparsityPattern< MemorySpace > &input_sparsity_pattern)
SparsityPattern< MemorySpace > & operator=(const SparsityPattern< MemorySpace > &input_sparsity_pattern)
void print(std::ostream &out, const bool write_extended_trilinos_info=false) const
unsigned int max_entries_per_row() const
Teuchos::RCP< TpetraTypes::GraphType< MemorySpace > > nonlocal_graph
void print_gnuplot(std::ostream &out) const
MPI_Comm get_mpi_communicator() const
Teuchos::RCP< const TpetraTypes::MapType< MemorySpace > > domain_partitioner() const
std::uint64_t n_nonzero_elements() const
size_type bandwidth() const
bool exists(const size_type i, const size_type j) const
Teuchos::RCP< TpetraTypes::MapType< MemorySpace > > column_space_map
unsigned int local_size() const
virtual void add_row_entries(const ::types::global_dof_index &row, const ArrayView< const ::types::global_dof_index > &columns, const bool indices_are_sorted=false) override
size_type row_length(const size_type row) const
std::size_t memory_consumption() const
bool row_is_stored_locally(const size_type i) const
void reinit(const size_type m, const size_type n, const size_type n_entries_per_row=0)
virtual void resize(const size_type rows, const size_type cols)
#define DEAL_II_NAMESPACE_OPEN
constexpr bool running_in_debug_mode()
#define DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NOT_IMPLEMENTED()
static ::ExceptionBase & ExcIO()
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
IndexSet complete_index_set(const IndexSet::size_type N)
void reinit_sp(const Teuchos::RCP< TpetraTypes::MapType< MemorySpace > > &row_map, const Teuchos::RCP< TpetraTypes::MapType< MemorySpace > > &col_map, const size_type< MemorySpace > n_entries_per_row, Teuchos::RCP< TpetraTypes::MapType< MemorySpace > > &column_space_map, Teuchos::RCP< TpetraTypes::GraphType< MemorySpace > > &graph, Teuchos::RCP< TpetraTypes::GraphType< MemorySpace > > &nonlocal_graph)
typename SparsityPattern< MemorySpace >::size_type size_type
Tpetra::CrsGraph< LO, GO, NodeType< MemorySpace > > GraphType
Tpetra::Map< LO, GO, NodeType< MemorySpace > > MapType
unsigned int n_mpi_processes(const MPI_Comm mpi_communicator)
T max(const T &t, const MPI_Comm mpi_communicator)
Teuchos::RCP< T > make_rcp(Args &&...args)
MPI_Comm teuchos_comm_to_mpi_comm(const Teuchos::RCP< const Teuchos::Comm< int > > &teuchos_comm)
const Teuchos::RCP< const Teuchos::Comm< int > > & tpetra_comm_self()
::VectorizedArray< Number, width > abs(const ::VectorizedArray< Number, width > &)