All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
trilinos_vector.h
1 //---------------------------------------------------------------------------
2 // @f$Id: trilinos_vector.h 29753 2013-06-05 13:50:55Z turcksin @f$
3 //
4 // Copyright (C) 2008, 2009, 2010, 2012, 2013 by the deal.II authors
5 //
6 // This file is subject to QPL and may not be distributed
7 // without copyright and license information. Please refer
8 // to the file deal.II/doc/license.html for the text and
9 // further information on this license.
10 //
11 //---------------------------------------------------------------------------
12 #ifndef __deal2__trilinos_vector_h
13 #define __deal2__trilinos_vector_h
14 
15 
16 #include <deal.II/base/config.h>
17 
18 #ifdef DEAL_II_WITH_TRILINOS
19 
20 # include <deal.II/base/std_cxx1x/shared_ptr.h>
21 # include <deal.II/base/subscriptor.h>
22 # include <deal.II/base/index_set.h>
23 # include <deal.II/base/utilities.h>
24 # include <deal.II/lac/exceptions.h>
25 # include <deal.II/lac/vector.h>
26 # include <deal.II/lac/trilinos_vector_base.h>
27 
28 # include "Epetra_Map.h"
29 # include "Epetra_LocalMap.h"
30 
31 DEAL_II_NAMESPACE_OPEN
32 
33 
34 // forward declaration
35 template <typename> class Vector;
36 
41 namespace TrilinosWrappers
42 {
43  class SparseMatrix;
44 
45  namespace
46  {
47 #ifndef DEAL_II_USE_LARGE_INDEX_TYPE
48  // define a helper function that queries the global ID of local ID of
49  // an Epetra_BlockMap object by calling either the 32- or 64-bit
50  // function necessary.
51  int gid(const Epetra_BlockMap &map, int i)
52  {
53  return map.GID(i);
54  }
55 #else
56  // define a helper function that queries the global ID of local ID of
57  // an Epetra_BlockMap object by calling either the 32- or 64-bit
58  // function necessary.
59  long long int gid(const Epetra_BlockMap &map, int i)
60  {
61  return map.GID64(i);
62  }
63 #endif
64  }
65 
74  namespace MPI
75  {
76  class BlockVector;
77 
188  class Vector : public VectorBase
189  {
190  public:
194  typedef ::types::global_dof_index size_type;
195 
207  static const bool supports_distributed_data = true;
208 
223  Vector ();
224 
229  Vector (const Vector &V);
230 
234  ~Vector ();
235 
271  void reinit (const VectorBase &v,
272  const bool fast = false,
273  const bool allow_different_maps = false);
274 
275  void reinit (const BlockVector &v,
276  const bool import_data = false);
277 
288  Vector &operator = (const TrilinosScalar s);
289 
298  Vector &
299  operator = (const Vector &V);
300 
314  Vector &
315  operator = (const ::TrilinosWrappers::Vector &V);
316 
334  template <typename Number>
335  Vector &
336  operator = (const ::Vector<Number> &v);
337 
379  (const ::TrilinosWrappers::SparseMatrix &matrix,
380  const Vector &vector);
382 
397  Vector (const Epetra_Map &parallel_partitioning);
398 
410  explicit Vector (const Epetra_Map &parallel_partitioning,
411  const VectorBase &v);
412 
418  template <typename number>
419  void reinit (const Epetra_Map &parallel_partitioner,
420  const ::Vector<number> &v);
421 
429  void reinit (const Epetra_Map &parallel_partitioning,
430  const bool fast = false);
431 
438  template <typename Number>
439  explicit Vector (const Epetra_Map &parallel_partitioning,
440  const ::Vector<Number> &v);
442 
456  Vector (const IndexSet &parallel_partitioning,
457  const MPI_Comm &communicator = MPI_COMM_WORLD);
458 
470  explicit Vector (const IndexSet &parallel_partitioning,
471  const VectorBase &v,
472  const MPI_Comm &communicator = MPI_COMM_WORLD);
473 
480  template <typename Number>
481  explicit Vector (const IndexSet &parallel_partitioning,
482  const ::Vector<Number> &v,
483  const MPI_Comm &communicator = MPI_COMM_WORLD);
484 
494  void reinit (const IndexSet &parallel_partitioning,
495  const MPI_Comm &communicator = MPI_COMM_WORLD,
496  const bool fast = false);
498  };
499 
500 
501 
502 
503 // ------------------- inline and template functions --------------
504 
505 
514  inline
515  void swap (Vector &u, Vector &v)
516  {
517  u.swap (v);
518  }
519 
520 
521 #ifndef DOXYGEN
522 
523  template <typename number>
524  Vector::Vector (const Epetra_Map &input_map,
525  const ::Vector<number> &v)
526  {
527  reinit (input_map, v);
528  }
529 
530 
531 
532  template <typename number>
533  Vector::Vector (const IndexSet &parallel_partitioner,
534  const ::Vector<number> &v,
535  const MPI_Comm &communicator)
536  {
537  *this = Vector(parallel_partitioner.make_trilinos_map (communicator, true),
538  v);
539  }
540 
541 
542 
543 
544  template <typename number>
545  void Vector::reinit (const Epetra_Map &parallel_partitioner,
546  const ::Vector<number> &v)
547  {
548  if (vector.get() != 0 && vector->Map().SameAs(parallel_partitioner))
549  vector.reset (new Epetra_FEVector(parallel_partitioner));
550 
551  has_ghosts = vector->Map().UniqueGIDs()==false;
552 
553  const int size = parallel_partitioner.NumMyElements();
554 
555  // Need to copy out values, since the
556  // deal.II might not use doubles, so
557  // that a direct access is not possible.
558  for (int i=0; i<size; ++i)
559  (*vector)[0][i] = v(gid(parallel_partitioner,i));
560  }
561 
562 
563  inline
564  Vector &
565  Vector::operator = (const TrilinosScalar s)
566  {
568 
569  return *this;
570  }
571 
572 
573  template <typename Number>
574  Vector &
575  Vector::operator = (const ::Vector<Number> &v)
576  {
577  if (size() != v.size())
578  {
579  vector.reset (new Epetra_FEVector(Epetra_Map
580  (static_cast<TrilinosWrappers::types::int_type>(v.size()), 0,
581 #ifdef DEAL_II_WITH_MPI
582  Epetra_MpiComm(MPI_COMM_SELF)
583 #else
584  Epetra_SerialComm()
585 #endif
586  )));
587  }
588 
589  reinit (vector_partitioner(), v);
590  return *this;
591  }
592 
593 
594 #endif
595 
596  } /* end of namespace MPI */
597 
598 
599 
611  class Vector : public VectorBase
612  {
613  public:
617  typedef ::types::global_dof_index size_type;
618 
632  static const bool supports_distributed_data = false;
633 
642  Vector ();
643 
649  Vector (const size_type n);
650 
663  Vector (const Epetra_Map &partitioning);
664 
676  Vector (const IndexSet &partitioning,
677  const MPI_Comm &communicator = MPI_COMM_WORLD);
678 
686  explicit Vector (const VectorBase &V);
687 
694  template <typename Number>
695  explicit Vector (const ::Vector<Number> &v);
696 
702  void reinit (const size_type n,
703  const bool fast = false);
704 
721  void reinit (const Epetra_Map &input_map,
722  const bool fast = false);
723 
740  void reinit (const IndexSet &input_map,
741  const MPI_Comm &communicator = MPI_COMM_WORLD,
742  const bool fast = false);
743 
750  void reinit (const VectorBase &V,
751  const bool fast = false,
752  const bool allow_different_maps = false);
753 
764  Vector &operator = (const TrilinosScalar s);
765 
772  Vector &
773  operator = (const MPI::Vector &V);
774 
779  template <typename Number>
780  Vector &
781  operator = (const ::Vector<Number> &V);
782 
788  Vector &
789  operator = (const Vector &V);
790 
806  void update_ghost_values () const;
807  };
808 
809 
810 
811 // ------------------- inline and template functions --------------
812 
813 
822  inline
823  void swap (Vector &u, Vector &v)
824  {
825  u.swap (v);
826  }
827 
828 
829 #ifndef DOXYGEN
830 
831  template <typename number>
832  Vector::Vector (const ::Vector<number> &v)
833  {
834  Epetra_LocalMap map ((TrilinosWrappers::types::int_type)v.size(), 0, Utilities::Trilinos::comm_self());
835  vector.reset (new Epetra_FEVector(map));
836  *this = v;
837  }
838 
839 
840 
841  inline
842  Vector &
843  Vector::operator = (const TrilinosScalar s)
844  {
846 
847  return *this;
848  }
849 
850 
851 
852  template <typename Number>
853  Vector &
854  Vector::operator = (const ::Vector<Number> &v)
855  {
856  if (size() != v.size())
857  {
858  vector.reset();
859 
860  Epetra_LocalMap map ((TrilinosWrappers::types::int_type)v.size(), 0,
862  vector.reset (new Epetra_FEVector(map));
863  }
864 
865  const Epetra_Map &map = vector_partitioner();
866  const TrilinosWrappers::types::int_type size = map.NumMyElements();
867 
868  Assert (map.MaxLID() == size-1,
869  ExcDimensionMismatch(map.MaxLID(), size-1));
870 
871  // Need to copy out values, since the
872  // deal.II might not use doubles, so
873  // that a direct access is not possible.
874  for (TrilinosWrappers::types::int_type i=0; i<size; ++i)
875  (*vector)[0][i] = v(i);
876 
877  return *this;
878  }
879 
880 
881 
882  inline
883  void
885  {}
886 
887 
888 #endif
889 
890 
891 }
892 
893 
896 DEAL_II_NAMESPACE_CLOSE
897 
898 #endif // DEAL_II_WITH_TRILINOS
899 
900 /*---------------------------- trilinos_vector.h ---------------------------*/
901 
902 #endif
903 /*---------------------------- trilinos_vector.h ---------------------------*/
904 

deal.II
generated on Thu Jun 20 2013 13:43:44 by doxygen 1.8.3.1.
Hosting by IWR Universität Heidelberg