Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
trilinos_vector.h
1 // ---------------------------------------------------------------------
2 // @f$Id: trilinos_vector.h 30040 2013-07-18 17:06:48Z maier @f$
3 //
4 // Copyright (C) 2008 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__trilinos_vector_h
18 #define __deal2__trilinos_vector_h
19 
20 
21 #include <deal.II/base/config.h>
22 
23 #ifdef DEAL_II_WITH_TRILINOS
24 
25 # include <deal.II/base/std_cxx1x/shared_ptr.h>
26 # include <deal.II/base/subscriptor.h>
27 # include <deal.II/base/index_set.h>
28 # include <deal.II/base/utilities.h>
29 # include <deal.II/lac/exceptions.h>
30 # include <deal.II/lac/vector.h>
31 # include <deal.II/lac/trilinos_vector_base.h>
32 
33 # include "Epetra_Map.h"
34 # include "Epetra_LocalMap.h"
35 
36 DEAL_II_NAMESPACE_OPEN
37 
38 
39 // forward declaration
40 template <typename> class Vector;
41 
46 namespace TrilinosWrappers
47 {
48  class SparseMatrix;
49 
50  namespace
51  {
52 #ifndef DEAL_II_USE_LARGE_INDEX_TYPE
53  // define a helper function that queries the global ID of local ID of
54  // an Epetra_BlockMap object by calling either the 32- or 64-bit
55  // function necessary.
56  int gid(const Epetra_BlockMap &map, int i)
57  {
58  return map.GID(i);
59  }
60 #else
61  // define a helper function that queries the global ID of local ID of
62  // an Epetra_BlockMap object by calling either the 32- or 64-bit
63  // function necessary.
64  long long int gid(const Epetra_BlockMap &map, int i)
65  {
66  return map.GID64(i);
67  }
68 #endif
69  }
70 
79  namespace MPI
80  {
81  class BlockVector;
82 
193  class Vector : public VectorBase
194  {
195  public:
199  typedef ::types::global_dof_index size_type;
200 
212  static const bool supports_distributed_data = true;
213 
228  Vector ();
229 
234  Vector (const Vector &V);
235 
239  ~Vector ();
240 
276  void reinit (const VectorBase &v,
277  const bool fast = false,
278  const bool allow_different_maps = false);
279 
280  void reinit (const BlockVector &v,
281  const bool import_data = false);
282 
289  void reinit (const IndexSet &local,
290  const IndexSet &ghost,
291  const MPI_Comm &communicator = MPI_COMM_WORLD);
292 
303  Vector &operator = (const TrilinosScalar s);
304 
313  Vector &
314  operator = (const Vector &V);
315 
329  Vector &
330  operator = (const ::TrilinosWrappers::Vector &V);
331 
349  template <typename Number>
350  Vector &
351  operator = (const ::Vector<Number> &v);
352 
394  (const ::TrilinosWrappers::SparseMatrix &matrix,
395  const Vector &vector);
397 
412  explicit Vector (const Epetra_Map &parallel_partitioning);
413 
425  Vector (const Epetra_Map &parallel_partitioning,
426  const VectorBase &v);
427 
433  template <typename number>
434  void reinit (const Epetra_Map &parallel_partitioner,
435  const ::Vector<number> &v);
436 
444  void reinit (const Epetra_Map &parallel_partitioning,
445  const bool fast = false);
446 
453  template <typename Number>
454  Vector (const Epetra_Map &parallel_partitioning,
455  const ::Vector<Number> &v);
457 
471  explicit Vector (const IndexSet &parallel_partitioning,
472  const MPI_Comm &communicator = MPI_COMM_WORLD);
473 
477  Vector (const IndexSet &local,
478  const IndexSet &ghost,
479  const MPI_Comm &communicator = MPI_COMM_WORLD);
480 
492  Vector (const IndexSet &parallel_partitioning,
493  const VectorBase &v,
494  const MPI_Comm &communicator = MPI_COMM_WORLD);
495 
502  template <typename Number>
503  Vector (const IndexSet &parallel_partitioning,
504  const ::Vector<Number> &v,
505  const MPI_Comm &communicator = MPI_COMM_WORLD);
506 
516  void reinit (const IndexSet &parallel_partitioning,
517  const MPI_Comm &communicator = MPI_COMM_WORLD,
518  const bool fast = false);
520  };
521 
522 
523 
524 
525 // ------------------- inline and template functions --------------
526 
527 
536  inline
537  void swap (Vector &u, Vector &v)
538  {
539  u.swap (v);
540  }
541 
542 
543 #ifndef DOXYGEN
544 
545  template <typename number>
546  Vector::Vector (const Epetra_Map &input_map,
547  const ::Vector<number> &v)
548  {
549  reinit (input_map, v);
550  }
551 
552 
553 
554  template <typename number>
555  Vector::Vector (const IndexSet &parallel_partitioner,
556  const ::Vector<number> &v,
557  const MPI_Comm &communicator)
558  {
559  *this = Vector(parallel_partitioner.make_trilinos_map (communicator, true),
560  v);
561  }
562 
563 
564 
565 
566  template <typename number>
567  void Vector::reinit (const Epetra_Map &parallel_partitioner,
568  const ::Vector<number> &v)
569  {
570  if (vector.get() != 0 && vector->Map().SameAs(parallel_partitioner))
571  vector.reset (new Epetra_FEVector(parallel_partitioner));
572 
573  has_ghosts = vector->Map().UniqueGIDs()==false;
574 
575  const int size = parallel_partitioner.NumMyElements();
576 
577  // Need to copy out values, since the
578  // deal.II might not use doubles, so
579  // that a direct access is not possible.
580  for (int i=0; i<size; ++i)
581  (*vector)[0][i] = v(gid(parallel_partitioner,i));
582  }
583 
584 
585  inline
586  Vector &
587  Vector::operator = (const TrilinosScalar s)
588  {
590 
591  return *this;
592  }
593 
594 
595  template <typename Number>
596  Vector &
597  Vector::operator = (const ::Vector<Number> &v)
598  {
599  if (size() != v.size())
600  {
601  vector.reset (new Epetra_FEVector(Epetra_Map
602  (static_cast<TrilinosWrappers::types::int_type>(v.size()), 0,
603 #ifdef DEAL_II_WITH_MPI
604  Epetra_MpiComm(MPI_COMM_SELF)
605 #else
606  Epetra_SerialComm()
607 #endif
608  )));
609  }
610 
611  reinit (vector_partitioner(), v);
612  return *this;
613  }
614 
615 
616 #endif
617 
618  } /* end of namespace MPI */
619 
620 
621 
633  class Vector : public VectorBase
634  {
635  public:
639  typedef ::types::global_dof_index size_type;
640 
654  static const bool supports_distributed_data = false;
655 
664  Vector ();
665 
671  explicit Vector (const size_type n);
672 
685  explicit Vector (const Epetra_Map &partitioning);
686 
698  explicit Vector (const IndexSet &partitioning,
699  const MPI_Comm &communicator = MPI_COMM_WORLD);
700 
708  explicit Vector (const VectorBase &V);
709 
716  template <typename Number>
717  explicit Vector (const ::Vector<Number> &v);
718 
724  void reinit (const size_type n,
725  const bool fast = false);
726 
743  void reinit (const Epetra_Map &input_map,
744  const bool fast = false);
745 
762  void reinit (const IndexSet &input_map,
763  const MPI_Comm &communicator = MPI_COMM_WORLD,
764  const bool fast = false);
765 
772  void reinit (const VectorBase &V,
773  const bool fast = false,
774  const bool allow_different_maps = false);
775 
786  Vector &operator = (const TrilinosScalar s);
787 
794  Vector &
795  operator = (const MPI::Vector &V);
796 
801  template <typename Number>
802  Vector &
803  operator = (const ::Vector<Number> &V);
804 
810  Vector &
811  operator = (const Vector &V);
812 
828  void update_ghost_values () const;
829  };
830 
831 
832 
833 // ------------------- inline and template functions --------------
834 
835 
844  inline
845  void swap (Vector &u, Vector &v)
846  {
847  u.swap (v);
848  }
849 
850 
851 #ifndef DOXYGEN
852 
853  template <typename number>
854  Vector::Vector (const ::Vector<number> &v)
855  {
856  Epetra_LocalMap map ((TrilinosWrappers::types::int_type)v.size(), 0, Utilities::Trilinos::comm_self());
857  vector.reset (new Epetra_FEVector(map));
858  *this = v;
859  }
860 
861 
862 
863  inline
864  Vector &
865  Vector::operator = (const TrilinosScalar s)
866  {
868 
869  return *this;
870  }
871 
872 
873 
874  template <typename Number>
875  Vector &
876  Vector::operator = (const ::Vector<Number> &v)
877  {
878  if (size() != v.size())
879  {
880  vector.reset();
881 
882  Epetra_LocalMap map ((TrilinosWrappers::types::int_type)v.size(), 0,
884  vector.reset (new Epetra_FEVector(map));
885  }
886 
887  const Epetra_Map &map = vector_partitioner();
888  const TrilinosWrappers::types::int_type size = map.NumMyElements();
889 
890  Assert (map.MaxLID() == size-1,
891  ExcDimensionMismatch(map.MaxLID(), size-1));
892 
893  // Need to copy out values, since the
894  // deal.II might not use doubles, so
895  // that a direct access is not possible.
896  for (TrilinosWrappers::types::int_type i=0; i<size; ++i)
897  (*vector)[0][i] = v(i);
898 
899  return *this;
900  }
901 
902 
903 
904  inline
905  void
907  {}
908 
909 
910 #endif
911 
912 
913 }
914 
915 
918 DEAL_II_NAMESPACE_CLOSE
919 
920 #endif // DEAL_II_WITH_TRILINOS
921 
922 /*---------------------------- trilinos_vector.h ---------------------------*/
923 
924 #endif
925 /*---------------------------- 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)
#define Assert(cond, exc)
Definition: exceptions.h:299
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
std_cxx1x::shared_ptr< Epetra_FEVector > vector
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