Reference documentation for deal.II version 9.4.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
trilinos_vector.cc
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2008 - 2022 by the deal.II authors
4//
5// This file is part of the deal.II library.
6//
7// The deal.II library is free software; you can use it, redistribute
8// it, and/or modify it under the terms of the GNU Lesser General
9// Public License as published by the Free Software Foundation; either
10// version 2.1 of the License, or (at your option) any later version.
11// The full text of the license can be found in the file LICENSE.md at
12// the top level directory of deal.II.
13//
14// ---------------------------------------------------------------------
15
17
18#ifdef DEAL_II_WITH_TRILINOS
19
20# include <deal.II/base/mpi.h>
21
26
27# include <boost/io/ios_state.hpp>
28
29# include <Epetra_Export.h>
30# include <Epetra_Import.h>
31# include <Epetra_Vector.h>
32
33# include <cmath>
34# include <memory>
35
36
38
39namespace TrilinosWrappers
40{
41# ifndef DOXYGEN
42 namespace internal
43 {
44 VectorReference::operator TrilinosScalar() const
45 {
46 AssertIndexRange(index, vector.size());
47
48 // Trilinos allows for vectors to be referenced by the [] or ()
49 // operators but only () checks index bounds. We check these bounds by
50 // ourselves, so we can use []. Note that we can only get local values.
51
52 const TrilinosWrappers::types::int_type local_index =
53 vector.vector->Map().LID(
54 static_cast<TrilinosWrappers::types::int_type>(index));
55 Assert(local_index >= 0,
57 index,
58 vector.local_size(),
59 vector.vector->Map().MinMyGID(),
60 vector.vector->Map().MaxMyGID()));
61
62
63 return (*(vector.vector))[0][local_index];
64 }
65 } // namespace internal
66# endif
67
68 namespace MPI
69 {
71 : Subscriptor()
72 , last_action(Zero)
73 , compressed(true)
74 , has_ghosts(false)
75 , vector(new Epetra_FEVector(
76 Epetra_Map(0, 0, 0, Utilities::Trilinos::comm_self())))
77 {}
78
79
80
81 Vector::Vector(const IndexSet &parallel_partitioning,
82 const MPI_Comm &communicator)
83 : Vector()
84 {
85 reinit(parallel_partitioning, communicator);
86 }
87
88
89
91 : Vector()
92 {
94 vector = std::make_unique<Epetra_FEVector>(*v.vector);
96 }
97
98
99
100 Vector::Vector(Vector &&v) noexcept
101 : Vector()
102 {
103 // initialize a minimal, valid object and swap
104 swap(v);
105 }
106
107
108
109 Vector::Vector(const IndexSet &parallel_partitioner,
110 const Vector & v,
111 const MPI_Comm &communicator)
112 : Vector()
113 {
114 AssertThrow(parallel_partitioner.size() ==
115 static_cast<size_type>(
117 ExcDimensionMismatch(parallel_partitioner.size(),
119 v.vector->Map())));
120
121 vector = std::make_unique<Epetra_FEVector>(
122 parallel_partitioner.make_trilinos_map(communicator, true));
123 reinit(v, false, true);
124 }
125
126
127
129 const IndexSet &ghost,
130 const MPI_Comm &communicator)
131 : Vector()
132 {
133 reinit(local, ghost, communicator, false);
134 }
135
136
137
138 void
140 {
141 // When we clear the vector, reset the pointer and generate an empty
142 // vector.
143 Epetra_Map map(0, 0, Epetra_MpiComm(MPI_COMM_SELF));
144
145 has_ghosts = false;
146 vector = std::make_unique<Epetra_FEVector>(map);
147 last_action = Zero;
148 }
149
150
151
152 void
153 Vector::reinit(const IndexSet &parallel_partitioner,
154 const MPI_Comm &communicator,
155 const bool /*omit_zeroing_entries*/)
156 {
157 nonlocal_vector.reset();
158
159 const bool overlapping =
160 !parallel_partitioner.is_ascending_and_one_to_one(communicator);
161
162 Epetra_Map map =
163 parallel_partitioner.make_trilinos_map(communicator, overlapping);
164
165 vector = std::make_unique<Epetra_FEVector>(map);
166
167 has_ghosts = vector->Map().UniqueGIDs() == false;
168
169 // If the IndexSets are overlapping, we don't really know
170 // which process owns what. So we decide that no process
171 // owns anything in that case. In particular asking for
172 // the locally owned elements is not allowed.
173 if (has_ghosts)
174 {
177 }
178 else
179 owned_elements = parallel_partitioner;
180
181# ifdef DEBUG
182 const size_type n_elements_global =
184
185 Assert(has_ghosts || n_elements_global == size(), ExcInternalError());
186# endif
187
188 last_action = Zero;
189 }
190
191
192
193 void
195 const bool omit_zeroing_entries,
196 const bool allow_different_maps)
197 {
198 nonlocal_vector.reset();
199
200 // In case we do not allow to have different maps, this call means that
201 // we have to reset the vector. So clear the vector, initialize our map
202 // with the map in v, and generate the vector.
203 if (allow_different_maps == false)
204 {
205 // check equality for MPI communicators: We can only choose the fast
206 // version in case the underlying Epetra_MpiComm object is the same,
207 // otherwise we might access an MPI_Comm object that has been
208 // deleted
209 const Epetra_MpiComm *my_comm =
210 dynamic_cast<const Epetra_MpiComm *>(&vector->Comm());
211 const Epetra_MpiComm *v_comm =
212 dynamic_cast<const Epetra_MpiComm *>(&v.vector->Comm());
213 const bool same_communicators =
214 my_comm != nullptr && v_comm != nullptr &&
215 my_comm->DataPtr() == v_comm->DataPtr();
216 if (!same_communicators ||
217 vector->Map().SameAs(v.vector->Map()) == false)
218 {
219 vector = std::make_unique<Epetra_FEVector>(v.vector->Map());
221 last_action = Zero;
223 }
224 else if (omit_zeroing_entries == false)
225 {
226 // old and new vectors have exactly the same map, i.e. size and
227 // parallel distribution
228 int ierr;
229 ierr = vector->GlobalAssemble(last_action);
230 (void)ierr;
231 Assert(ierr == 0, ExcTrilinosError(ierr));
232
233 ierr = vector->PutScalar(0.0);
234 Assert(ierr == 0, ExcTrilinosError(ierr));
235
236 last_action = Zero;
237 }
238 }
239
240 // Otherwise, we have to check that the two vectors are already of the
241 // same size, create an object for the data exchange and then insert all
242 // the data. The first assertion is only a check whether the user knows
243 // what they are doing.
244 else
245 {
246 Assert(omit_zeroing_entries == false,
248 "It is not possible to exchange data with the "
249 "option 'omit_zeroing_entries' set, which would not write "
250 "elements."));
251
252 AssertThrow(size() == v.size(),
254
255 Epetra_Import data_exchange(vector->Map(), v.vector->Map());
256
257 const int ierr = vector->Import(*v.vector, data_exchange, Insert);
258 AssertThrow(ierr == 0, ExcTrilinosError(ierr));
259
260 last_action = Insert;
261 }
262# ifdef DEBUG
263 const Epetra_MpiComm *comm_ptr =
264 dynamic_cast<const Epetra_MpiComm *>(&(v.vector->Comm()));
265 Assert(comm_ptr != nullptr, ExcInternalError());
266 const size_type n_elements_global =
267 Utilities::MPI::sum(owned_elements.n_elements(), comm_ptr->Comm());
268 Assert(has_ghosts || n_elements_global == size(), ExcInternalError());
269# endif
270 }
271
272
273
274 void
275 Vector::reinit(const MPI::BlockVector &v, const bool import_data)
276 {
277 nonlocal_vector.reset();
280
281 // In case we do not allow to have different maps, this call means that
282 // we have to reset the vector. So clear the vector, initialize our map
283 // with the map in v, and generate the vector.
284 if (v.n_blocks() == 0)
285 return;
286
287 // create a vector that holds all the elements contained in the block
288 // vector. need to manually create an Epetra_Map.
289 size_type n_elements = 0, added_elements = 0, block_offset = 0;
290 for (size_type block = 0; block < v.n_blocks(); ++block)
291 n_elements += v.block(block).local_size();
292 std::vector<TrilinosWrappers::types::int_type> global_ids(n_elements, -1);
293 for (size_type block = 0; block < v.n_blocks(); ++block)
294 {
295 TrilinosWrappers::types::int_type *glob_elements =
297 v.block(block).trilinos_partitioner());
298 for (size_type i = 0; i < v.block(block).local_size(); ++i)
299 global_ids[added_elements++] = glob_elements[i] + block_offset;
300 owned_elements.add_indices(v.block(block).owned_elements,
301 block_offset);
302 block_offset += v.block(block).size();
303 }
304
305 Assert(n_elements == added_elements, ExcInternalError());
306 Epetra_Map new_map(v.size(),
307 n_elements,
308 global_ids.data(),
309 0,
310 v.block(0).trilinos_partitioner().Comm());
311
312 auto actual_vec = std::make_unique<Epetra_FEVector>(new_map);
313
314 TrilinosScalar *entries = (*actual_vec)[0];
315 for (size_type block = 0; block < v.n_blocks(); ++block)
316 {
317 v.block(block).trilinos_vector().ExtractCopy(entries, 0);
318 entries += v.block(block).local_size();
319 }
320
321 if (import_data == true)
322 {
324 *actual_vec)) == v.size(),
326 *actual_vec),
327 v.size()));
328
329 Epetra_Import data_exchange(vector->Map(), actual_vec->Map());
330
331 const int ierr = vector->Import(*actual_vec, data_exchange, Insert);
332 AssertThrow(ierr == 0, ExcTrilinosError(ierr));
333
334 last_action = Insert;
335 }
336 else
337 vector = std::move(actual_vec);
338# ifdef DEBUG
339 const Epetra_MpiComm *comm_ptr =
340 dynamic_cast<const Epetra_MpiComm *>(&(vector->Comm()));
341 Assert(comm_ptr != nullptr, ExcInternalError());
342 const size_type n_elements_global =
343 Utilities::MPI::sum(owned_elements.n_elements(), comm_ptr->Comm());
344
345 Assert(has_ghosts || n_elements_global == size(), ExcInternalError());
346# endif
347 }
348
349
350
351 void
352 Vector::reinit(const IndexSet &locally_owned_entries,
353 const IndexSet &ghost_entries,
354 const MPI_Comm &communicator,
355 const bool vector_writable)
356 {
357 nonlocal_vector.reset();
358 owned_elements = locally_owned_entries;
359 if (vector_writable == false)
360 {
361 IndexSet parallel_partitioner = locally_owned_entries;
362 parallel_partitioner.add_indices(ghost_entries);
363 Epetra_Map map =
364 parallel_partitioner.make_trilinos_map(communicator, true);
365 vector = std::make_unique<Epetra_FEVector>(map);
366 }
367 else
368 {
369 Epetra_Map map =
370 locally_owned_entries.make_trilinos_map(communicator, true);
371 Assert(map.IsOneToOne(),
372 ExcMessage("A writable vector must not have ghost entries in "
373 "its parallel partitioning"));
374
375 if (vector->Map().SameAs(map) == false)
376 vector = std::make_unique<Epetra_FEVector>(map);
377 else
378 {
379 const int ierr = vector->PutScalar(0.);
380 (void)ierr;
381 Assert(ierr == 0, ExcTrilinosError(ierr));
382 }
383
384 IndexSet nonlocal_entries(ghost_entries);
385 nonlocal_entries.subtract_set(locally_owned_entries);
386 if (Utilities::MPI::n_mpi_processes(communicator) > 1)
387 {
388 Epetra_Map nonlocal_map =
389 nonlocal_entries.make_trilinos_map(communicator, true);
391 std::make_unique<Epetra_MultiVector>(nonlocal_map, 1);
392 }
393 }
394
395 has_ghosts = vector->Map().UniqueGIDs() == false;
396
397 last_action = Zero;
398
399# ifdef DEBUG
400 const size_type n_elements_global =
402
403 Assert(has_ghosts || n_elements_global == size(), ExcInternalError());
404# endif
405 }
406
407
408
409 Vector &
411 {
412 Assert(vector.get() != nullptr,
413 ExcMessage("Vector is not constructed properly."));
414
415 // check equality for MPI communicators to avoid accessing a possibly
416 // invalid MPI_Comm object
417 const Epetra_MpiComm *my_comm =
418 dynamic_cast<const Epetra_MpiComm *>(&vector->Comm());
419 const Epetra_MpiComm *v_comm =
420 dynamic_cast<const Epetra_MpiComm *>(&v.vector->Comm());
421 const bool same_communicators = my_comm != nullptr && v_comm != nullptr &&
422 my_comm->DataPtr() == v_comm->DataPtr();
423 // Need to ask MPI whether the communicators are the same. We would like
424 // to use the following checks but currently we cannot make sure the
425 // memory of my_comm is not stale from some MPI_Comm_free
426 // somewhere. This can happen when a vector lives in GrowingVectorMemory
427 // data structures. Thus, the following code is commented out.
428 //
429 // if (my_comm != nullptr &&
430 // v_comm != nullptr &&
431 // my_comm->DataPtr() != v_comm->DataPtr())
432 // {
433 // int communicators_same = 0;
434 // const int ierr = MPI_Comm_compare (my_comm->GetMpiComm(),
435 // v_comm->GetMpiComm(),
436 // &communicators_same);
437 // AssertThrowMPI(ierr);
438 // if (!(communicators_same == MPI_IDENT ||
439 // communicators_same == MPI_CONGRUENT))
440 // same_communicators = false;
441 // else
442 // same_communicators = true;
443 // }
444
445 // distinguish three cases. First case: both vectors have the same
446 // layout (just need to copy the local data, not reset the memory and
447 // the underlying Epetra_Map). The third case means that we have to
448 // rebuild the calling vector.
449 if (same_communicators && v.vector->Map().SameAs(vector->Map()))
450 {
451 *vector = *v.vector;
452 if (v.nonlocal_vector.get() != nullptr)
454 std::make_unique<Epetra_MultiVector>(v.nonlocal_vector->Map(), 1);
455 last_action = Zero;
456 }
457 // Second case: vectors have the same global
458 // size, but different parallel layouts (and
459 // one of them a one-to-one mapping). Then we
460 // can call the import/export functionality.
461 else if (size() == v.size() &&
462 (v.vector->Map().UniqueGIDs() || vector->Map().UniqueGIDs()))
463 {
464 reinit(v, false, true);
465 }
466 // Third case: Vectors do not have the same
467 // size.
468 else
469 {
470 vector = std::make_unique<Epetra_FEVector>(*v.vector);
471 last_action = Zero;
474 }
475
476 if (v.nonlocal_vector.get() != nullptr)
478 std::make_unique<Epetra_MultiVector>(v.nonlocal_vector->Map(), 1);
479
480 return *this;
481 }
482
483
484
485 Vector &
487 {
488 swap(v);
489 return *this;
490 }
491
492
493
494 template <typename number>
495 Vector &
496 Vector::operator=(const ::Vector<number> &v)
497 {
498 Assert(size() == v.size(), ExcDimensionMismatch(size(), v.size()));
499
500 // this is probably not very efficient but works. in particular, we could
501 // do better if we know that number==TrilinosScalar because then we could
502 // elide the copying of elements
503 //
504 // let's hope this isn't a particularly frequent operation
505 std::pair<size_type, size_type> local_range = this->local_range();
506 for (size_type i = local_range.first; i < local_range.second; ++i)
507 (*vector)[0][i - local_range.first] = v(i);
508
509 return *this;
510 }
511
512
513
514 void
516 const Vector & v)
517 {
518 Assert(m.trilinos_matrix().Filled() == true,
519 ExcMessage("Matrix is not compressed. "
520 "Cannot find exchange information!"));
521 Assert(v.vector->Map().UniqueGIDs() == true,
522 ExcMessage("The input vector has overlapping data, "
523 "which is not allowed."));
524
525 if (vector->Map().SameAs(m.trilinos_matrix().ColMap()) == false)
526 vector =
527 std::make_unique<Epetra_FEVector>(m.trilinos_matrix().ColMap());
528
529 Epetra_Import data_exchange(vector->Map(), v.vector->Map());
530 const int ierr = vector->Import(*v.vector, data_exchange, Insert);
531
532 AssertThrow(ierr == 0, ExcTrilinosError(ierr));
533
534 last_action = Insert;
535 }
536
537
538 void
540 const VectorOperation::values operation)
541 {
542 Assert(
543 this->size() == rwv.size(),
545 "Both vectors need to have the same size for import() to work!"));
546 // TODO: a generic import() function should handle any kind of data layout
547 // in ReadWriteVector, but this function is of limited use as this class
548 // will (hopefully) be retired eventually.
551
552 if (operation == VectorOperation::insert)
553 {
554 for (const auto idx : this->locally_owned_elements())
555 (*this)[idx] = rwv[idx];
556 }
557 else if (operation == VectorOperation::add)
558 {
559 for (const auto idx : this->locally_owned_elements())
560 (*this)[idx] += rwv[idx];
561 }
562 else
564
565 this->compress(operation);
566 }
567
568
569 void
571 {
572 // Select which mode to send to Trilinos. Note that we use last_action if
573 // available and ignore what the user tells us to detect wrongly mixed
574 // operations. Typically given_last_action is only used on machines that
575 // do not execute an operation (because they have no own cells for
576 // example).
577 Epetra_CombineMode mode = last_action;
578 if (last_action == Zero)
579 {
580 if (given_last_action == ::VectorOperation::add)
581 mode = Add;
582 else if (given_last_action == ::VectorOperation::insert)
583 mode = Insert;
584 else
585 Assert(
586 false,
588 "compress() can only be called with VectorOperation add, insert, or unknown"));
589 }
590 else
591 {
592 Assert(
593 ((last_action == Add) &&
594 (given_last_action == ::VectorOperation::add)) ||
595 ((last_action == Insert) &&
596 (given_last_action == ::VectorOperation::insert)),
598 "The last operation on the Vector and the given last action in the compress() call do not agree!"));
599 }
600
601
602# ifdef DEBUG
603 // check that every process has decided to use the same mode. This will
604 // otherwise result in undefined behavior in the call to
605 // GlobalAssemble().
606 const double double_mode = mode;
607 const Epetra_MpiComm *comm_ptr =
608 dynamic_cast<const Epetra_MpiComm *>(&(trilinos_partitioner().Comm()));
609 Assert(comm_ptr != nullptr, ExcInternalError());
610
611 const Utilities::MPI::MinMaxAvg result =
612 Utilities::MPI::min_max_avg(double_mode, comm_ptr->GetMpiComm());
613 Assert(result.max == result.min,
615 "Not all processors agree whether the last operation on "
616 "this vector was an addition or a set operation. This will "
617 "prevent the compress() operation from succeeding."));
618
619# endif
620
621 // Now pass over the information about what we did last to the vector.
622 int ierr = 0;
623 if (nonlocal_vector.get() == nullptr || mode != Add)
624 ierr = vector->GlobalAssemble(mode);
625 else
626 {
627 Epetra_Export exporter(nonlocal_vector->Map(), vector->Map());
628 ierr = vector->Export(*nonlocal_vector, exporter, mode);
629 AssertThrow(ierr == 0, ExcTrilinosError(ierr));
630 ierr = nonlocal_vector->PutScalar(0.);
631 }
632 AssertThrow(ierr == 0, ExcTrilinosError(ierr));
633 last_action = Zero;
634
635 compressed = true;
636 }
637
638
639
641 Vector::operator()(const size_type index) const
642 {
643 // Extract local indices in the vector.
644 TrilinosWrappers::types::int_type trilinos_i = vector->Map().LID(
645 static_cast<TrilinosWrappers::types::int_type>(index));
646 TrilinosScalar value = 0.;
647
648 // If the element is not present on the current processor, we can't
649 // continue. This is the main difference to the el() function.
650 if (trilinos_i == -1)
651 {
652 Assert(false,
654 local_size(),
655 vector->Map().MinMyGID(),
656 vector->Map().MaxMyGID()));
657 }
658 else
659 value = (*vector)[0][trilinos_i];
660
661 return value;
662 }
663
664
665
666 void
667 Vector::add(const Vector &v, const bool allow_different_maps)
668 {
669 if (allow_different_maps == false)
670 *this += v;
671 else
672 {
674 AssertThrow(size() == v.size(),
676
677# if DEAL_II_TRILINOS_VERSION_GTE(11, 11, 0)
678 Epetra_Import data_exchange(vector->Map(), v.vector->Map());
679 int ierr =
680 vector->Import(*v.vector, data_exchange, Epetra_AddLocalAlso);
681 AssertThrow(ierr == 0, ExcTrilinosError(ierr));
682 last_action = Add;
683# else
684 // In versions older than 11.11 the Import function is broken for
685 // adding Hence, we provide a workaround in this case
686
687 Epetra_MultiVector dummy(vector->Map(), 1, false);
688 Epetra_Import data_exchange(dummy.Map(), v.vector->Map());
689
690 int ierr = dummy.Import(*v.vector, data_exchange, Insert);
691 AssertThrow(ierr == 0, ExcTrilinosError(ierr));
692
693 ierr = vector->Update(1.0, dummy, 1.0);
694 AssertThrow(ierr == 0, ExcTrilinosError(ierr));
695# endif
696 }
697 }
698
699
700
701 bool
703 {
704 Assert(size() == v.size(), ExcDimensionMismatch(size(), v.size()));
705 if (local_size() != v.local_size())
706 return false;
707
708 size_type i;
709 for (i = 0; i < local_size(); ++i)
710 if ((*(v.vector))[0][i] != (*vector)[0][i])
711 return false;
712
713 return true;
714 }
715
716
717
718 bool
720 {
721 Assert(size() == v.size(), ExcDimensionMismatch(size(), v.size()));
722
723 return (!(*this == v));
724 }
725
726
727
728 bool
730 {
731 // get a representation of the vector and
732 // loop over all the elements
733 TrilinosScalar * start_ptr = (*vector)[0];
734 const TrilinosScalar *ptr = start_ptr, *eptr = start_ptr + local_size();
735 unsigned int flag = 0;
736 while (ptr != eptr)
737 {
738 if (*ptr != 0)
739 {
740 flag = 1;
741 break;
742 }
743 ++ptr;
744 }
745
746 // in parallel, check that the vector
747 // is zero on _all_ processors.
748 const Epetra_MpiComm *mpi_comm =
749 dynamic_cast<const Epetra_MpiComm *>(&vector->Map().Comm());
750 Assert(mpi_comm != nullptr, ExcInternalError());
751 unsigned int num_nonzero = Utilities::MPI::sum(flag, mpi_comm->Comm());
752 return num_nonzero == 0;
753 }
754
755
756
757 bool
759 {
760 // get a representation of the vector and
761 // loop over all the elements
762 TrilinosScalar * start_ptr = (*vector)[0];
763 const TrilinosScalar *ptr = start_ptr, *eptr = start_ptr + local_size();
764 unsigned int flag = 0;
765 while (ptr != eptr)
766 {
767 if (*ptr < 0.0)
768 {
769 flag = 1;
770 break;
771 }
772 ++ptr;
773 }
774
775 // in parallel, check that the vector
776 // is zero on _all_ processors.
777 const auto max_n_negative =
779 return max_n_negative == 0;
780 }
781
782
783
784 void
785 Vector::print(std::ostream & out,
786 const unsigned int precision,
787 const bool scientific,
788 const bool across) const
789 {
790 AssertThrow(out.fail() == false, ExcIO());
791 boost::io::ios_flags_saver restore_flags(out);
792
793
794 out.precision(precision);
795 if (scientific)
796 out.setf(std::ios::scientific, std::ios::floatfield);
797 else
798 out.setf(std::ios::fixed, std::ios::floatfield);
799
800 if (size() != local_size())
801 {
802 auto global_id = [&](const size_type index) {
803 return gid(vector->Map(), index);
804 };
805 out << "size:" << size() << " local_size:" << local_size() << " :"
806 << std::endl;
807 for (size_type i = 0; i < local_size(); ++i)
808 out << "[" << global_id(i) << "]: " << (*(vector))[0][i]
809 << std::endl;
810 }
811 else
812 {
813 TrilinosScalar *val;
814 int leading_dimension;
815 int ierr = vector->ExtractView(&val, &leading_dimension);
816
817 AssertThrow(ierr == 0, ExcTrilinosError(ierr));
818 if (across)
819 for (size_type i = 0; i < size(); ++i)
820 out << static_cast<double>(val[i]) << ' ';
821 else
822 for (size_type i = 0; i < size(); ++i)
823 out << static_cast<double>(val[i]) << std::endl;
824 out << std::endl;
825 }
826
827 AssertThrow(out.fail() == false, ExcIO());
828 }
829
830
831
832 void
834 {
841 }
842
843
844
845 std::size_t
847 {
848 // TODO[TH]: No accurate memory
849 // consumption for Trilinos vectors
850 // yet. This is a rough approximation with
851 // one index and the value per local
852 // entry.
853 return sizeof(*this) +
854 this->local_size() *
855 (sizeof(double) + sizeof(TrilinosWrappers::types::int_type));
856 }
857
858 // explicit instantiations
859# ifndef DOXYGEN
860# include "trilinos_vector.inst"
861# endif
862 } // namespace MPI
863} // namespace TrilinosWrappers
864
866
867#endif // DEAL_II_WITH_TRILINOS
size_type size() const
Definition: index_set.h:1636
size_type n_elements() const
Definition: index_set.h:1834
void set_size(const size_type size)
Definition: index_set.h:1624
void clear()
Definition: index_set.h:1612
void subtract_set(const IndexSet &other)
Definition: index_set.cc:266
Epetra_Map make_trilinos_map(const MPI_Comm &communicator=MPI_COMM_WORLD, const bool overlapping=false) const
Definition: index_set.cc:601
bool is_ascending_and_one_to_one(const MPI_Comm &communicator) const
Definition: index_set.cc:671
void add_indices(const ForwardIterator &begin, const ForwardIterator &end)
Definition: index_set.h:1705
const Epetra_CrsMatrix & trilinos_matrix() const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:442
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:443
static ::ExceptionBase & ExcIO()
static ::ExceptionBase & ExcGhostsPresent()
unsigned int n_blocks() const
static ::ExceptionBase & ExcNotImplemented()
std::size_t size() const
#define Assert(cond, exc)
Definition: exceptions.h:1473
#define AssertIndexRange(index, range)
Definition: exceptions.h:1732
static ::ExceptionBase & ExcInternalError()
BlockType & block(const unsigned int i)
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
Definition: exceptions.h:1583
void add(const std::vector< size_type > &indices, const std::vector< TrilinosScalar > &values)
const MPI_Comm & get_mpi_communicator() const
std::unique_ptr< Epetra_MultiVector > nonlocal_vector
const Epetra_BlockMap & trilinos_partitioner() const
int gid(const Epetra_BlockMap &map, int i)
reference operator()(const size_type index)
void import_nonlocal_data_for_fe(const ::TrilinosWrappers::SparseMatrix &matrix, const Vector &vector)
std::unique_ptr< Epetra_FEVector > vector
size_type local_size() const
void print(std::ostream &out, const unsigned int precision=3, const bool scientific=true, const bool across=true) const
IndexSet locally_owned_elements() const
void swap(BlockVector &u, BlockVector &v)
static ::ExceptionBase & ExcAccessToNonLocalElement(size_type arg1, size_type arg2, size_type arg3, size_type arg4)
void reinit(const Vector &v, const bool omit_zeroing_entries=false, const bool allow_different_maps=false)
std::pair< size_type, size_type > local_range() const
void compress(::VectorOperation::values operation)
bool operator!=(const Vector &v) const
bool operator==(const Vector &v) const
static ::ExceptionBase & ExcTrilinosError(int arg1)
Vector & operator=(const TrilinosScalar s)
void import(const LinearAlgebra::ReadWriteVector< double > &rwv, const VectorOperation::values operation)
std::size_t memory_consumption() const
const IndexSet & get_stored_elements() const
void swap(MemorySpaceData< Number, MemorySpace > &, MemorySpaceData< Number, MemorySpace > &)
TrilinosWrappers::types::int_type global_length(const Epetra_MultiVector &vector)
TrilinosWrappers::types::int_type * my_global_elements(const Epetra_BlockMap &map)
TrilinosWrappers::types::int64_type n_global_elements(const Epetra_BlockMap &map)
T sum(const T &t, const MPI_Comm &mpi_communicator)
unsigned int n_mpi_processes(const MPI_Comm &mpi_communicator)
Definition: mpi.cc:140
MinMaxAvg min_max_avg(const double my_value, const MPI_Comm &mpi_communicator)
Definition: mpi.cc:114
T max(const T &t, const MPI_Comm &mpi_communicator)
const Epetra_Comm & comm_self()
Definition: utilities.cc:1110
double TrilinosScalar
Definition: types.h:163