Reference documentation for deal.II version Git 06e7726 2015-03-27 10:02:22 +0100
 All Classes Namespaces 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 
249  class Vector : public VectorBase
250  {
251  public:
255  typedef ::types::global_dof_index size_type;
256 
267  static const bool supports_distributed_data = true;
268 
278  Vector ();
279 
283  Vector (const Vector &V);
284 
288  ~Vector ();
289 
309  void reinit (const VectorBase &v,
310  const bool fast = false,
311  const bool allow_different_maps = false);
312 
316  void reinit (const BlockVector &v,
317  const bool import_data = false);
318 
325  Vector &operator = (const TrilinosScalar s);
326 
332  Vector &
333  operator = (const Vector &V);
334 
342  Vector &
343  operator = (const ::TrilinosWrappers::Vector &V);
344 
353  template <typename Number>
354  Vector &
355  operator = (const ::Vector<Number> &v);
356 
375  (const ::TrilinosWrappers::SparseMatrix &matrix,
376  const Vector &vector);
378 
396  explicit Vector (const Epetra_Map &parallel_partitioning);
397 
412  Vector (const Epetra_Map &parallel_partitioning,
413  const VectorBase &v);
414 
427  template <typename number>
428  void reinit (const Epetra_Map &parallel_partitioner,
429  const ::Vector<number> &v);
430 
443  void reinit (const Epetra_Map &parallel_partitioning,
444  const bool fast = false);
445 
458  template <typename Number>
459  Vector (const Epetra_Map &parallel_partitioning,
460  const ::Vector<Number> &v);
462 
480  explicit Vector (const IndexSet &parallel_partitioning,
481  const MPI_Comm &communicator = MPI_COMM_WORLD);
482 
494  Vector (const IndexSet &local,
495  const IndexSet &ghost,
496  const MPI_Comm &communicator = MPI_COMM_WORLD);
497 
512  Vector (const IndexSet &parallel_partitioning,
513  const VectorBase &v,
514  const MPI_Comm &communicator = MPI_COMM_WORLD);
515 
528  template <typename Number>
529  Vector (const IndexSet &parallel_partitioning,
530  const ::Vector<Number> &v,
531  const MPI_Comm &communicator = MPI_COMM_WORLD);
532 
548  void reinit (const IndexSet &parallel_partitioning,
549  const MPI_Comm &communicator = MPI_COMM_WORLD,
550  const bool fast = false);
551 
577  void reinit (const IndexSet &locally_owned_entries,
578  const IndexSet &ghost_entries,
579  const MPI_Comm &communicator = MPI_COMM_WORLD,
580  const bool vector_writable = false);
582  };
583 
584 
585 
586 
587 // ------------------- inline and template functions --------------
588 
589 
598  inline
599  void swap (Vector &u, Vector &v)
600  {
601  u.swap (v);
602  }
603 
604 
605 #ifndef DOXYGEN
606 
607  template <typename number>
608  Vector::Vector (const Epetra_Map &input_map,
609  const ::Vector<number> &v)
610  {
611  reinit (input_map, v);
612  }
613 
614 
615 
616  template <typename number>
617  Vector::Vector (const IndexSet &parallel_partitioner,
618  const ::Vector<number> &v,
619  const MPI_Comm &communicator)
620  {
621  *this = Vector(parallel_partitioner.make_trilinos_map (communicator, true),
622  v);
623  }
624 
625 
626 
627 
628  template <typename number>
629  void Vector::reinit (const Epetra_Map &parallel_partitioner,
630  const ::Vector<number> &v)
631  {
632  if (vector.get() == 0 || vector->Map().SameAs(parallel_partitioner) == false)
633  vector.reset (new Epetra_FEVector(parallel_partitioner));
634 
635  has_ghosts = vector->Map().UniqueGIDs()==false;
636 
637  const int size = parallel_partitioner.NumMyElements();
638 
639  // Need to copy out values, since the deal.II might not use doubles, so
640  // that a direct access is not possible.
641  for (int i=0; i<size; ++i)
642  (*vector)[0][i] = v(gid(parallel_partitioner,i));
643  }
644 
645 
646  inline
647  Vector &
648  Vector::operator = (const TrilinosScalar s)
649  {
651 
652  return *this;
653  }
654 
655 
656  template <typename Number>
657  Vector &
658  Vector::operator = (const ::Vector<Number> &v)
659  {
660  if (size() != v.size())
661  {
662  vector.reset (new Epetra_FEVector(Epetra_Map
663  (static_cast<TrilinosWrappers::types::int_type>(v.size()), 0,
664 #ifdef DEAL_II_WITH_MPI
665  Epetra_MpiComm(MPI_COMM_SELF)
666 #else
667  Epetra_SerialComm()
668 #endif
669  )));
670  }
671 
672  reinit (vector_partitioner(), v);
673  return *this;
674  }
675 
676 
677 #endif
678 
679  } /* end of namespace MPI */
680 
681 
682 
694  class Vector : public VectorBase
695  {
696  public:
700  typedef ::types::global_dof_index size_type;
701 
712  static const bool supports_distributed_data = false;
713 
719  Vector ();
720 
724  explicit Vector (const size_type n);
725 
735  explicit Vector (const Epetra_Map &partitioning);
736 
746  explicit Vector (const IndexSet &partitioning,
747  const MPI_Comm &communicator = MPI_COMM_WORLD);
748 
753  explicit Vector (const VectorBase &V);
754 
759  template <typename Number>
760  explicit Vector (const ::Vector<Number> &v);
761 
766  void reinit (const size_type n,
767  const bool fast = false);
768 
782  void reinit (const Epetra_Map &input_map,
783  const bool fast = false);
784 
798  void reinit (const IndexSet &input_map,
799  const MPI_Comm &communicator = MPI_COMM_WORLD,
800  const bool fast = false);
801 
806  void reinit (const VectorBase &V,
807  const bool fast = false,
808  const bool allow_different_maps = false);
809 
816  Vector &operator = (const TrilinosScalar s);
817 
822  Vector &
823  operator = (const MPI::Vector &V);
824 
828  template <typename Number>
829  Vector &
830  operator = (const ::Vector<Number> &V);
831 
836  Vector &
837  operator = (const Vector &V);
838 
850  void update_ghost_values () const;
851  };
852 
853 
854 
855 // ------------------- inline and template functions --------------
856 
857 
866  inline
867  void swap (Vector &u, Vector &v)
868  {
869  u.swap (v);
870  }
871 
872 
873 #ifndef DOXYGEN
874 
875  template <typename number>
876  Vector::Vector (const ::Vector<number> &v)
877  {
878  Epetra_LocalMap map ((TrilinosWrappers::types::int_type)v.size(), 0, Utilities::Trilinos::comm_self());
879  vector.reset (new Epetra_FEVector(map));
880  *this = v;
881  }
882 
883 
884 
885  inline
886  Vector &
887  Vector::operator = (const TrilinosScalar s)
888  {
890 
891  return *this;
892  }
893 
894 
895 
896  template <typename Number>
897  Vector &
898  Vector::operator = (const ::Vector<Number> &v)
899  {
900  if (size() != v.size())
901  {
902  vector.reset();
903 
904  Epetra_LocalMap map ((TrilinosWrappers::types::int_type)v.size(), 0,
906  vector.reset (new Epetra_FEVector(map));
907  }
908 
909  const Epetra_Map &map = vector_partitioner();
910  const TrilinosWrappers::types::int_type size = map.NumMyElements();
911 
912  Assert (map.MaxLID() == size-1,
913  ExcDimensionMismatch(map.MaxLID(), size-1));
914 
915  // Need to copy out values, since the
916  // deal.II might not use doubles, so
917  // that a direct access is not possible.
918  for (TrilinosWrappers::types::int_type i=0; i<size; ++i)
919  (*vector)[0][i] = v(i);
920 
921  return *this;
922  }
923 
924 
925 
926  inline
927  void
929  {}
930 
931 
932 #endif
933 
934 
935 }
936 
937 
940 DEAL_II_NAMESPACE_CLOSE
941 
942 #endif // DEAL_II_WITH_TRILINOS
943 
944 /*---------------------------- trilinos_vector.h ---------------------------*/
945 
946 #endif
947 /*---------------------------- 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
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