Reference documentation for deal.II version Git a64c6c9 2017-01-19 10:07:38 -0700
trilinos_vector.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2008 - 2016 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 at
12 // the top level of the deal.II distribution.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii__trilinos_vector_h
17 #define dealii__trilinos_vector_h
18 
19 
20 #include <deal.II/base/config.h>
21 
22 #ifdef DEAL_II_WITH_TRILINOS
23 
24 # include <deal.II/base/std_cxx11/shared_ptr.h>
25 # include <deal.II/base/subscriptor.h>
26 # include <deal.II/base/index_set.h>
27 # include <deal.II/base/utilities.h>
28 # include <deal.II/lac/exceptions.h>
29 # include <deal.II/lac/vector.h>
30 # include <deal.II/lac/trilinos_vector_base.h>
31 
32 DEAL_II_DISABLE_EXTRA_DIAGNOSTICS
33 # include "Epetra_Map.h"
34 # include "Epetra_LocalMap.h"
35 DEAL_II_ENABLE_EXTRA_DIAGNOSTICS
36 
37 DEAL_II_NAMESPACE_OPEN
38 
39 
40 // forward declaration
41 template <typename> class Vector;
42 
47 namespace TrilinosWrappers
48 {
49  class SparseMatrix;
50 
51  namespace
52  {
53 #ifndef DEAL_II_WITH_64BIT_INDICES
54  // define a helper function that queries the global ID of local ID of
55  // an Epetra_BlockMap object by calling either the 32- or 64-bit
56  // function necessary.
57  inline
58  int gid(const Epetra_BlockMap &map, int i)
59  {
60  return map.GID(i);
61  }
62 #else
63  // define a helper function that queries the global ID of local ID of
64  // an Epetra_BlockMap object by calling either the 32- or 64-bit
65  // function necessary.
66  inline
67  long long int gid(const Epetra_BlockMap &map, int i)
68  {
69  return map.GID64(i);
70  }
71 #endif
72  }
73 
82  namespace MPI
83  {
84  class BlockVector;
85 
253  class Vector : public VectorBase
254  {
255  public:
259  typedef ::types::global_dof_index size_type;
260 
271  static const bool supports_distributed_data = true;
272 
282  Vector ();
283 
287  Vector (const Vector &v);
288 
289 #ifdef DEAL_II_WITH_CXX11
290 
297  Vector (Vector &&v);
298 #endif
299 
303  ~Vector ();
304 
328  void reinit (const VectorBase &v,
329  const bool omit_zeroing_entries = false,
330  const bool allow_different_maps = false);
331 
335  void reinit (const BlockVector &v,
336  const bool import_data = false);
337 
344  Vector &operator= (const TrilinosScalar s);
345 
351  Vector &operator= (const Vector &v);
352 
353 #ifdef DEAL_II_WITH_CXX11
354 
361  Vector &operator= (Vector &&v);
362 #endif
363 
371  Vector &operator= (const ::TrilinosWrappers::Vector &v);
372 
381  template <typename Number>
382  Vector &operator= (const ::Vector<Number> &v);
383 
402  (const ::TrilinosWrappers::SparseMatrix &matrix,
403  const Vector &vector);
405 
425  explicit Vector (const Epetra_Map &parallel_partitioning) DEAL_II_DEPRECATED;
426 
443  Vector (const Epetra_Map &parallel_partitioning,
444  const VectorBase &v) DEAL_II_DEPRECATED;
445 
460  template <typename number>
461  void reinit (const Epetra_Map &parallel_partitioner,
462  const ::Vector<number> &v) DEAL_II_DEPRECATED;
463 
478  void reinit (const Epetra_Map &parallel_partitioning,
479  const bool omit_zeroing_entries = false) DEAL_II_DEPRECATED;
480 
495  template <typename Number>
496  Vector (const Epetra_Map &parallel_partitioning,
497  const ::Vector<Number> &v) DEAL_II_DEPRECATED;
499 
521  explicit Vector (const IndexSet &parallel_partitioning,
522  const MPI_Comm &communicator = MPI_COMM_WORLD);
523 
535  Vector (const IndexSet &local,
536  const IndexSet &ghost,
537  const MPI_Comm &communicator = MPI_COMM_WORLD);
538 
553  Vector (const IndexSet &parallel_partitioning,
554  const VectorBase &v,
555  const MPI_Comm &communicator = MPI_COMM_WORLD);
556 
569  template <typename Number>
570  Vector (const IndexSet &parallel_partitioning,
571  const ::Vector<Number> &v,
572  const MPI_Comm &communicator = MPI_COMM_WORLD);
573 
597  void reinit (const IndexSet &parallel_partitioning,
598  const MPI_Comm &communicator = MPI_COMM_WORLD,
599  const bool omit_zeroing_entries = false);
600 
626  void reinit (const IndexSet &locally_owned_entries,
627  const IndexSet &ghost_entries,
628  const MPI_Comm &communicator = MPI_COMM_WORLD,
629  const bool vector_writable = false);
631  };
632 
633 
634 
635 
636 // ------------------- inline and template functions --------------
637 
638 
647  inline
648  void swap (Vector &u, Vector &v)
649  {
650  u.swap (v);
651  }
652 
653 
654 #ifndef DOXYGEN
655 
656  template <typename number>
657  Vector::Vector (const Epetra_Map &input_map,
658  const ::Vector<number> &v)
659  {
660  reinit (input_map, v);
661  }
662 
663 
664 
665  template <typename number>
666  Vector::Vector (const IndexSet &parallel_partitioner,
667  const ::Vector<number> &v,
668  const MPI_Comm &communicator)
669  {
670  *this = Vector(parallel_partitioner.make_trilinos_map (communicator, true),
671  v);
672  owned_elements = parallel_partitioner;
673  }
674 
675 
676 
677 
678  template <typename number>
679  void Vector::reinit (const Epetra_Map &parallel_partitioner,
680  const ::Vector<number> &v)
681  {
682  if (vector.get() == 0 || vector->Map().SameAs(parallel_partitioner) == false)
683  reinit(parallel_partitioner);
684 
685  const int size = parallel_partitioner.NumMyElements();
686 
687  // Need to copy out values, since the deal.II might not use doubles, so
688  // that a direct access is not possible.
689  for (int i=0; i<size; ++i)
690  (*vector)[0][i] = v(gid(parallel_partitioner,i));
691 
692  }
693 
694 
695  inline
696  Vector &
697  Vector::operator= (const TrilinosScalar s)
698  {
700 
701  return *this;
702  }
703 
704 
705  template <typename Number>
706  Vector &
707  Vector::operator= (const ::Vector<Number> &v)
708  {
709  if (size() != v.size())
710  {
711  vector.reset (new Epetra_FEVector(Epetra_Map
712  (static_cast<TrilinosWrappers::types::int_type>(v.size()), 0,
713 #ifdef DEAL_II_WITH_MPI
714  Epetra_MpiComm(MPI_COMM_SELF)
715 #else
716  Epetra_SerialComm()
717 #endif
718  )));
719  }
720 
721  const Epetra_Map &map = vector_partitioner();
722  const int size = map.NumMyElements();
723 
724  // Need to copy out values, since the deal.II might not use doubles, so
725  // that a direct access is not possible.
726  for (int i=0; i<size; ++i)
727  (*vector)[0][i] = v(gid(map,i));
728 
729  return *this;
730  }
731 
732 
733 #endif /* DOXYGEN */
734 
735  } /* end of namespace MPI */
736 
737 
738 
752  class Vector : public VectorBase
753  {
754  public:
758  typedef ::types::global_dof_index size_type;
759 
770  static const bool supports_distributed_data = false;
771 
777  Vector () DEAL_II_DEPRECATED;
778 
782  explicit Vector (const size_type n) DEAL_II_DEPRECATED;
783 
793  explicit Vector (const Epetra_Map &partitioning) DEAL_II_DEPRECATED;
794 
804  explicit Vector (const IndexSet &partitioning,
805  const MPI_Comm &communicator = MPI_COMM_WORLD) DEAL_II_DEPRECATED;
806 
811  explicit Vector (const VectorBase &V) DEAL_II_DEPRECATED;
812 
817  template <typename Number>
818  explicit Vector (const ::Vector<Number> &v) DEAL_II_DEPRECATED;
819 
829  void reinit (const size_type n,
830  const bool omit_zeroing_entries = false);
831 
849  void reinit (const Epetra_Map &input_map,
850  const bool omit_zeroing_entries = false);
851 
869  void reinit (const IndexSet &input_map,
870  const MPI_Comm &communicator = MPI_COMM_WORLD,
871  const bool omit_zeroing_entries = false);
872 
886  void reinit (const VectorBase &V,
887  const bool omit_zeroing_entries = false,
888  const bool allow_different_maps = false);
889 
896  Vector &operator= (const TrilinosScalar s);
897 
902  Vector &operator= (const MPI::Vector &v);
903 
907  template <typename Number>
908  Vector &operator= (const ::Vector<Number> &v);
909 
914  Vector &operator= (const Vector &v);
915 
927  void update_ghost_values () const;
928  };
929 
930 
931 
932 // ------------------- inline and template functions --------------
933 
934 
943  inline
944  void swap (Vector &u, Vector &v)
945  {
946  u.swap (v);
947  }
948 
949 
950 #ifndef DOXYGEN
951 
952  template <typename number>
953  Vector::Vector (const ::Vector<number> &v)
954  {
955  Epetra_LocalMap map ((TrilinosWrappers::types::int_type)v.size(), 0, Utilities::Trilinos::comm_self());
956  vector.reset (new Epetra_FEVector(map));
957  *this = v;
958  }
959 
960 
961 
962  inline
963  Vector &
964  Vector::operator= (const TrilinosScalar s)
965  {
967 
968  return *this;
969  }
970 
971 
972 
973  template <typename Number>
974  Vector &
975  Vector::operator= (const ::Vector<Number> &v)
976  {
977  if (size() != v.size())
978  {
979  vector.reset();
980 
981  Epetra_LocalMap map ((TrilinosWrappers::types::int_type)v.size(), 0,
983  vector.reset (new Epetra_FEVector(map));
984  owned_elements = v.locally_owned_elements();
985  }
986 
987  const Epetra_Map &map = vector_partitioner();
988  const TrilinosWrappers::types::int_type size = map.NumMyElements();
989 
990  Assert (map.MaxLID() == size-1,
991  ExcDimensionMismatch(map.MaxLID(), size-1));
992 
993  // Need to copy out values, since the
994  // deal.II might not use doubles, so
995  // that a direct access is not possible.
996  for (TrilinosWrappers::types::int_type i=0; i<size; ++i)
997  (*vector)[0][i] = v(i);
998 
999  return *this;
1000  }
1001 
1002 
1003 
1004  inline
1005  void
1007  {}
1008 
1009 
1010 #endif
1011 
1012 } /* namespace TrilinosWrappers */
1013 
1017 namespace internal
1018 {
1019  namespace LinearOperator
1020  {
1021  template <typename> class ReinitHelper;
1022 
1027  template<>
1028  class ReinitHelper<TrilinosWrappers::MPI::Vector>
1029  {
1030  public:
1031  template <typename Matrix>
1032  static
1033  void reinit_range_vector (const Matrix &matrix,
1035  bool omit_zeroing_entries)
1036  {
1037  v.reinit(matrix.locally_owned_range_indices(), matrix.get_mpi_communicator(), omit_zeroing_entries);
1038  }
1039 
1040  template <typename Matrix>
1041  static
1042  void reinit_domain_vector(const Matrix &matrix,
1044  bool omit_zeroing_entries)
1045  {
1046  v.reinit(matrix.locally_owned_domain_indices(), matrix.get_mpi_communicator(), omit_zeroing_entries);
1047  }
1048  };
1049 
1054  template<>
1056  {
1057  public:
1058  template <typename Matrix>
1059  static
1060  void reinit_range_vector (const Matrix &matrix,
1062  bool omit_zeroing_entries)
1063  {
1064  v.reinit(matrix.locally_owned_range_indices(),
1065  matrix.get_mpi_communicator(),
1066  omit_zeroing_entries);
1067  }
1068 
1069  template <typename Matrix>
1070  static
1071  void reinit_domain_vector(const Matrix &matrix,
1073  bool omit_zeroing_entries)
1074  {
1075  v.reinit(matrix.locally_owned_domain_indices(),
1076  matrix.get_mpi_communicator(),
1077  omit_zeroing_entries);
1078  }
1079  };
1080 
1081  } /* namespace LinearOperator */
1082 } /* namespace internal */
1083 
1084 
1085 DEAL_II_NAMESPACE_CLOSE
1086 
1087 #endif // DEAL_II_WITH_TRILINOS
1088 
1089 /*---------------------------- trilinos_vector.h ---------------------------*/
1090 
1091 #endif
1092 /*---------------------------- trilinos_vector.h ---------------------------*/
Vector & operator=(const TrilinosScalar s)
::types::global_dof_index size_type
const Epetra_Comm & comm_self()
Definition: utilities.cc:734
void reinit(const VectorBase &v, const bool omit_zeroing_entries=false, const bool allow_different_maps=false)
std_cxx11::shared_ptr< Epetra_FEVector > vector
#define Assert(cond, exc)
Definition: exceptions.h:313
static::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
void swap(Vector &u, Vector &v)
const Epetra_Map & vector_partitioner() const
void import_nonlocal_data_for_fe(const ::TrilinosWrappers::SparseMatrix &matrix, const Vector &vector)
Epetra_Map make_trilinos_map(const MPI_Comm &communicator=MPI_COMM_WORLD, const bool overlapping=false) const
Definition: index_set.cc:535
VectorBase & operator=(const TrilinosScalar s)
void reinit(const size_type n, const bool omit_zeroing_entries=false)
void update_ghost_values() const
::types::global_dof_index size_type
size_type size() const
static const bool supports_distributed_data