Reference documentation for deal.II version Git f2670cf 2014-12-18 09:41:14 +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 
395  explicit Vector (const Epetra_Map &parallel_partitioning);
396 
410  Vector (const Epetra_Map &parallel_partitioning,
411  const VectorBase &v);
412 
424  template <typename number>
425  void reinit (const Epetra_Map &parallel_partitioner,
426  const ::Vector<number> &v);
427 
439  void reinit (const Epetra_Map &parallel_partitioning,
440  const bool fast = false);
441 
453  template <typename Number>
454  Vector (const Epetra_Map &parallel_partitioning,
455  const ::Vector<Number> &v);
457 
474  explicit Vector (const IndexSet &parallel_partitioning,
475  const MPI_Comm &communicator = MPI_COMM_WORLD);
476 
487  Vector (const IndexSet &local,
488  const IndexSet &ghost,
489  const MPI_Comm &communicator = MPI_COMM_WORLD);
490 
504  Vector (const IndexSet &parallel_partitioning,
505  const VectorBase &v,
506  const MPI_Comm &communicator = MPI_COMM_WORLD);
507 
519  template <typename Number>
520  Vector (const IndexSet &parallel_partitioning,
521  const ::Vector<Number> &v,
522  const MPI_Comm &communicator = MPI_COMM_WORLD);
523 
538  void reinit (const IndexSet &parallel_partitioning,
539  const MPI_Comm &communicator = MPI_COMM_WORLD,
540  const bool fast = false);
541 
566  void reinit (const IndexSet &locally_owned_entries,
567  const IndexSet &ghost_entries,
568  const MPI_Comm &communicator = MPI_COMM_WORLD,
569  const bool vector_writable = false);
571  };
572 
573 
574 
575 
576 // ------------------- inline and template functions --------------
577 
578 
587  inline
588  void swap (Vector &u, Vector &v)
589  {
590  u.swap (v);
591  }
592 
593 
594 #ifndef DOXYGEN
595 
596  template <typename number>
597  Vector::Vector (const Epetra_Map &input_map,
598  const ::Vector<number> &v)
599  {
600  reinit (input_map, v);
601  }
602 
603 
604 
605  template <typename number>
606  Vector::Vector (const IndexSet &parallel_partitioner,
607  const ::Vector<number> &v,
608  const MPI_Comm &communicator)
609  {
610  *this = Vector(parallel_partitioner.make_trilinos_map (communicator, true),
611  v);
612  }
613 
614 
615 
616 
617  template <typename number>
618  void Vector::reinit (const Epetra_Map &parallel_partitioner,
619  const ::Vector<number> &v)
620  {
621  if (vector.get() == 0 || vector->Map().SameAs(parallel_partitioner) == false)
622  vector.reset (new Epetra_FEVector(parallel_partitioner));
623 
624  has_ghosts = vector->Map().UniqueGIDs()==false;
625 
626  const int size = parallel_partitioner.NumMyElements();
627 
628  // Need to copy out values, since the deal.II might not use doubles, so
629  // that a direct access is not possible.
630  for (int i=0; i<size; ++i)
631  (*vector)[0][i] = v(gid(parallel_partitioner,i));
632  }
633 
634 
635  inline
636  Vector &
637  Vector::operator = (const TrilinosScalar s)
638  {
640 
641  return *this;
642  }
643 
644 
645  template <typename Number>
646  Vector &
647  Vector::operator = (const ::Vector<Number> &v)
648  {
649  if (size() != v.size())
650  {
651  vector.reset (new Epetra_FEVector(Epetra_Map
652  (static_cast<TrilinosWrappers::types::int_type>(v.size()), 0,
653 #ifdef DEAL_II_WITH_MPI
654  Epetra_MpiComm(MPI_COMM_SELF)
655 #else
656  Epetra_SerialComm()
657 #endif
658  )));
659  }
660 
661  reinit (vector_partitioner(), v);
662  return *this;
663  }
664 
665 
666 #endif
667 
668  } /* end of namespace MPI */
669 
670 
671 
683  class Vector : public VectorBase
684  {
685  public:
689  typedef ::types::global_dof_index size_type;
690 
701  static const bool supports_distributed_data = false;
702 
708  Vector ();
709 
713  explicit Vector (const size_type n);
714 
724  explicit Vector (const Epetra_Map &partitioning);
725 
735  explicit Vector (const IndexSet &partitioning,
736  const MPI_Comm &communicator = MPI_COMM_WORLD);
737 
742  explicit Vector (const VectorBase &V);
743 
748  template <typename Number>
749  explicit Vector (const ::Vector<Number> &v);
750 
755  void reinit (const size_type n,
756  const bool fast = false);
757 
771  void reinit (const Epetra_Map &input_map,
772  const bool fast = false);
773 
787  void reinit (const IndexSet &input_map,
788  const MPI_Comm &communicator = MPI_COMM_WORLD,
789  const bool fast = false);
790 
795  void reinit (const VectorBase &V,
796  const bool fast = false,
797  const bool allow_different_maps = false);
798 
805  Vector &operator = (const TrilinosScalar s);
806 
811  Vector &
812  operator = (const MPI::Vector &V);
813 
817  template <typename Number>
818  Vector &
819  operator = (const ::Vector<Number> &V);
820 
825  Vector &
826  operator = (const Vector &V);
827 
839  void update_ghost_values () const;
840  };
841 
842 
843 
844 // ------------------- inline and template functions --------------
845 
846 
855  inline
856  void swap (Vector &u, Vector &v)
857  {
858  u.swap (v);
859  }
860 
861 
862 #ifndef DOXYGEN
863 
864  template <typename number>
865  Vector::Vector (const ::Vector<number> &v)
866  {
867  Epetra_LocalMap map ((TrilinosWrappers::types::int_type)v.size(), 0, Utilities::Trilinos::comm_self());
868  vector.reset (new Epetra_FEVector(map));
869  *this = v;
870  }
871 
872 
873 
874  inline
875  Vector &
876  Vector::operator = (const TrilinosScalar s)
877  {
879 
880  return *this;
881  }
882 
883 
884 
885  template <typename Number>
886  Vector &
887  Vector::operator = (const ::Vector<Number> &v)
888  {
889  if (size() != v.size())
890  {
891  vector.reset();
892 
893  Epetra_LocalMap map ((TrilinosWrappers::types::int_type)v.size(), 0,
895  vector.reset (new Epetra_FEVector(map));
896  }
897 
898  const Epetra_Map &map = vector_partitioner();
899  const TrilinosWrappers::types::int_type size = map.NumMyElements();
900 
901  Assert (map.MaxLID() == size-1,
902  ExcDimensionMismatch(map.MaxLID(), size-1));
903 
904  // Need to copy out values, since the
905  // deal.II might not use doubles, so
906  // that a direct access is not possible.
907  for (TrilinosWrappers::types::int_type i=0; i<size; ++i)
908  (*vector)[0][i] = v(i);
909 
910  return *this;
911  }
912 
913 
914 
915  inline
916  void
918  {}
919 
920 
921 #endif
922 
923 
924 }
925 
926 
929 DEAL_II_NAMESPACE_CLOSE
930 
931 #endif // DEAL_II_WITH_TRILINOS
932 
933 /*---------------------------- trilinos_vector.h ---------------------------*/
934 
935 #endif
936 /*---------------------------- 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:291
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