Reference documentation for deal.II version Git a670ecb 2017-04-28 20:25:58 -0400
trilinos_vector.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2008 - 2017 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/subscriptor.h>
25 # include <deal.II/base/index_set.h>
26 # include <deal.II/base/utilities.h>
27 # include <deal.II/lac/exceptions.h>
28 # include <deal.II/lac/vector.h>
29 # include <deal.II/lac/trilinos_vector_base.h>
30 # include <deal.II/lac/vector_type_traits.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 #include <memory>
38 
39 DEAL_II_NAMESPACE_OPEN
40 
41 
42 // forward declaration
43 template <typename> class Vector;
44 
49 namespace TrilinosWrappers
50 {
51  class SparseMatrix;
52 
53  namespace
54  {
55 #ifndef DEAL_II_WITH_64BIT_INDICES
56  // define a helper function that queries the global ID of local ID of
57  // an Epetra_BlockMap object by calling either the 32- or 64-bit
58  // function necessary.
59  inline
60  int gid(const Epetra_BlockMap &map, int i)
61  {
62  return map.GID(i);
63  }
64 #else
65  // define a helper function that queries the global ID of local ID of
66  // an Epetra_BlockMap object by calling either the 32- or 64-bit
67  // function necessary.
68  inline
69  long long int gid(const Epetra_BlockMap &map, int i)
70  {
71  return map.GID64(i);
72  }
73 #endif
74  }
75 
84  namespace MPI
85  {
86  class BlockVector;
87 
255  class Vector : public VectorBase
256  {
257  public:
261  typedef ::types::global_dof_index size_type;
262 
276  static const bool supports_distributed_data DEAL_II_DEPRECATED = true;
277 
287  Vector ();
288 
292  Vector (const Vector &v);
293 
294 #ifdef DEAL_II_WITH_CXX11
295 
302  Vector (Vector &&v);
303 #endif
304 
308  ~Vector ();
309 
333  void reinit (const VectorBase &v,
334  const bool omit_zeroing_entries = false,
335  const bool allow_different_maps = false);
336 
340  void reinit (const BlockVector &v,
341  const bool import_data = false);
342 
349  Vector &operator= (const TrilinosScalar s);
350 
356  Vector &operator= (const Vector &v);
357 
358 #ifdef DEAL_II_WITH_CXX11
359 
366  Vector &operator= (Vector &&v);
367 #endif
368 
376  Vector &operator= (const ::TrilinosWrappers::Vector &v);
377 
386  template <typename Number>
387  Vector &operator= (const ::Vector<Number> &v);
388 
407  (const ::TrilinosWrappers::SparseMatrix &matrix,
408  const Vector &vector);
410 
430  explicit Vector (const Epetra_Map &parallel_partitioning) DEAL_II_DEPRECATED;
431 
448  Vector (const Epetra_Map &parallel_partitioning,
449  const VectorBase &v) DEAL_II_DEPRECATED;
450 
465  template <typename number>
466  void reinit (const Epetra_Map &parallel_partitioner,
467  const ::Vector<number> &v) DEAL_II_DEPRECATED;
468 
483  void reinit (const Epetra_Map &parallel_partitioning,
484  const bool omit_zeroing_entries = false) DEAL_II_DEPRECATED;
485 
500  template <typename Number>
501  Vector (const Epetra_Map &parallel_partitioning,
502  const ::Vector<Number> &v) DEAL_II_DEPRECATED;
504 
526  explicit Vector (const IndexSet &parallel_partitioning,
527  const MPI_Comm &communicator = MPI_COMM_WORLD);
528 
540  Vector (const IndexSet &local,
541  const IndexSet &ghost,
542  const MPI_Comm &communicator = MPI_COMM_WORLD);
543 
558  Vector (const IndexSet &parallel_partitioning,
559  const VectorBase &v,
560  const MPI_Comm &communicator = MPI_COMM_WORLD);
561 
574  template <typename Number>
575  Vector (const IndexSet &parallel_partitioning,
576  const ::Vector<Number> &v,
577  const MPI_Comm &communicator = MPI_COMM_WORLD);
578 
602  void reinit (const IndexSet &parallel_partitioning,
603  const MPI_Comm &communicator = MPI_COMM_WORLD,
604  const bool omit_zeroing_entries = false);
605 
631  void reinit (const IndexSet &locally_owned_entries,
632  const IndexSet &ghost_entries,
633  const MPI_Comm &communicator = MPI_COMM_WORLD,
634  const bool vector_writable = false);
636  };
637 
638 
639 
640 
641 // ------------------- inline and template functions --------------
642 
643 
652  inline
653  void swap (Vector &u, Vector &v)
654  {
655  u.swap (v);
656  }
657 
658 
659 #ifndef DOXYGEN
660 
661  template <typename number>
662  Vector::Vector (const Epetra_Map &input_map,
663  const ::Vector<number> &v)
664  {
665  reinit (input_map, v);
666  }
667 
668 
669 
670  template <typename number>
671  Vector::Vector (const IndexSet &parallel_partitioner,
672  const ::Vector<number> &v,
673  const MPI_Comm &communicator)
674  {
675  *this = Vector(parallel_partitioner.make_trilinos_map (communicator, true),
676  v);
677  owned_elements = parallel_partitioner;
678  }
679 
680 
681 
682 
683  template <typename number>
684  void Vector::reinit (const Epetra_Map &parallel_partitioner,
685  const ::Vector<number> &v)
686  {
687  if (vector.get() == nullptr || vector->Map().SameAs(parallel_partitioner) == false)
688  reinit(parallel_partitioner);
689 
690  const int size = parallel_partitioner.NumMyElements();
691 
692  // Need to copy out values, since the deal.II might not use doubles, so
693  // that a direct access is not possible.
694  for (int i=0; i<size; ++i)
695  (*vector)[0][i] = v(gid(parallel_partitioner,i));
696 
697  }
698 
699 
700  inline
701  Vector &
702  Vector::operator= (const TrilinosScalar s)
703  {
705 
706  return *this;
707  }
708 
709 
710  template <typename Number>
711  Vector &
712  Vector::operator= (const ::Vector<Number> &v)
713  {
714  if (size() != v.size())
715  {
716  vector.reset (new Epetra_FEVector(Epetra_Map
717  (static_cast<TrilinosWrappers::types::int_type>(v.size()), 0,
718 #ifdef DEAL_II_WITH_MPI
719  Epetra_MpiComm(MPI_COMM_SELF)
720 #else
721  Epetra_SerialComm()
722 #endif
723  )));
724  }
725 
726  const Epetra_Map &map = vector_partitioner();
727  const int size = map.NumMyElements();
728 
729  // Need to copy out values, since the deal.II might not use doubles, so
730  // that a direct access is not possible.
731  for (int i=0; i<size; ++i)
732  (*vector)[0][i] = v(gid(map,i));
733 
734  return *this;
735  }
736 
737 
738 #endif /* DOXYGEN */
739 
740  } /* end of namespace MPI */
741 
742 
743 
757  class Vector : public VectorBase
758  {
759  public:
763  typedef ::types::global_dof_index size_type;
764 
778  static const bool supports_distributed_data DEAL_II_DEPRECATED = false;
779 
785  Vector () DEAL_II_DEPRECATED;
786 
790  explicit Vector (const size_type n) DEAL_II_DEPRECATED;
791 
801  explicit Vector (const Epetra_Map &partitioning) DEAL_II_DEPRECATED;
802 
812  explicit Vector (const IndexSet &partitioning,
813  const MPI_Comm &communicator = MPI_COMM_WORLD) DEAL_II_DEPRECATED;
814 
819  explicit Vector (const VectorBase &V) DEAL_II_DEPRECATED;
820 
825  template <typename Number>
826  explicit Vector (const ::Vector<Number> &v) DEAL_II_DEPRECATED;
827 
837  void reinit (const size_type n,
838  const bool omit_zeroing_entries = false);
839 
857  void reinit (const Epetra_Map &input_map,
858  const bool omit_zeroing_entries = false);
859 
877  void reinit (const IndexSet &input_map,
878  const MPI_Comm &communicator = MPI_COMM_WORLD,
879  const bool omit_zeroing_entries = false);
880 
894  void reinit (const VectorBase &V,
895  const bool omit_zeroing_entries = false,
896  const bool allow_different_maps = false);
897 
904  Vector &operator= (const TrilinosScalar s);
905 
910  Vector &operator= (const MPI::Vector &v);
911 
915  template <typename Number>
916  Vector &operator= (const ::Vector<Number> &v);
917 
922  Vector &operator= (const Vector &v);
923 
935  void update_ghost_values () const;
936  };
937 
938 
939 // ------------------- inline and template functions --------------
940 
941 
950  inline
951  void swap (Vector &u, Vector &v)
952  {
953  u.swap (v);
954  }
955 
956 
957 #ifndef DOXYGEN
958 
959  template <typename number>
960  Vector::Vector (const ::Vector<number> &v)
961  {
962  Epetra_LocalMap map ((TrilinosWrappers::types::int_type)v.size(), 0, Utilities::Trilinos::comm_self());
963  vector.reset (new Epetra_FEVector(map));
964  *this = v;
965  }
966 
967 
968 
969  inline
970  Vector &
971  Vector::operator= (const TrilinosScalar s)
972  {
974 
975  return *this;
976  }
977 
978 
979 
980  template <typename Number>
981  Vector &
982  Vector::operator= (const ::Vector<Number> &v)
983  {
984  if (size() != v.size())
985  {
986  vector.reset();
987 
988  Epetra_LocalMap map ((TrilinosWrappers::types::int_type)v.size(), 0,
990  vector.reset (new Epetra_FEVector(map));
991  owned_elements = v.locally_owned_elements();
992  }
993 
994  const Epetra_Map &map = vector_partitioner();
995  const TrilinosWrappers::types::int_type size = map.NumMyElements();
996 
997  Assert (map.MaxLID() == size-1,
998  ExcDimensionMismatch(map.MaxLID(), size-1));
999 
1000  // Need to copy out values, since the
1001  // deal.II might not use doubles, so
1002  // that a direct access is not possible.
1003  for (TrilinosWrappers::types::int_type i=0; i<size; ++i)
1004  (*vector)[0][i] = v(i);
1005 
1006  return *this;
1007  }
1008 
1009 
1010 
1011  inline
1012  void
1014  {}
1015 
1016 
1017 #endif
1018 
1019 } /* namespace TrilinosWrappers */
1020 
1024 namespace internal
1025 {
1026  namespace LinearOperator
1027  {
1028  template <typename> class ReinitHelper;
1029 
1034  template<>
1035  class ReinitHelper<TrilinosWrappers::MPI::Vector>
1036  {
1037  public:
1038  template <typename Matrix>
1039  static
1040  void reinit_range_vector (const Matrix &matrix,
1042  bool omit_zeroing_entries)
1043  {
1044  v.reinit(matrix.locally_owned_range_indices(), matrix.get_mpi_communicator(), omit_zeroing_entries);
1045  }
1046 
1047  template <typename Matrix>
1048  static
1049  void reinit_domain_vector(const Matrix &matrix,
1051  bool omit_zeroing_entries)
1052  {
1053  v.reinit(matrix.locally_owned_domain_indices(), matrix.get_mpi_communicator(), omit_zeroing_entries);
1054  }
1055  };
1056 
1061  template<>
1063  {
1064  public:
1065  template <typename Matrix>
1066  static
1067  void reinit_range_vector (const Matrix &matrix,
1069  bool omit_zeroing_entries)
1070  {
1071  v.reinit(matrix.locally_owned_range_indices(),
1072  matrix.get_mpi_communicator(),
1073  omit_zeroing_entries);
1074  }
1075 
1076  template <typename Matrix>
1077  static
1078  void reinit_domain_vector(const Matrix &matrix,
1080  bool omit_zeroing_entries)
1081  {
1082  v.reinit(matrix.locally_owned_domain_indices(),
1083  matrix.get_mpi_communicator(),
1084  omit_zeroing_entries);
1085  }
1086  };
1087 
1088  } /* namespace LinearOperator */
1089 } /* namespace internal */
1090 
1091 
1097 template <>
1098 struct is_serial_vector< TrilinosWrappers::Vector > : std::true_type
1099 {
1100 };
1101 
1102 
1108 template <>
1109 struct is_serial_vector< TrilinosWrappers::MPI::Vector > : std::false_type
1110 {
1111 };
1112 
1113 
1114 DEAL_II_NAMESPACE_CLOSE
1115 
1116 #endif // DEAL_II_WITH_TRILINOS
1117 
1118 /*---------------------------- trilinos_vector.h ---------------------------*/
1119 
1120 #endif
1121 /*---------------------------- trilinos_vector.h ---------------------------*/
Vector & operator=(const TrilinosScalar s)
::types::global_dof_index size_type
std::shared_ptr< Epetra_FEVector > vector
const Epetra_Comm & comm_self()
Definition: utilities.cc:740
void reinit(const VectorBase &v, const bool omit_zeroing_entries=false, const bool allow_different_maps=false)
#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