Reference documentation for deal.II version Git 437fb4a 2015-08-02 17:50:50 -0500
trilinos_vector.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2008 - 2015 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 
324  void reinit (const VectorBase &v,
325  const bool fast = false,
326  const bool allow_different_maps = false);
327 
331  void reinit (const BlockVector &v,
332  const bool import_data = false);
333 
340  Vector &operator= (const TrilinosScalar s);
341 
347  Vector &operator= (const Vector &v);
348 
349 #ifdef DEAL_II_WITH_CXX11
350 
357  Vector &operator= (Vector &&v);
358 #endif
359 
367  Vector &operator= (const ::TrilinosWrappers::Vector &v);
368 
377  template <typename Number>
378  Vector &operator= (const ::Vector<Number> &v);
379 
398  (const ::TrilinosWrappers::SparseMatrix &matrix,
399  const Vector &vector);
401 
421  explicit Vector (const Epetra_Map &parallel_partitioning) DEAL_II_DEPRECATED;
422 
439  Vector (const Epetra_Map &parallel_partitioning,
440  const VectorBase &v) DEAL_II_DEPRECATED;
441 
456  template <typename number>
457  void reinit (const Epetra_Map &parallel_partitioner,
458  const ::Vector<number> &v) DEAL_II_DEPRECATED;
459 
474  void reinit (const Epetra_Map &parallel_partitioning,
475  const bool fast = false) DEAL_II_DEPRECATED;
476 
491  template <typename Number>
492  Vector (const Epetra_Map &parallel_partitioning,
493  const ::Vector<Number> &v) DEAL_II_DEPRECATED;
495 
513  explicit Vector (const IndexSet &parallel_partitioning,
514  const MPI_Comm &communicator = MPI_COMM_WORLD);
515 
527  Vector (const IndexSet &local,
528  const IndexSet &ghost,
529  const MPI_Comm &communicator = MPI_COMM_WORLD);
530 
545  Vector (const IndexSet &parallel_partitioning,
546  const VectorBase &v,
547  const MPI_Comm &communicator = MPI_COMM_WORLD);
548 
561  template <typename Number>
562  Vector (const IndexSet &parallel_partitioning,
563  const ::Vector<Number> &v,
564  const MPI_Comm &communicator = MPI_COMM_WORLD);
565 
581  void reinit (const IndexSet &parallel_partitioning,
582  const MPI_Comm &communicator = MPI_COMM_WORLD,
583  const bool fast = false);
584 
610  void reinit (const IndexSet &locally_owned_entries,
611  const IndexSet &ghost_entries,
612  const MPI_Comm &communicator = MPI_COMM_WORLD,
613  const bool vector_writable = false);
615  };
616 
617 
618 
619 
620 // ------------------- inline and template functions --------------
621 
622 
631  inline
632  void swap (Vector &u, Vector &v)
633  {
634  u.swap (v);
635  }
636 
637 
638 #ifndef DOXYGEN
639 
640  template <typename number>
641  Vector::Vector (const Epetra_Map &input_map,
642  const ::Vector<number> &v)
643  {
644  reinit (input_map, v);
645  }
646 
647 
648 
649  template <typename number>
650  Vector::Vector (const IndexSet &parallel_partitioner,
651  const ::Vector<number> &v,
652  const MPI_Comm &communicator)
653  {
654  *this = Vector(parallel_partitioner.make_trilinos_map (communicator, true),
655  v);
656  }
657 
658 
659 
660 
661  template <typename number>
662  void Vector::reinit (const Epetra_Map &parallel_partitioner,
663  const ::Vector<number> &v)
664  {
665  if (vector.get() == 0 || vector->Map().SameAs(parallel_partitioner) == false)
666  vector.reset (new Epetra_FEVector(parallel_partitioner));
667 
668  has_ghosts = vector->Map().UniqueGIDs()==false;
669 
670  const int size = parallel_partitioner.NumMyElements();
671 
672  // Need to copy out values, since the deal.II might not use doubles, so
673  // that a direct access is not possible.
674  for (int i=0; i<size; ++i)
675  (*vector)[0][i] = v(gid(parallel_partitioner,i));
676  }
677 
678 
679  inline
680  Vector &
681  Vector::operator= (const TrilinosScalar s)
682  {
684 
685  return *this;
686  }
687 
688 
689  template <typename Number>
690  Vector &
691  Vector::operator= (const ::Vector<Number> &v)
692  {
693  if (size() != v.size())
694  {
695  vector.reset (new Epetra_FEVector(Epetra_Map
696  (static_cast<TrilinosWrappers::types::int_type>(v.size()), 0,
697 #ifdef DEAL_II_WITH_MPI
698  Epetra_MpiComm(MPI_COMM_SELF)
699 #else
700  Epetra_SerialComm()
701 #endif
702  )));
703  }
704 
705  reinit (vector_partitioner(), v);
706  return *this;
707  }
708 
709 
710 #endif /* DOXYGEN */
711 
712  } /* end of namespace MPI */
713 
714 
715 
729  class Vector : public VectorBase
730  {
731  public:
735  typedef ::types::global_dof_index size_type;
736 
747  static const bool supports_distributed_data = false;
748 
754  Vector () DEAL_II_DEPRECATED;
755 
759  explicit Vector (const size_type n) DEAL_II_DEPRECATED;
760 
770  explicit Vector (const Epetra_Map &partitioning) DEAL_II_DEPRECATED;
771 
781  explicit Vector (const IndexSet &partitioning,
782  const MPI_Comm &communicator = MPI_COMM_WORLD) DEAL_II_DEPRECATED;
783 
788  explicit Vector (const VectorBase &V) DEAL_II_DEPRECATED;
789 
794  template <typename Number>
795  explicit Vector (const ::Vector<Number> &v) DEAL_II_DEPRECATED;
796 
801  void reinit (const size_type n,
802  const bool fast = false);
803 
817  void reinit (const Epetra_Map &input_map,
818  const bool fast = false);
819 
833  void reinit (const IndexSet &input_map,
834  const MPI_Comm &communicator = MPI_COMM_WORLD,
835  const bool fast = false);
836 
841  void reinit (const VectorBase &V,
842  const bool fast = false,
843  const bool allow_different_maps = false);
844 
851  Vector &operator= (const TrilinosScalar s);
852 
857  Vector &operator= (const MPI::Vector &v);
858 
862  template <typename Number>
863  Vector &operator= (const ::Vector<Number> &v);
864 
869  Vector &operator= (const Vector &v);
870 
882  void update_ghost_values () const;
883  };
884 
885 
886 
887 // ------------------- inline and template functions --------------
888 
889 
898  inline
899  void swap (Vector &u, Vector &v)
900  {
901  u.swap (v);
902  }
903 
904 
905 #ifndef DOXYGEN
906 
907  template <typename number>
908  Vector::Vector (const ::Vector<number> &v)
909  {
910  Epetra_LocalMap map ((TrilinosWrappers::types::int_type)v.size(), 0, Utilities::Trilinos::comm_self());
911  vector.reset (new Epetra_FEVector(map));
912  *this = v;
913  }
914 
915 
916 
917  inline
918  Vector &
919  Vector::operator= (const TrilinosScalar s)
920  {
922 
923  return *this;
924  }
925 
926 
927 
928  template <typename Number>
929  Vector &
930  Vector::operator= (const ::Vector<Number> &v)
931  {
932  if (size() != v.size())
933  {
934  vector.reset();
935 
936  Epetra_LocalMap map ((TrilinosWrappers::types::int_type)v.size(), 0,
938  vector.reset (new Epetra_FEVector(map));
939  }
940 
941  const Epetra_Map &map = vector_partitioner();
942  const TrilinosWrappers::types::int_type size = map.NumMyElements();
943 
944  Assert (map.MaxLID() == size-1,
945  ExcDimensionMismatch(map.MaxLID(), size-1));
946 
947  // Need to copy out values, since the
948  // deal.II might not use doubles, so
949  // that a direct access is not possible.
950  for (TrilinosWrappers::types::int_type i=0; i<size; ++i)
951  (*vector)[0][i] = v(i);
952 
953  return *this;
954  }
955 
956 
957 
958  inline
959  void
961  {}
962 
963 
964 #endif
965 
966 } /* namespace TrilinosWrappers */
967 
971 namespace internal
972 {
973  namespace LinearOperator
974  {
975  template <typename> class ReinitHelper;
976 
981  template<>
982  class ReinitHelper<TrilinosWrappers::MPI::Vector>
983  {
984  public:
985  template <typename Matrix>
986  static
987  void reinit_range_vector (const Matrix &matrix,
989  bool fast)
990  {
991  v.reinit(matrix.locally_owned_range_indices(), matrix.get_mpi_communicator(), fast);
992  }
993 
994  template <typename Matrix>
995  static
996  void reinit_domain_vector(const Matrix &matrix,
998  bool fast)
999  {
1000  v.reinit(matrix.locally_owned_domain_indices(), matrix.get_mpi_communicator(), fast);
1001  }
1002  };
1003 
1008  template<>
1010  {
1011  public:
1012  template <typename Matrix>
1013  static
1014  void reinit_range_vector (const Matrix &matrix,
1016  bool fast)
1017  {
1018  v.reinit(matrix.locally_owned_range_indices(), matrix.get_mpi_communicator(), fast);
1019  }
1020 
1021  template <typename Matrix>
1022  static
1023  void reinit_domain_vector(const Matrix &matrix,
1025  bool fast)
1026  {
1027  v.reinit(matrix.locally_owned_domain_indices(), matrix.get_mpi_communicator(), fast);
1028  }
1029  };
1030 
1031  } /* namespace LinearOperator */
1032 } /* namespace internal */
1033 
1034 
1035 DEAL_II_NAMESPACE_CLOSE
1036 
1037 #endif // DEAL_II_WITH_TRILINOS
1038 
1039 /*---------------------------- trilinos_vector.h ---------------------------*/
1040 
1041 #endif
1042 /*---------------------------- trilinos_vector.h ---------------------------*/
Vector & operator=(const TrilinosScalar s)
void reinit(const size_type n, const bool fast=false)
::types::global_dof_index size_type
void import_nonlocal_data_for_fe(const ::TrilinosWrappers::SparseMatrix &matrix, const Vector &vector)
#define DEAL_II_DEPRECATED
Definition: config.h:88
std_cxx11::shared_ptr< Epetra_FEVector > vector
#define Assert(cond, exc)
Definition: exceptions.h:293
const Epetra_Comm & comm_self()
void swap(Vector &u, Vector &v)
const Epetra_Map & vector_partitioner() const
Epetra_Map make_trilinos_map(const MPI_Comm &communicator=MPI_COMM_WORLD, const bool overlapping=false) const
VectorBase & operator=(const TrilinosScalar s)
void reinit(const VectorBase &v, const bool fast=false, const bool allow_different_maps=false)
void update_ghost_values() const
::types::global_dof_index size_type
size_type size() const
static const bool supports_distributed_data