Reference documentation for deal.II version Git 0183d66 2014-10-17 10:54:07 +0200
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
trilinos_vector.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2008 - 2014 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 # include "Epetra_Map.h"
33 # include "Epetra_LocalMap.h"
34 
35 DEAL_II_NAMESPACE_OPEN
36 
37 
38 // forward declaration
39 template <typename> class Vector;
40 
45 namespace TrilinosWrappers
46 {
47  class SparseMatrix;
48 
49  namespace
50  {
51 #ifndef DEAL_II_WITH_64BIT_INDICES
52  // define a helper function that queries the global ID of local ID of
53  // an Epetra_BlockMap object by calling either the 32- or 64-bit
54  // function necessary.
55  int gid(const Epetra_BlockMap &map, int i)
56  {
57  return map.GID(i);
58  }
59 #else
60  // define a helper function that queries the global ID of local ID of
61  // an Epetra_BlockMap object by calling either the 32- or 64-bit
62  // function necessary.
63  long long int gid(const Epetra_BlockMap &map, int i)
64  {
65  return map.GID64(i);
66  }
67 #endif
68  }
69 
78  namespace MPI
79  {
80  class BlockVector;
81 
250  class Vector : public VectorBase
251  {
252  public:
256  typedef ::types::global_dof_index size_type;
257 
268  static const bool supports_distributed_data = true;
269 
279  Vector ();
280 
284  Vector (const Vector &V);
285 
289  ~Vector ();
290 
310  void reinit (const VectorBase &v,
311  const bool fast = false,
312  const bool allow_different_maps = false);
313 
317  void reinit (const BlockVector &v,
318  const bool import_data = false);
319 
326  Vector &operator = (const TrilinosScalar s);
327 
333  Vector &
334  operator = (const Vector &V);
335 
343  Vector &
344  operator = (const ::TrilinosWrappers::Vector &V);
345 
354  template <typename Number>
355  Vector &
356  operator = (const ::Vector<Number> &v);
357 
376  (const ::TrilinosWrappers::SparseMatrix &matrix,
377  const Vector &vector);
379 
396  explicit Vector (const Epetra_Map &parallel_partitioning);
397 
411  Vector (const Epetra_Map &parallel_partitioning,
412  const VectorBase &v);
413 
425  template <typename number>
426  void reinit (const Epetra_Map &parallel_partitioner,
427  const ::Vector<number> &v);
428 
440  void reinit (const Epetra_Map &parallel_partitioning,
441  const bool fast = false);
442 
454  template <typename Number>
455  Vector (const Epetra_Map &parallel_partitioning,
456  const ::Vector<Number> &v);
458 
475  explicit Vector (const IndexSet &parallel_partitioning,
476  const MPI_Comm &communicator = MPI_COMM_WORLD);
477 
488  Vector (const IndexSet &local,
489  const IndexSet &ghost,
490  const MPI_Comm &communicator = MPI_COMM_WORLD);
491 
505  Vector (const IndexSet &parallel_partitioning,
506  const VectorBase &v,
507  const MPI_Comm &communicator = MPI_COMM_WORLD);
508 
520  template <typename Number>
521  Vector (const IndexSet &parallel_partitioning,
522  const ::Vector<Number> &v,
523  const MPI_Comm &communicator = MPI_COMM_WORLD);
524 
539  void reinit (const IndexSet &parallel_partitioning,
540  const MPI_Comm &communicator = MPI_COMM_WORLD,
541  const bool fast = false);
542 
568  void reinit (const IndexSet &locally_owned_entries,
569  const IndexSet &ghost_entries,
570  const MPI_Comm &communicator = MPI_COMM_WORLD,
571  const bool vector_writable = false);
573  };
574 
575 
576 
577 
578 // ------------------- inline and template functions --------------
579 
580 
589  inline
590  void swap (Vector &u, Vector &v)
591  {
592  u.swap (v);
593  }
594 
595 
596 #ifndef DOXYGEN
597 
598  template <typename number>
599  Vector::Vector (const Epetra_Map &input_map,
600  const ::Vector<number> &v)
601  {
602  reinit (input_map, v);
603  }
604 
605 
606 
607  template <typename number>
608  Vector::Vector (const IndexSet &parallel_partitioner,
609  const ::Vector<number> &v,
610  const MPI_Comm &communicator)
611  {
612  *this = Vector(parallel_partitioner.make_trilinos_map (communicator, true),
613  v);
614  }
615 
616 
617 
618 
619  template <typename number>
620  void Vector::reinit (const Epetra_Map &parallel_partitioner,
621  const ::Vector<number> &v)
622  {
623  if (vector.get() == 0 || vector->Map().SameAs(parallel_partitioner) == false)
624  vector.reset (new Epetra_FEVector(parallel_partitioner));
625 
626  has_ghosts = vector->Map().UniqueGIDs()==false;
627 
628  const int size = parallel_partitioner.NumMyElements();
629 
630  // Need to copy out values, since the deal.II might not use doubles, so
631  // that a direct access is not possible.
632  for (int i=0; i<size; ++i)
633  (*vector)[0][i] = v(gid(parallel_partitioner,i));
634  }
635 
636 
637  inline
638  Vector &
639  Vector::operator = (const TrilinosScalar s)
640  {
642 
643  return *this;
644  }
645 
646 
647  template <typename Number>
648  Vector &
649  Vector::operator = (const ::Vector<Number> &v)
650  {
651  if (size() != v.size())
652  {
653  vector.reset (new Epetra_FEVector(Epetra_Map
654  (static_cast<TrilinosWrappers::types::int_type>(v.size()), 0,
655 #ifdef DEAL_II_WITH_MPI
656  Epetra_MpiComm(MPI_COMM_SELF)
657 #else
658  Epetra_SerialComm()
659 #endif
660  )));
661  }
662 
663  reinit (vector_partitioner(), v);
664  return *this;
665  }
666 
667 
668 #endif
669 
670  } /* end of namespace MPI */
671 
672 
673 
685  class Vector : public VectorBase
686  {
687  public:
691  typedef ::types::global_dof_index size_type;
692 
703  static const bool supports_distributed_data = false;
704 
710  Vector ();
711 
715  explicit Vector (const size_type n);
716 
726  explicit Vector (const Epetra_Map &partitioning);
727 
737  explicit Vector (const IndexSet &partitioning,
738  const MPI_Comm &communicator = MPI_COMM_WORLD);
739 
744  explicit Vector (const VectorBase &V);
745 
750  template <typename Number>
751  explicit Vector (const ::Vector<Number> &v);
752 
757  void reinit (const size_type n,
758  const bool fast = false);
759 
773  void reinit (const Epetra_Map &input_map,
774  const bool fast = false);
775 
789  void reinit (const IndexSet &input_map,
790  const MPI_Comm &communicator = MPI_COMM_WORLD,
791  const bool fast = false);
792 
797  void reinit (const VectorBase &V,
798  const bool fast = false,
799  const bool allow_different_maps = false);
800 
807  Vector &operator = (const TrilinosScalar s);
808 
813  Vector &
814  operator = (const MPI::Vector &V);
815 
819  template <typename Number>
820  Vector &
821  operator = (const ::Vector<Number> &V);
822 
827  Vector &
828  operator = (const Vector &V);
829 
841  void update_ghost_values () const;
842  };
843 
844 
845 
846 // ------------------- inline and template functions --------------
847 
848 
857  inline
858  void swap (Vector &u, Vector &v)
859  {
860  u.swap (v);
861  }
862 
863 
864 #ifndef DOXYGEN
865 
866  template <typename number>
867  Vector::Vector (const ::Vector<number> &v)
868  {
869  Epetra_LocalMap map ((TrilinosWrappers::types::int_type)v.size(), 0, Utilities::Trilinos::comm_self());
870  vector.reset (new Epetra_FEVector(map));
871  *this = v;
872  }
873 
874 
875 
876  inline
877  Vector &
878  Vector::operator = (const TrilinosScalar s)
879  {
881 
882  return *this;
883  }
884 
885 
886 
887  template <typename Number>
888  Vector &
889  Vector::operator = (const ::Vector<Number> &v)
890  {
891  if (size() != v.size())
892  {
893  vector.reset();
894 
895  Epetra_LocalMap map ((TrilinosWrappers::types::int_type)v.size(), 0,
897  vector.reset (new Epetra_FEVector(map));
898  }
899 
900  const Epetra_Map &map = vector_partitioner();
901  const TrilinosWrappers::types::int_type size = map.NumMyElements();
902 
903  Assert (map.MaxLID() == size-1,
904  ExcDimensionMismatch(map.MaxLID(), size-1));
905 
906  // Need to copy out values, since the
907  // deal.II might not use doubles, so
908  // that a direct access is not possible.
909  for (TrilinosWrappers::types::int_type i=0; i<size; ++i)
910  (*vector)[0][i] = v(i);
911 
912  return *this;
913  }
914 
915 
916 
917  inline
918  void
920  {}
921 
922 
923 #endif
924 
925 
926 }
927 
928 
931 DEAL_II_NAMESPACE_CLOSE
932 
933 #endif // DEAL_II_WITH_TRILINOS
934 
935 /*---------------------------- trilinos_vector.h ---------------------------*/
936 
937 #endif
938 /*---------------------------- 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:298
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)
::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
void update_ghost_values() const
::types::global_dof_index size_type
size_type size() const
static const bool supports_distributed_data