22 #ifdef DEAL_II_WITH_TRILINOS
23 # ifdef DEAL_II_WITH_MPI
24 # include <Epetra_MpiComm.h>
26 # include <Epetra_Map.h>
27 # include <Epetra_SerialComm.h>
34 #ifdef DEAL_II_WITH_TRILINOS
39 # ifdef DEAL_II_WITH_64BIT_INDICES
43 , index_space_size(1 + map.MaxAllGID64())
46 Assert(map.MinAllGID64() == 0,
48 "The Epetra_BlockMap does not contain the global index 0, "
49 "which means some entries are not present on any processor."));
56 const size_type n_indices = map.NumMyElements();
58 reinterpret_cast<size_type *
>(map.MyGlobalElements64());
59 add_indices(indices, indices + n_indices);
70 , index_space_size(1 + map.MaxAllGID())
73 Assert(map.MinAllGID() == 0,
75 "The Epetra_BlockMap does not contain the global index 0, "
76 "which means some entries are not present on any processor."));
83 const size_type n_indices = map.NumMyElements();
84 unsigned int * indices =
85 reinterpret_cast<unsigned int *
>(map.MyGlobalElements());
111 std::vector<Range>::iterator store =
ranges.begin();
112 for (std::vector<Range>::iterator i =
ranges.begin(); i !=
ranges.end();)
114 std::vector<Range>::iterator next = i;
121 while (next !=
ranges.end() && (next->begin <= last_index))
123 last_index =
std::max(last_index, next->end);
129 *store =
Range(first_index, last_index);
133 if (store !=
ranges.end())
135 std::vector<Range> new_ranges(
ranges.begin(), store);
140 size_type next_index = 0, largest_range_size = 0;
141 for (std::vector<Range>::iterator i =
ranges.begin(); i !=
ranges.end(); ++i)
145 i->nth_index_in_set = next_index;
146 next_index += (i->end - i->begin);
147 if (i->end - i->begin > largest_range_size)
149 largest_range_size = i->end - i->begin;
171 std::vector<Range>::const_iterator r1 =
ranges.begin(),
179 if (r1->end <= r2->begin)
181 else if (r2->end <= r1->begin)
186 Assert(((r1->begin <= r2->begin) && (r1->end > r2->begin)) ||
187 ((r2->begin <= r1->begin) && (r2->end > r1->begin)),
191 result.add_range(
std::max(r1->begin, r2->begin),
197 if (r1->end <= r2->end)
215 ExcMessage(
"End index needs to be larger or equal to begin index!"));
219 std::vector<Range>::const_iterator r1 =
ranges.begin();
221 while (r1 !=
ranges.end())
223 if ((r1->end >
begin) && (r1->begin <
end))
228 else if (r1->begin >=
end)
238 std::vector<IndexSet>
240 const std::vector<types::global_dof_index> &n_indices_per_block)
const
242 std::vector<IndexSet> partitioned;
243 const unsigned int n_blocks = n_indices_per_block.size();
247 for (
const auto n_block_indices : n_indices_per_block)
249 partitioned.push_back(this->
get_view(start, start + n_block_indices));
250 start += n_block_indices;
255 for (
const auto &
partition : partitioned)
275 std::vector<Range> new_ranges;
277 std::vector<Range>::iterator own_it =
ranges.begin();
278 std::vector<Range>::iterator other_it = other.
ranges.begin();
280 while (own_it !=
ranges.end() && other_it != other.
ranges.end())
283 if (own_it->end <= other_it->begin)
285 new_ranges.push_back(*own_it);
290 if (own_it->begin >= other_it->end)
298 if (own_it->begin < other_it->begin)
300 Range r(own_it->begin, other_it->begin);
302 new_ranges.push_back(r);
307 own_it->begin = other_it->end;
308 if (own_it->begin > own_it->end)
310 own_it->begin = own_it->end;
319 for (; own_it !=
ranges.end(); ++own_it)
320 new_ranges.push_back(*own_it);
325 const std::vector<Range>::iterator
end = new_ranges.end();
326 for (std::vector<Range>::iterator it = new_ranges.begin(); it !=
end; ++it)
338 for (
const auto el : *
this)
339 set.add_indices(other, el * other.
size());
351 "pop_back() failed, because this IndexSet contains no entries."));
369 "pop_front() failed, because this IndexSet contains no entries."));
397 boost::container::small_vector<std::pair<size_type, size_type>, 200>
399 const bool ranges_are_sorted)
401 if (!ranges_are_sorted)
402 std::sort(tmp_ranges.begin(), tmp_ranges.end());
412 if (tmp_ranges.size() > 9)
415 tmp_set.
ranges.reserve(tmp_ranges.size());
416 for (
const auto &i : tmp_ranges)
421 if (this->
ranges.size() <= 1)
423 if (this->
ranges.size() == 1)
431 for (
const auto &i : tmp_ranges)
440 if ((
this == &other) && (offset == 0))
443 if (other.
ranges.size() != 0)
451 std::vector<Range>::const_iterator r1 =
ranges.begin(),
452 r2 = other.
ranges.begin();
454 std::vector<Range> new_ranges;
461 if (r2 == other.
ranges.end() ||
462 (r1 !=
ranges.end() && r1->end < (r2->begin + offset)))
464 new_ranges.push_back(*r1);
467 else if (r1 ==
ranges.end() || (r2->end + offset) < r1->begin)
469 new_ranges.emplace_back(r2->begin + offset, r2->end + offset);
477 std::max(r1->end, r2->end + offset));
478 new_ranges.push_back(next);
495 out <<
size() <<
" ";
496 out <<
ranges.size() << std::endl;
497 std::vector<Range>::const_iterator r =
ranges.begin();
498 for (; r !=
ranges.end(); ++r)
500 out << r->begin <<
" " << r->end << std::endl;
512 unsigned int n_ranges;
517 for (
unsigned int i = 0; i < n_ranges; ++i)
534 std::size_t n_ranges =
ranges.size();
535 out.write(
reinterpret_cast<const char *
>(&n_ranges),
sizeof(n_ranges));
536 if (
ranges.empty() ==
false)
537 out.write(
reinterpret_cast<const char *
>(&*
ranges.begin()),
546 std::size_t n_ranges;
547 in.read(
reinterpret_cast<char *
>(&
size),
sizeof(
size));
548 in.read(
reinterpret_cast<char *
>(&n_ranges),
sizeof(n_ranges));
554 in.read(
reinterpret_cast<char *
>(&*
ranges.begin()),
577 std::vector<Range>::const_iterator p = std::upper_bound(
585 return ((index >= p->begin) && (index < p->
end));
593 return (p->end > index);
610 return p->begin + (n - p->nth_index_in_set);
622 std::vector<Range>::const_iterator p =
626 if (p ==
ranges.end() || p->end == n || p->begin > n)
632 return (n - p->begin) + p->nth_index_in_set;
646 std::vector<Range>::const_iterator main_range =
652 std::vector<Range>::const_iterator range_begin, range_end;
653 if (global_index < main_range->
begin)
655 range_begin =
ranges.begin();
656 range_end = main_range;
660 range_begin = main_range;
666 const std::vector<Range>::const_iterator p =
679 if (global_index < p->
begin)
695 for (
const auto &range :
ranges)
696 for (
size_type entry = range.begin; entry < range.end; ++entry)
697 indices.push_back(entry);
704 #ifdef DEAL_II_WITH_TRILINOS
705 # ifdef DEAL_II_TRILINOS_WITH_TPETRA
707 Tpetra::Map<int, types::global_dof_index>
709 const bool overlapping)
const
720 ExcMessage(
"You are trying to create an Tpetra::Map object "
721 "that partitions elements of an index set "
722 "between processors. However, the union of the "
723 "index sets on different processors does not "
724 "contain all indices exactly once: the sum of "
725 "the number of entries the various processors "
726 "want to store locally is " +
728 " whereas the total size of the object to be "
731 ". In other words, there are "
732 "either indices that are not spoken for "
733 "by any processor, or there are indices that are "
734 "claimed by multiple processors."));
743 return Tpetra::Map<int, types::global_dof_index>(
747 # ifdef DEAL_II_WITH_MPI
748 Teuchos::rcp(
new Teuchos::MpiComm<int>(communicator))
750 Teuchos::rcp(
new Teuchos::Comm<int>())
755 std::vector<size_type> indices;
757 std::vector<types::global_dof_index> int_indices(indices.size());
758 std::copy(indices.begin(), indices.end(), int_indices.begin());
759 const Teuchos::ArrayView<types::global_dof_index> arr_view(int_indices);
760 return Tpetra::Map<int, types::global_dof_index>(
764 # ifdef DEAL_II_WITH_MPI
765 Teuchos::rcp(
new Teuchos::MpiComm<int>(communicator))
767 Teuchos::rcp(
new Teuchos::Comm<int>())
778 const bool overlapping)
const
789 ExcMessage(
"You are trying to create an Epetra_Map object "
790 "that partitions elements of an index set "
791 "between processors. However, the union of the "
792 "index sets on different processors does not "
793 "contain all indices exactly once: the sum of "
794 "the number of entries the various processors "
795 "want to store locally is " +
797 " whereas the total size of the object to be "
800 ". In other words, there are "
801 "either indices that are not spoken for "
802 "by any processor, or there are indices that are "
803 "claimed by multiple processors."));
816 # ifdef DEAL_II_WITH_MPI
817 Epetra_MpiComm(communicator)
824 std::vector<size_type> indices;
834 # ifdef DEAL_II_WITH_MPI
835 Epetra_MpiComm(communicator)
859 #ifdef DEAL_II_WITH_MPI
861 const bool all_contiguous =
866 bool is_globally_ascending =
true;
873 const std::vector<types::global_dof_index> global_dofs =
883 for (; index < global_dofs.size(); ++index)
888 if (new_dof <= old_dof)
890 is_globally_ascending =
false;
900 int is_ascending = is_globally_ascending ? 1 : 0;
901 int ierr = MPI_Bcast(&is_ascending, 1, MPI_INT, 0, communicator);
904 return (is_ascending == 1);
ElementIterator begin() const
bool is_element_binary_search(const size_type local_index) const
size_type index_within_set_binary_search(const size_type global_index) const
bool is_contiguous() const
ElementIterator at(const size_type global_index) const
void fill_index_vector(std::vector< size_type > &indices) const
std::vector< IndexSet > split_by_block(const std::vector< types::global_dof_index > &n_indices_per_block) const
size_type n_elements() const
void add_range_lower_bound(const Range &range)
ElementIterator begin() const
void set_size(const size_type size)
void read(std::istream &in)
Tpetra::Map< int, types::global_dof_index > make_tpetra_map(const MPI_Comm &communicator=MPI_COMM_WORLD, const bool overlapping=false) const
IndexSet tensor_product(const IndexSet &other) const
void write(std::ostream &out) const
void block_read(std::istream &in)
void add_ranges_internal(boost::container::small_vector< std::pair< size_type, size_type >, 200 > &tmp_ranges, const bool ranges_are_sorted)
IntervalIterator begin_intervals() const
std::vector< Range > ranges
void subtract_set(const IndexSet &other)
ElementIterator end() const
Threads::Mutex compress_mutex
size_type index_space_size
Epetra_Map make_trilinos_map(const MPI_Comm &communicator=MPI_COMM_WORLD, const bool overlapping=false) const
void block_write(std::ostream &out) const
IndexSet get_view(const size_type begin, const size_type end) const
void add_range(const size_type begin, const size_type end)
bool is_ascending_and_one_to_one(const MPI_Comm &communicator) const
std::size_t memory_consumption() const
size_type nth_index_in_set_binary_search(const size_type local_index) const
types::global_dof_index size_type
void add_indices(const ForwardIterator &begin, const ForwardIterator &end)
IndexSet operator&(const IndexSet &is) const
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
__global__ void set(Number *val, const Number s, const size_type N)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define Assert(cond, exc)
std::string to_string(const T &t)
#define AssertDimension(dim1, dim2)
#define AssertThrowMPI(error_code)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcIO()
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
types::global_dof_index size_type
std::enable_if< IsBlockVector< VectorType >::value, unsigned int >::type n_blocks(const VectorType &vector)
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
void swap(MemorySpaceData< Number, MemorySpace > &, MemorySpaceData< Number, MemorySpace > &)
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
VectorType::value_type * begin(VectorType &V)
VectorType::value_type * end(VectorType &V)
TrilinosWrappers::types::int_type global_index(const Epetra_BlockMap &map, const ::types::global_dof_index i)
TrilinosWrappers::types::int64_type n_global_elements(const Epetra_BlockMap &map)
std::vector< T > gather(const MPI_Comm &comm, const T &object_to_send, const unsigned int root_process=0)
unsigned int this_mpi_process(const MPI_Comm &mpi_communicator)
T min(const T &t, const MPI_Comm &mpi_communicator)
T sum(const T &t, const MPI_Comm &mpi_communicator)
std::string compress(const std::string &input)
Iterator lower_bound(Iterator first, Iterator last, const T &val)
void copy(const T *begin, const T *end, U *dest)
static const unsigned int invalid_unsigned_int
const types::global_dof_index invalid_dof_index
unsigned int global_dof_index
static bool end_compare(const IndexSet::Range &x, const IndexSet::Range &y)
static bool nth_index_compare(const IndexSet::Range &x, const IndexSet::Range &y)
size_type nth_index_in_set