00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012 #ifndef __deal2__trilinos_vector_h
00013 #define __deal2__trilinos_vector_h
00014
00015
00016 #include <deal.II/base/config.h>
00017
00018 #ifdef DEAL_II_USE_TRILINOS
00019
00020 # include <deal.II/base/std_cxx1x/shared_ptr.h>
00021 # include <deal.II/base/subscriptor.h>
00022 # include <deal.II/base/index_set.h>
00023 # include <deal.II/base/utilities.h>
00024 # include <deal.II/lac/exceptions.h>
00025 # include <deal.II/lac/vector.h>
00026 # include <deal.II/lac/trilinos_vector_base.h>
00027
00028 # include "Epetra_Map.h"
00029 # include "Epetra_LocalMap.h"
00030
00031 DEAL_II_NAMESPACE_OPEN
00032
00033
00034
00035 template <typename> class Vector;
00036
00041 namespace TrilinosWrappers
00042 {
00043 class SparseMatrix;
00052 namespace MPI
00053 {
00054 class BlockVector;
00055
00166 class Vector : public VectorBase
00167 {
00168 public:
00183 Vector ();
00184
00189 Vector (const Vector &V);
00190
00194 ~Vector ();
00195
00231 void reinit (const VectorBase &v,
00232 const bool fast = false,
00233 const bool allow_different_maps = false);
00234
00235 void reinit (const BlockVector &v,
00236 const bool import_data = false);
00237
00248 Vector & operator = (const TrilinosScalar s);
00249
00258 Vector &
00259 operator = (const Vector &V);
00260
00274 Vector &
00275 operator = (const ::TrilinosWrappers::Vector &V);
00276
00294 template <typename Number>
00295 Vector &
00296 operator = (const ::Vector<Number> &v);
00297
00338 void import_nonlocal_data_for_fe
00339 (const ::TrilinosWrappers::SparseMatrix &matrix,
00340 const Vector &vector);
00342
00357 Vector (const Epetra_Map ¶llel_partitioning);
00358
00370 explicit Vector (const Epetra_Map ¶llel_partitioning,
00371 const VectorBase &v);
00372
00378 template <typename number>
00379 void reinit (const Epetra_Map ¶llel_partitioner,
00380 const ::Vector<number> &v);
00381
00389 void reinit (const Epetra_Map ¶llel_partitioning,
00390 const bool fast = false);
00391
00398 template <typename Number>
00399 explicit Vector (const Epetra_Map ¶llel_partitioning,
00400 const ::Vector<Number> &v);
00402
00416 Vector (const IndexSet ¶llel_partitioning,
00417 const MPI_Comm &communicator = MPI_COMM_WORLD);
00418
00430 explicit Vector (const IndexSet ¶llel_partitioning,
00431 const VectorBase &v,
00432 const MPI_Comm &communicator = MPI_COMM_WORLD);
00433
00440 template <typename Number>
00441 explicit Vector (const IndexSet ¶llel_partitioning,
00442 const ::Vector<Number> &v,
00443 const MPI_Comm &communicator = MPI_COMM_WORLD);
00444
00451 void reinit (const IndexSet ¶llel_partitioning,
00452 const MPI_Comm &communicator = MPI_COMM_WORLD,
00453 const bool fast = false);
00455 };
00456
00457
00458
00459
00460
00461
00462
00471 inline
00472 void swap (Vector &u, Vector &v)
00473 {
00474 u.swap (v);
00475 }
00476
00477
00478 #ifndef DOXYGEN
00479
00480 template <typename number>
00481 Vector::Vector (const Epetra_Map &input_map,
00482 const ::Vector<number> &v)
00483 {
00484 reinit (input_map, v);
00485 }
00486
00487
00488
00489 template <typename number>
00490 Vector::Vector (const IndexSet ¶llel_partitioner,
00491 const ::Vector<number> &v,
00492 const MPI_Comm &communicator)
00493 {
00494 *this = Vector(parallel_partitioner.make_trilinos_map (communicator, true),
00495 v);
00496 }
00497
00498
00499
00500
00501 template <typename number>
00502 void Vector::reinit (const Epetra_Map ¶llel_partitioner,
00503 const ::Vector<number> &v)
00504 {
00505 if (vector.get() != 0 && vector->Map().SameAs(parallel_partitioner))
00506 vector.reset (new Epetra_FEVector(parallel_partitioner));
00507
00508 const int size = parallel_partitioner.NumMyElements();
00509
00510
00511
00512
00513 for (int i=0; i<size; ++i)
00514 (*vector)[0][i] = v(parallel_partitioner.GID(i));
00515 }
00516
00517
00518 inline
00519 Vector &
00520 Vector::operator = (const TrilinosScalar s)
00521 {
00522 VectorBase::operator = (s);
00523
00524 return *this;
00525 }
00526
00527
00528 template <typename Number>
00529 Vector &
00530 Vector::operator = (const ::Vector<Number> &v)
00531 {
00532 if (size() != v.size())
00533 {
00534 vector.reset (new Epetra_FEVector(Epetra_Map
00535 (v.size(), 0,
00536 #ifdef DEAL_II_COMPILER_SUPPORTS_MPI
00537 Epetra_MpiComm(MPI_COMM_SELF)
00538 #else
00539 Epetra_SerialComm()
00540 #endif
00541 )));
00542 }
00543
00544 reinit (vector_partitioner(), v);
00545 return *this;
00546 }
00547
00548
00549 #endif
00550
00551 }
00552
00553
00554
00566 class Vector : public VectorBase
00567 {
00568 public:
00577 Vector ();
00578
00584 Vector (const unsigned int n);
00585
00598 Vector (const Epetra_Map &partitioning);
00599
00611 Vector (const IndexSet &partitioning,
00612 const MPI_Comm &communicator = MPI_COMM_WORLD);
00613
00621 explicit Vector (const VectorBase &V);
00622
00629 template <typename Number>
00630 explicit Vector (const ::Vector<Number> &v);
00631
00637 void reinit (const unsigned int n,
00638 const bool fast = false);
00639
00656 void reinit (const Epetra_Map &input_map,
00657 const bool fast = false);
00658
00675 void reinit (const IndexSet &input_map,
00676 const MPI_Comm &communicator = MPI_COMM_WORLD,
00677 const bool fast = false);
00678
00685 void reinit (const VectorBase &V,
00686 const bool fast = false,
00687 const bool allow_different_maps = false);
00688
00699 Vector & operator = (const TrilinosScalar s);
00700
00707 Vector &
00708 operator = (const MPI::Vector &V);
00709
00714 template <typename Number>
00715 Vector &
00716 operator = (const ::Vector<Number> &V);
00717
00723 Vector &
00724 operator = (const Vector &V);
00725
00741 void update_ghost_values () const;
00742 };
00743
00744
00745
00746
00747
00748
00757 inline
00758 void swap (Vector &u, Vector &v)
00759 {
00760 u.swap (v);
00761 }
00762
00763
00764 #ifndef DOXYGEN
00765
00766 template <typename number>
00767 Vector::Vector (const ::Vector<number> &v)
00768 {
00769 Epetra_LocalMap map ((int)v.size(), 0, Utilities::Trilinos::comm_self());
00770 vector.reset (new Epetra_FEVector(map));
00771 *this = v;
00772 }
00773
00774
00775
00776 inline
00777 Vector &
00778 Vector::operator = (const TrilinosScalar s)
00779 {
00780 VectorBase::operator = (s);
00781
00782 return *this;
00783 }
00784
00785
00786
00787 template <typename Number>
00788 Vector &
00789 Vector::operator = (const ::Vector<Number> &v)
00790 {
00791 if (size() != v.size())
00792 {
00793 vector.reset();
00794
00795 Epetra_LocalMap map ((int)v.size(), 0,
00796 Utilities::Trilinos::comm_self());
00797 vector.reset (new Epetra_FEVector(map));
00798 }
00799
00800 const Epetra_Map & map = vector_partitioner();
00801 const int size = map.NumMyElements();
00802
00803 Assert (map.MaxLID() == size-1,
00804 ExcDimensionMismatch(map.MaxLID(), size-1));
00805
00806
00807
00808
00809 for (int i=0; i<size; ++i)
00810 (*vector)[0][i] = v(i);
00811
00812 return *this;
00813 }
00814
00815
00816
00817 inline
00818 void
00819 Vector::update_ghost_values () const
00820 {}
00821
00822
00823 #endif
00824
00825
00826 }
00827
00828
00831 DEAL_II_NAMESPACE_CLOSE
00832
00833 #endif // DEAL_II_USE_TRILINOS
00834
00835
00836
00837 #endif
00838
00839