Reference documentation for deal.II version Git ca90e6c 2015-05-21 13:56:49 -0400
 All Classes Namespaces Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
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 __deal2__trilinos_vector_h
17 #define __deal2__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 
419  explicit Vector (const Epetra_Map &parallel_partitioning);
420 
435  Vector (const Epetra_Map &parallel_partitioning,
436  const VectorBase &v);
437 
450  template <typename number>
451  void reinit (const Epetra_Map &parallel_partitioner,
452  const ::Vector<number> &v);
453 
466  void reinit (const Epetra_Map &parallel_partitioning,
467  const bool fast = false);
468 
481  template <typename Number>
482  Vector (const Epetra_Map &parallel_partitioning,
483  const ::Vector<Number> &v);
485 
503  explicit Vector (const IndexSet &parallel_partitioning,
504  const MPI_Comm &communicator = MPI_COMM_WORLD);
505 
517  Vector (const IndexSet &local,
518  const IndexSet &ghost,
519  const MPI_Comm &communicator = MPI_COMM_WORLD);
520 
535  Vector (const IndexSet &parallel_partitioning,
536  const VectorBase &v,
537  const MPI_Comm &communicator = MPI_COMM_WORLD);
538 
551  template <typename Number>
552  Vector (const IndexSet &parallel_partitioning,
553  const ::Vector<Number> &v,
554  const MPI_Comm &communicator = MPI_COMM_WORLD);
555 
571  void reinit (const IndexSet &parallel_partitioning,
572  const MPI_Comm &communicator = MPI_COMM_WORLD,
573  const bool fast = false);
574 
600  void reinit (const IndexSet &locally_owned_entries,
601  const IndexSet &ghost_entries,
602  const MPI_Comm &communicator = MPI_COMM_WORLD,
603  const bool vector_writable = false);
605  };
606 
607 
608 
609 
610 // ------------------- inline and template functions --------------
611 
612 
621  inline
622  void swap (Vector &u, Vector &v)
623  {
624  u.swap (v);
625  }
626 
627 
628 #ifndef DOXYGEN
629 
630  template <typename number>
631  Vector::Vector (const Epetra_Map &input_map,
632  const ::Vector<number> &v)
633  {
634  reinit (input_map, v);
635  }
636 
637 
638 
639  template <typename number>
640  Vector::Vector (const IndexSet &parallel_partitioner,
641  const ::Vector<number> &v,
642  const MPI_Comm &communicator)
643  {
644  *this = Vector(parallel_partitioner.make_trilinos_map (communicator, true),
645  v);
646  }
647 
648 
649 
650 
651  template <typename number>
652  void Vector::reinit (const Epetra_Map &parallel_partitioner,
653  const ::Vector<number> &v)
654  {
655  if (vector.get() == 0 || vector->Map().SameAs(parallel_partitioner) == false)
656  vector.reset (new Epetra_FEVector(parallel_partitioner));
657 
658  has_ghosts = vector->Map().UniqueGIDs()==false;
659 
660  const int size = parallel_partitioner.NumMyElements();
661 
662  // Need to copy out values, since the deal.II might not use doubles, so
663  // that a direct access is not possible.
664  for (int i=0; i<size; ++i)
665  (*vector)[0][i] = v(gid(parallel_partitioner,i));
666  }
667 
668 
669  inline
670  Vector &
671  Vector::operator= (const TrilinosScalar s)
672  {
674 
675  return *this;
676  }
677 
678 
679  template <typename Number>
680  Vector &
681  Vector::operator= (const ::Vector<Number> &v)
682  {
683  if (size() != v.size())
684  {
685  vector.reset (new Epetra_FEVector(Epetra_Map
686  (static_cast<TrilinosWrappers::types::int_type>(v.size()), 0,
687 #ifdef DEAL_II_WITH_MPI
688  Epetra_MpiComm(MPI_COMM_SELF)
689 #else
690  Epetra_SerialComm()
691 #endif
692  )));
693  }
694 
695  reinit (vector_partitioner(), v);
696  return *this;
697  }
698 
699 
700 #endif /* DOXYGEN */
701 
702  } /* end of namespace MPI */
703 
704 
705 
719  class Vector : public VectorBase
720  {
721  public:
725  typedef ::types::global_dof_index size_type;
726 
737  static const bool supports_distributed_data = false;
738 
745 
749  explicit Vector (const size_type n) DEAL_II_DEPRECATED;
750 
760  explicit Vector (const Epetra_Map &partitioning) DEAL_II_DEPRECATED;
761 
771  explicit Vector (const IndexSet &partitioning,
772  const MPI_Comm &communicator = MPI_COMM_WORLD) DEAL_II_DEPRECATED;
773 
778  explicit Vector (const VectorBase &V) DEAL_II_DEPRECATED;
779 
784  template <typename Number>
785  explicit Vector (const ::Vector<Number> &v) DEAL_II_DEPRECATED;
786 
791  void reinit (const size_type n,
792  const bool fast = false);
793 
807  void reinit (const Epetra_Map &input_map,
808  const bool fast = false);
809 
823  void reinit (const IndexSet &input_map,
824  const MPI_Comm &communicator = MPI_COMM_WORLD,
825  const bool fast = false);
826 
831  void reinit (const VectorBase &V,
832  const bool fast = false,
833  const bool allow_different_maps = false);
834 
841  Vector &operator= (const TrilinosScalar s);
842 
847  Vector &operator= (const MPI::Vector &v);
848 
852  template <typename Number>
853  Vector &operator= (const ::Vector<Number> &v);
854 
859  Vector &operator= (const Vector &v);
860 
872  void update_ghost_values () const;
873  };
874 
875 
876 
877 // ------------------- inline and template functions --------------
878 
879 
888  inline
889  void swap (Vector &u, Vector &v)
890  {
891  u.swap (v);
892  }
893 
894 
895 #ifndef DOXYGEN
896 
897  template <typename number>
898  Vector::Vector (const ::Vector<number> &v)
899  {
900  Epetra_LocalMap map ((TrilinosWrappers::types::int_type)v.size(), 0, Utilities::Trilinos::comm_self());
901  vector.reset (new Epetra_FEVector(map));
902  *this = v;
903  }
904 
905 
906 
907  inline
908  Vector &
909  Vector::operator= (const TrilinosScalar s)
910  {
912 
913  return *this;
914  }
915 
916 
917 
918  template <typename Number>
919  Vector &
920  Vector::operator= (const ::Vector<Number> &v)
921  {
922  if (size() != v.size())
923  {
924  vector.reset();
925 
926  Epetra_LocalMap map ((TrilinosWrappers::types::int_type)v.size(), 0,
928  vector.reset (new Epetra_FEVector(map));
929  }
930 
931  const Epetra_Map &map = vector_partitioner();
932  const TrilinosWrappers::types::int_type size = map.NumMyElements();
933 
934  Assert (map.MaxLID() == size-1,
935  ExcDimensionMismatch(map.MaxLID(), size-1));
936 
937  // Need to copy out values, since the
938  // deal.II might not use doubles, so
939  // that a direct access is not possible.
940  for (TrilinosWrappers::types::int_type i=0; i<size; ++i)
941  (*vector)[0][i] = v(i);
942 
943  return *this;
944  }
945 
946 
947 
948  inline
949  void
951  {}
952 
953 
954 #endif
955 
956 } /* namespace TrilinosWrappers */
957 
961 namespace internal
962 {
963  namespace LinearOperator
964  {
965  template <typename> class ReinitHelper;
966 
971  template<>
972  class ReinitHelper<TrilinosWrappers::MPI::Vector>
973  {
974  public:
975  template <typename Matrix>
976  static
977  void reinit_range_vector (const Matrix &matrix,
979  bool fast)
980  {
981  v.reinit(matrix.range_partitioner(), fast);
982  }
983 
984  template <typename Matrix>
985  static
986  void reinit_domain_vector(const Matrix &matrix,
988  bool fast)
989  {
990  v.reinit(matrix.domain_partitioner(), fast);
991  }
992  };
993 
998  template<>
999  class ReinitHelper<TrilinosWrappers::Vector>
1000  {
1001  public:
1002  template <typename Matrix>
1003  static
1004  void reinit_range_vector (const Matrix &matrix,
1006  bool fast)
1007  {
1008  v.reinit(matrix.range_partitioner(), fast);
1009  }
1010 
1011  template <typename Matrix>
1012  static
1013  void reinit_domain_vector(const Matrix &matrix,
1015  bool fast)
1016  {
1017  v.reinit(matrix.domain_partitioner(), fast);
1018  }
1019  };
1020 
1021  } /* namespace LinearOperator */
1022 } /* namespace internal */
1023 
1024 
1025 DEAL_II_NAMESPACE_CLOSE
1026 
1027 #endif // DEAL_II_WITH_TRILINOS
1028 
1029 /*---------------------------- trilinos_vector.h ---------------------------*/
1030 
1031 #endif
1032 /*---------------------------- trilinos_vector.h ---------------------------*/
Vector & operator=(const TrilinosScalar s)
static const bool supports_distributed_data
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)
Vector & operator=(const TrilinosScalar s)
void swap(Vector &u, Vector &v)
BlockDynamicSparsityPattern BlockCompressedSparsityPattern DEAL_II_DEPRECATED
Vector() DEAL_II_DEPRECATED
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