Reference documentation for deal.II version Git 672aa0c 2015-04-25 06:51:17 -0500
 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 
292  ~Vector ();
293 
313  void reinit (const VectorBase &v,
314  const bool fast = false,
315  const bool allow_different_maps = false);
316 
320  void reinit (const BlockVector &v,
321  const bool import_data = false);
322 
329  Vector &operator = (const TrilinosScalar s);
330 
336  Vector &
337  operator = (const Vector &V);
338 
346  Vector &
347  operator = (const ::TrilinosWrappers::Vector &V);
348 
357  template <typename Number>
358  Vector &
359  operator = (const ::Vector<Number> &v);
360 
379  (const ::TrilinosWrappers::SparseMatrix &matrix,
380  const Vector &vector);
382 
400  explicit Vector (const Epetra_Map &parallel_partitioning);
401 
416  Vector (const Epetra_Map &parallel_partitioning,
417  const VectorBase &v);
418 
431  template <typename number>
432  void reinit (const Epetra_Map &parallel_partitioner,
433  const ::Vector<number> &v);
434 
447  void reinit (const Epetra_Map &parallel_partitioning,
448  const bool fast = false);
449 
462  template <typename Number>
463  Vector (const Epetra_Map &parallel_partitioning,
464  const ::Vector<Number> &v);
466 
484  explicit Vector (const IndexSet &parallel_partitioning,
485  const MPI_Comm &communicator = MPI_COMM_WORLD);
486 
498  Vector (const IndexSet &local,
499  const IndexSet &ghost,
500  const MPI_Comm &communicator = MPI_COMM_WORLD);
501 
516  Vector (const IndexSet &parallel_partitioning,
517  const VectorBase &v,
518  const MPI_Comm &communicator = MPI_COMM_WORLD);
519 
532  template <typename Number>
533  Vector (const IndexSet &parallel_partitioning,
534  const ::Vector<Number> &v,
535  const MPI_Comm &communicator = MPI_COMM_WORLD);
536 
552  void reinit (const IndexSet &parallel_partitioning,
553  const MPI_Comm &communicator = MPI_COMM_WORLD,
554  const bool fast = false);
555 
581  void reinit (const IndexSet &locally_owned_entries,
582  const IndexSet &ghost_entries,
583  const MPI_Comm &communicator = MPI_COMM_WORLD,
584  const bool vector_writable = false);
586  };
587 
588 
589 
590 
591 // ------------------- inline and template functions --------------
592 
593 
602  inline
603  void swap (Vector &u, Vector &v)
604  {
605  u.swap (v);
606  }
607 
608 
609 #ifndef DOXYGEN
610 
611  template <typename number>
612  Vector::Vector (const Epetra_Map &input_map,
613  const ::Vector<number> &v)
614  {
615  reinit (input_map, v);
616  }
617 
618 
619 
620  template <typename number>
621  Vector::Vector (const IndexSet &parallel_partitioner,
622  const ::Vector<number> &v,
623  const MPI_Comm &communicator)
624  {
625  *this = Vector(parallel_partitioner.make_trilinos_map (communicator, true),
626  v);
627  }
628 
629 
630 
631 
632  template <typename number>
633  void Vector::reinit (const Epetra_Map &parallel_partitioner,
634  const ::Vector<number> &v)
635  {
636  if (vector.get() == 0 || vector->Map().SameAs(parallel_partitioner) == false)
637  vector.reset (new Epetra_FEVector(parallel_partitioner));
638 
639  has_ghosts = vector->Map().UniqueGIDs()==false;
640 
641  const int size = parallel_partitioner.NumMyElements();
642 
643  // Need to copy out values, since the deal.II might not use doubles, so
644  // that a direct access is not possible.
645  for (int i=0; i<size; ++i)
646  (*vector)[0][i] = v(gid(parallel_partitioner,i));
647  }
648 
649 
650  inline
651  Vector &
652  Vector::operator = (const TrilinosScalar s)
653  {
655 
656  return *this;
657  }
658 
659 
660  template <typename Number>
661  Vector &
662  Vector::operator = (const ::Vector<Number> &v)
663  {
664  if (size() != v.size())
665  {
666  vector.reset (new Epetra_FEVector(Epetra_Map
667  (static_cast<TrilinosWrappers::types::int_type>(v.size()), 0,
668 #ifdef DEAL_II_WITH_MPI
669  Epetra_MpiComm(MPI_COMM_SELF)
670 #else
671  Epetra_SerialComm()
672 #endif
673  )));
674  }
675 
676  reinit (vector_partitioner(), v);
677  return *this;
678  }
679 
680 
681 #endif
682 
683  } /* end of namespace MPI */
684 
685 
686 
698  class Vector : public VectorBase
699  {
700  public:
704  typedef ::types::global_dof_index size_type;
705 
716  static const bool supports_distributed_data = false;
717 
723  Vector ();
724 
728  explicit Vector (const size_type n);
729 
739  explicit Vector (const Epetra_Map &partitioning);
740 
750  explicit Vector (const IndexSet &partitioning,
751  const MPI_Comm &communicator = MPI_COMM_WORLD);
752 
757  explicit Vector (const VectorBase &V);
758 
763  template <typename Number>
764  explicit Vector (const ::Vector<Number> &v);
765 
770  void reinit (const size_type n,
771  const bool fast = false);
772 
786  void reinit (const Epetra_Map &input_map,
787  const bool fast = false);
788 
802  void reinit (const IndexSet &input_map,
803  const MPI_Comm &communicator = MPI_COMM_WORLD,
804  const bool fast = false);
805 
810  void reinit (const VectorBase &V,
811  const bool fast = false,
812  const bool allow_different_maps = false);
813 
820  Vector &operator = (const TrilinosScalar s);
821 
826  Vector &
827  operator = (const MPI::Vector &V);
828 
832  template <typename Number>
833  Vector &
834  operator = (const ::Vector<Number> &V);
835 
840  Vector &
841  operator = (const Vector &V);
842 
854  void update_ghost_values () const;
855  };
856 
857 
858 
859 // ------------------- inline and template functions --------------
860 
861 
870  inline
871  void swap (Vector &u, Vector &v)
872  {
873  u.swap (v);
874  }
875 
876 
877 #ifndef DOXYGEN
878 
879  template <typename number>
880  Vector::Vector (const ::Vector<number> &v)
881  {
882  Epetra_LocalMap map ((TrilinosWrappers::types::int_type)v.size(), 0, Utilities::Trilinos::comm_self());
883  vector.reset (new Epetra_FEVector(map));
884  *this = v;
885  }
886 
887 
888 
889  inline
890  Vector &
891  Vector::operator = (const TrilinosScalar s)
892  {
894 
895  return *this;
896  }
897 
898 
899 
900  template <typename Number>
901  Vector &
902  Vector::operator = (const ::Vector<Number> &v)
903  {
904  if (size() != v.size())
905  {
906  vector.reset();
907 
908  Epetra_LocalMap map ((TrilinosWrappers::types::int_type)v.size(), 0,
910  vector.reset (new Epetra_FEVector(map));
911  }
912 
913  const Epetra_Map &map = vector_partitioner();
914  const TrilinosWrappers::types::int_type size = map.NumMyElements();
915 
916  Assert (map.MaxLID() == size-1,
917  ExcDimensionMismatch(map.MaxLID(), size-1));
918 
919  // Need to copy out values, since the
920  // deal.II might not use doubles, so
921  // that a direct access is not possible.
922  for (TrilinosWrappers::types::int_type i=0; i<size; ++i)
923  (*vector)[0][i] = v(i);
924 
925  return *this;
926  }
927 
928 
929 
930  inline
931  void
933  {}
934 
935 
936 #endif
937 
938 } /* namespace TrilinosWrappers */
939 
943 namespace internal
944 {
945  namespace LinearOperator
946  {
947  template <typename> class ReinitHelper;
948 
953  template<>
954  class ReinitHelper<TrilinosWrappers::MPI::Vector>
955  {
956  public:
957  template <typename Matrix>
958  static
959  void reinit_range_vector (const Matrix &matrix,
961  bool fast)
962  {
963  v.reinit(matrix.range_partitioner(), fast);
964  }
965 
966  template <typename Matrix>
967  static
968  void reinit_domain_vector(const Matrix &matrix,
970  bool fast)
971  {
972  v.reinit(matrix.domain_partitioner(), fast);
973  }
974  };
975 
980  template<>
981  class ReinitHelper<TrilinosWrappers::Vector>
982  {
983  public:
984  template <typename Matrix>
985  static
986  void reinit_range_vector (const Matrix &matrix,
988  bool fast)
989  {
990  v.reinit(matrix.range_partitioner(), fast);
991  }
992 
993  template <typename Matrix>
994  static
995  void reinit_domain_vector(const Matrix &matrix,
997  bool fast)
998  {
999  v.reinit(matrix.domain_partitioner(), fast);
1000  }
1001  };
1002 
1003  } /* namespace LinearOperator */
1004 } /* namespace internal */
1005 
1006 
1007 DEAL_II_NAMESPACE_CLOSE
1008 
1009 #endif // DEAL_II_WITH_TRILINOS
1010 
1011 /*---------------------------- trilinos_vector.h ---------------------------*/
1012 
1013 #endif
1014 /*---------------------------- 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)
std_cxx11::shared_ptr< Epetra_FEVector > vector
#define Assert(cond, exc)
Definition: exceptions.h:293
static void reinit_range_vector(const Matrix &matrix, Vector &v, bool fast)
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)
static void reinit_domain_vector(const Matrix &matrix, Vector &v, bool fast)
void update_ghost_values() const
::types::global_dof_index size_type
size_type size() const
static const bool supports_distributed_data