Reference documentation for deal.II version SVN Revision 32792
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
trilinos_vector.h
1 // ---------------------------------------------------------------------
2 // @f$Id: trilinos_vector.h 32300 2014-01-24 16:52:51Z kronbichler @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_WITH_64BIT_INDICES
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 
205  class Vector : public VectorBase
206  {
207  public:
211  typedef ::types::global_dof_index size_type;
212 
223  static const bool supports_distributed_data = true;
224 
234  Vector ();
235 
239  Vector (const Vector &V);
240 
244  ~Vector ();
245 
265  void reinit (const VectorBase &v,
266  const bool fast = false,
267  const bool allow_different_maps = false);
268 
272  void reinit (const BlockVector &v,
273  const bool import_data = false);
274 
281  Vector &operator = (const TrilinosScalar s);
282 
288  Vector &
289  operator = (const Vector &V);
290 
298  Vector &
299  operator = (const ::TrilinosWrappers::Vector &V);
300 
309  template <typename Number>
310  Vector &
311  operator = (const ::Vector<Number> &v);
312 
331  (const ::TrilinosWrappers::SparseMatrix &matrix,
332  const Vector &vector);
334 
344  explicit Vector (const Epetra_Map &parallel_partitioning);
345 
352  Vector (const Epetra_Map &parallel_partitioning,
353  const VectorBase &v);
354 
359  template <typename number>
360  void reinit (const Epetra_Map &parallel_partitioner,
361  const ::Vector<number> &v);
362 
367  void reinit (const Epetra_Map &parallel_partitioning,
368  const bool fast = false);
369 
374  template <typename Number>
375  Vector (const Epetra_Map &parallel_partitioning,
376  const ::Vector<Number> &v);
378 
388  explicit Vector (const IndexSet &parallel_partitioning,
389  const MPI_Comm &communicator = MPI_COMM_WORLD);
390 
394  Vector (const IndexSet &local,
395  const IndexSet &ghost,
396  const MPI_Comm &communicator = MPI_COMM_WORLD);
397 
404  Vector (const IndexSet &parallel_partitioning,
405  const VectorBase &v,
406  const MPI_Comm &communicator = MPI_COMM_WORLD);
407 
412  template <typename Number>
413  Vector (const IndexSet &parallel_partitioning,
414  const ::Vector<Number> &v,
415  const MPI_Comm &communicator = MPI_COMM_WORLD);
416 
423  void reinit (const IndexSet &parallel_partitioning,
424  const MPI_Comm &communicator = MPI_COMM_WORLD,
425  const bool fast = false);
426 
445  void reinit (const IndexSet &locally_owned_entries,
446  const IndexSet &ghost_entries,
447  const MPI_Comm &communicator = MPI_COMM_WORLD,
448  const bool vector_writable = false);
450  };
451 
452 
453 
454 
455 // ------------------- inline and template functions --------------
456 
457 
466  inline
467  void swap (Vector &u, Vector &v)
468  {
469  u.swap (v);
470  }
471 
472 
473 #ifndef DOXYGEN
474 
475  template <typename number>
476  Vector::Vector (const Epetra_Map &input_map,
477  const ::Vector<number> &v)
478  {
479  reinit (input_map, v);
480  }
481 
482 
483 
484  template <typename number>
485  Vector::Vector (const IndexSet &parallel_partitioner,
486  const ::Vector<number> &v,
487  const MPI_Comm &communicator)
488  {
489  *this = Vector(parallel_partitioner.make_trilinos_map (communicator, true),
490  v);
491  }
492 
493 
494 
495 
496  template <typename number>
497  void Vector::reinit (const Epetra_Map &parallel_partitioner,
498  const ::Vector<number> &v)
499  {
500  if (vector.get() == 0 || vector->Map().SameAs(parallel_partitioner) == false)
501  vector.reset (new Epetra_FEVector(parallel_partitioner));
502 
503  has_ghosts = vector->Map().UniqueGIDs()==false;
504 
505  const int size = parallel_partitioner.NumMyElements();
506 
507  // Need to copy out values, since the deal.II might not use doubles, so
508  // that a direct access is not possible.
509  for (int i=0; i<size; ++i)
510  (*vector)[0][i] = v(gid(parallel_partitioner,i));
511  }
512 
513 
514  inline
515  Vector &
516  Vector::operator = (const TrilinosScalar s)
517  {
519 
520  return *this;
521  }
522 
523 
524  template <typename Number>
525  Vector &
526  Vector::operator = (const ::Vector<Number> &v)
527  {
528  if (size() != v.size())
529  {
530  vector.reset (new Epetra_FEVector(Epetra_Map
531  (static_cast<TrilinosWrappers::types::int_type>(v.size()), 0,
532 #ifdef DEAL_II_WITH_MPI
533  Epetra_MpiComm(MPI_COMM_SELF)
534 #else
535  Epetra_SerialComm()
536 #endif
537  )));
538  }
539 
540  reinit (vector_partitioner(), v);
541  return *this;
542  }
543 
544 
545 #endif
546 
547  } /* end of namespace MPI */
548 
549 
550 
562  class Vector : public VectorBase
563  {
564  public:
568  typedef ::types::global_dof_index size_type;
569 
580  static const bool supports_distributed_data = false;
581 
587  Vector ();
588 
592  explicit Vector (const size_type n);
593 
600  explicit Vector (const Epetra_Map &partitioning);
601 
609  explicit Vector (const IndexSet &partitioning,
610  const MPI_Comm &communicator = MPI_COMM_WORLD);
611 
616  explicit Vector (const VectorBase &V);
617 
622  template <typename Number>
623  explicit Vector (const ::Vector<Number> &v);
624 
629  void reinit (const size_type n,
630  const bool fast = false);
631 
641  void reinit (const Epetra_Map &input_map,
642  const bool fast = false);
643 
653  void reinit (const IndexSet &input_map,
654  const MPI_Comm &communicator = MPI_COMM_WORLD,
655  const bool fast = false);
656 
661  void reinit (const VectorBase &V,
662  const bool fast = false,
663  const bool allow_different_maps = false);
664 
671  Vector &operator = (const TrilinosScalar s);
672 
677  Vector &
678  operator = (const MPI::Vector &V);
679 
683  template <typename Number>
684  Vector &
685  operator = (const ::Vector<Number> &V);
686 
691  Vector &
692  operator = (const Vector &V);
693 
705  void update_ghost_values () const;
706  };
707 
708 
709 
710 // ------------------- inline and template functions --------------
711 
712 
721  inline
722  void swap (Vector &u, Vector &v)
723  {
724  u.swap (v);
725  }
726 
727 
728 #ifndef DOXYGEN
729 
730  template <typename number>
731  Vector::Vector (const ::Vector<number> &v)
732  {
733  Epetra_LocalMap map ((TrilinosWrappers::types::int_type)v.size(), 0, Utilities::Trilinos::comm_self());
734  vector.reset (new Epetra_FEVector(map));
735  *this = v;
736  }
737 
738 
739 
740  inline
741  Vector &
742  Vector::operator = (const TrilinosScalar s)
743  {
745 
746  return *this;
747  }
748 
749 
750 
751  template <typename Number>
752  Vector &
753  Vector::operator = (const ::Vector<Number> &v)
754  {
755  if (size() != v.size())
756  {
757  vector.reset();
758 
759  Epetra_LocalMap map ((TrilinosWrappers::types::int_type)v.size(), 0,
761  vector.reset (new Epetra_FEVector(map));
762  }
763 
764  const Epetra_Map &map = vector_partitioner();
765  const TrilinosWrappers::types::int_type size = map.NumMyElements();
766 
767  Assert (map.MaxLID() == size-1,
768  ExcDimensionMismatch(map.MaxLID(), size-1));
769 
770  // Need to copy out values, since the
771  // deal.II might not use doubles, so
772  // that a direct access is not possible.
773  for (TrilinosWrappers::types::int_type i=0; i<size; ++i)
774  (*vector)[0][i] = v(i);
775 
776  return *this;
777  }
778 
779 
780 
781  inline
782  void
784  {}
785 
786 
787 #endif
788 
789 
790 }
791 
792 
795 DEAL_II_NAMESPACE_CLOSE
796 
797 #endif // DEAL_II_WITH_TRILINOS
798 
799 /*---------------------------- trilinos_vector.h ---------------------------*/
800 
801 #endif
802 /*---------------------------- 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