Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
trilinos_parallel_block_vector.h
1 // ---------------------------------------------------------------------
2 // @f$Id: trilinos_parallel_block_vector.h 31932 2013-12-08 02:15:54Z heister @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_parallel_block_vector_h
18 #define __deal2__trilinos_parallel_block_vector_h
19 
20 
21 #include <deal.II/base/config.h>
22 
23 #ifdef DEAL_II_WITH_TRILINOS
24 
25 # include <deal.II/lac/trilinos_vector.h>
26 # include <deal.II/lac/block_indices.h>
27 # include <deal.II/lac/block_vector_base.h>
28 # include <deal.II/lac/exceptions.h>
29 
30 DEAL_II_NAMESPACE_OPEN
31 
32 // forward declaration
33 template <typename Number> class BlockVector;
34 
39 namespace TrilinosWrappers
40 {
41  // forward declaration
42  namespace MPI
43  {
44  class BlockVector;
45  }
46  class BlockVector;
47  class BlockSparseMatrix;
48 
49 
50  namespace MPI
51  {
72  class BlockVector : public BlockVectorBase<Vector>
73  {
74  public:
80 
86 
91  typedef BaseClass::value_type value_type;
92  typedef BaseClass::pointer pointer;
93  typedef BaseClass::const_pointer const_pointer;
94  typedef BaseClass::reference reference;
95  typedef BaseClass::const_reference const_reference;
96  typedef BaseClass::size_type size_type;
99 
104  BlockVector ();
105 
115  explicit BlockVector (const std::vector<Epetra_Map> &parallel_partitioning);
116 
127  explicit BlockVector (const std::vector<IndexSet> &parallel_partitioning,
128  const MPI_Comm &communicator = MPI_COMM_WORLD);
129 
135  BlockVector (const std::vector<IndexSet> &parallel_partitioning,
136  const std::vector<IndexSet> &ghost_values,
137  const MPI_Comm &communicator);
138 
139 
146  BlockVector (const BlockVector &V);
147 
158  explicit BlockVector (const size_type num_blocks);
159 
163  ~BlockVector ();
164 
171  BlockVector &
172  operator = (const value_type s);
173 
178  BlockVector &
179  operator = (const BlockVector &V);
180 
186  BlockVector &
187  operator = (const ::TrilinosWrappers::BlockVector &V);
188 
209  template <typename Number>
210  BlockVector &
211  operator = (const ::BlockVector<Number> &V);
212 
225  void reinit (const std::vector<Epetra_Map> &parallel_partitioning,
226  const bool fast = false);
227 
240  void reinit (const std::vector<IndexSet> &parallel_partitioning,
241  const MPI_Comm &communicator = MPI_COMM_WORLD,
242  const bool fast = false);
247  void reinit (const std::vector<IndexSet> &partitioning,
248  const std::vector<IndexSet> &ghost_values,
249  const MPI_Comm &communicator = MPI_COMM_WORLD);
250 
251 
278  void reinit (const BlockVector &V,
279  const bool fast = false);
280 
292  void reinit (const size_type num_blocks);
293 
330  const BlockVector &v);
331 
332 
342  void compress (const Epetra_CombineMode last_action) DEAL_II_DEPRECATED;
343 
348 
349 
364  bool is_compressed () const;
365 
369  bool has_ghost_elements() const;
370 
401  void swap (BlockVector &v);
402 
406  void print (std::ostream &out,
407  const unsigned int precision = 3,
408  const bool scientific = true,
409  const bool across = true) const;
410 
414  DeclException0 (ExcIteratorRangeDoesNotMatchVectorSize);
415 
419  DeclException0 (ExcNonMatchingBlockVectors);
420  };
421 
422 
423 
424  /*----------------------- Inline functions ----------------------------------*/
425 
426 
427  inline
429  {}
430 
431 
432 
433  inline
434  BlockVector::BlockVector (const std::vector<Epetra_Map> &parallel_partitioning)
435  {
436  reinit (parallel_partitioning, false);
437  }
438 
439 
440 
441  inline
442  BlockVector::BlockVector (const std::vector<IndexSet> &parallel_partitioning,
443  const MPI_Comm &communicator)
444  {
445  reinit (parallel_partitioning, communicator, false);
446  }
447 
448 
449  inline
450  BlockVector::BlockVector (const std::vector<IndexSet> &parallel_partitioning,
451  const std::vector<IndexSet> &ghost_values,
452  const MPI_Comm &communicator)
453  {
454  reinit(parallel_partitioning, ghost_values, communicator);
455  }
456 
457 
458  inline
459  BlockVector::BlockVector (const size_type num_blocks)
460  {
461  reinit (num_blocks);
462  }
463 
464 
465 
466  inline
468  :
470  {
471  this->components.resize (v.n_blocks());
472  this->block_indices = v.block_indices;
473 
474  for (size_type i=0; i<this->n_blocks(); ++i)
475  this->components[i] = v.components[i];
476  }
477 
478 
479 
480  inline
481  bool
483  {
484  bool compressed = true;
485  for (unsigned int row=0; row<n_blocks(); ++row)
486  if (block(row).is_compressed() == false)
487  {
488  compressed = false;
489  break;
490  }
491 
492  return compressed;
493  }
494 
495 
496 
497  inline
498  void
499  BlockVector::compress (const Epetra_CombineMode last_action)
500  {
501  ::VectorOperation::values last_action_ =
502  ::VectorOperation::unknown;
503  if (last_action == Add)
504  last_action_ = ::VectorOperation::add;
505  else if (last_action == Insert)
506  last_action_ = ::VectorOperation::insert;
507  else
508  AssertThrow(false, ExcNotImplemented());
509 
510  this->compress(last_action_);
511  }
512 
513 
514 
515  template <typename Number>
516  BlockVector &
517  BlockVector::operator = (const ::BlockVector<Number> &v)
518  {
519  if (n_blocks() != v.n_blocks())
520  {
521  std::vector<size_type> block_sizes (v.n_blocks(), 0);
522  block_indices.reinit (block_sizes);
523  if (components.size() != n_blocks())
524  components.resize(n_blocks());
525  }
526 
527  for (size_type i=0; i<this->n_blocks(); ++i)
528  this->components[i] = v.block(i);
529 
530  collect_sizes();
531 
532  return *this;
533  }
534 
535 
536  inline
537  bool
539  {
540  bool ghosted=block(0).has_ghost_elements();
541 #ifdef DEBUG
542  for (unsigned int i=0; i<this->n_blocks(); ++i)
543  Assert(block(i).has_ghost_elements()==ghosted, ExcInternalError());
544 #endif
545  return ghosted;
546  }
547 
548  inline
549  void
551  {
552  Assert (n_blocks() == v.n_blocks(),
554 
555  for (unsigned int row=0; row<n_blocks(); ++row)
556  block(row).swap (v.block(row));
557  }
558 
559 
560 
569  inline
570  void swap (BlockVector &u,
571  BlockVector &v)
572  {
573  u.swap (v);
574  }
575 
576  } /* end of namespace MPI */
577 
578 }
579 
582 DEAL_II_NAMESPACE_CLOSE
583 
584 #endif // DEAL_II_WITH_TRILINOS
585 
586 #endif
void import_nonlocal_data_for_fe(const TrilinosWrappers::BlockSparseMatrix &m, const BlockVector &v)
#define AssertThrow(cond, exc)
Definition: exceptions.h:362
void compress() DEAL_II_DEPRECATED
#define Assert(cond, exc)
Definition: exceptions.h:299
BlockType & block(const unsigned int i)
void reinit(const unsigned int n_blocks, const size_type n_elements_per_block)
unsigned int n_blocks() const
BlockCompressedSparsityPattern CompressedBlockSparsityPattern DEAL_II_DEPRECATED
BlockVector & operator=(const value_type s)
void print(std::ostream &out, const unsigned int precision=3, const bool scientific=true, const bool across=true) const
::ExceptionBase & ExcNotImplemented()
::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
void reinit(const std::vector< Epetra_Map > &parallel_partitioning, const bool fast=false)
::ExceptionBase & ExcInternalError()
std::vector< VectorType > components
void swap(BlockVector &u, BlockVector &v)
DeclException0(ExcIteratorRangeDoesNotMatchVectorSize)