include/deal.II/lac/block_sparse_matrix.h

00001 //---------------------------------------------------------------------------
00002 //    @f$Id: block_sparse_matrix.h 25345 2012-03-31 08:37:04Z bangerth @f$
00003 //
00004 //    Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2009, 2010, 2012 by the deal.II authors
00005 //
00006 //    This file is subject to QPL and may not be  distributed
00007 //    without copyright and license information. Please refer
00008 //    to the file deal.II/doc/license.html for the  text  and
00009 //    further information on this license.
00010 //
00011 //---------------------------------------------------------------------------
00012 #ifndef __deal2__block_sparse_matrix_h
00013 #define __deal2__block_sparse_matrix_h
00014 
00015 
00016 #include <deal.II/base/config.h>
00017 #include <deal.II/base/table.h>
00018 #include <deal.II/lac/block_matrix_base.h>
00019 #include <deal.II/lac/block_vector.h>
00020 #include <deal.II/lac/sparse_matrix.h>
00021 #include <deal.II/lac/block_sparsity_pattern.h>
00022 #include <deal.II/lac/exceptions.h>
00023 
00024 #include <cmath>
00025 
00026 DEAL_II_NAMESPACE_OPEN
00027 
00028 
00044 template <typename number>
00045 class BlockSparseMatrix : public BlockMatrixBase<SparseMatrix<number> >
00046 {
00047   public:
00052     typedef BlockMatrixBase<SparseMatrix<number> > BaseClass;
00053 
00058     typedef typename BaseClass::BlockType  BlockType;
00059 
00064     typedef typename BaseClass::value_type      value_type;
00065     typedef typename BaseClass::pointer         pointer;
00066     typedef typename BaseClass::const_pointer   const_pointer;
00067     typedef typename BaseClass::reference       reference;
00068     typedef typename BaseClass::const_reference const_reference;
00069     typedef typename BaseClass::size_type       size_type;
00070     typedef typename BaseClass::iterator        iterator;
00071     typedef typename BaseClass::const_iterator  const_iterator;
00072 
00098     BlockSparseMatrix ();
00099 
00121     BlockSparseMatrix (const BlockSparsityPattern &sparsity);
00122 
00126     virtual ~BlockSparseMatrix ();
00127 
00128 
00129 
00135     BlockSparseMatrix &
00136     operator = (const BlockSparseMatrix &);
00137 
00153     BlockSparseMatrix &
00154     operator = (const double d);
00155 
00168     void clear ();
00169 
00187     virtual void reinit (const BlockSparsityPattern &sparsity);
00189 
00201     bool empty () const;
00202 
00207     unsigned int get_row_length (const unsigned int row) const;
00208 
00218     unsigned int n_nonzero_elements () const;
00219 
00227     unsigned int n_actually_nonzero_elements (const double threshold = 0.0) const;
00228 
00241     const BlockSparsityPattern &
00242     get_sparsity_pattern () const;
00243 
00249     std::size_t memory_consumption () const;
00251 
00261     template <typename block_number>
00262     void vmult (BlockVector<block_number>       &dst,
00263                 const BlockVector<block_number> &src) const;
00264 
00272     template <typename block_number,
00273               typename nonblock_number>
00274     void vmult (BlockVector<block_number>          &dst,
00275                 const Vector<nonblock_number> &src) const;
00276 
00284     template <typename block_number,
00285               typename nonblock_number>
00286     void vmult (Vector<nonblock_number>    &dst,
00287                 const BlockVector<block_number> &src) const;
00288 
00296     template <typename nonblock_number>
00297     void vmult (Vector<nonblock_number>       &dst,
00298                 const Vector<nonblock_number> &src) const;
00299 
00308     template <typename block_number>
00309     void Tvmult (BlockVector<block_number>       &dst,
00310                  const BlockVector<block_number> &src) const;
00311 
00319     template <typename block_number,
00320               typename nonblock_number>
00321     void Tvmult (BlockVector<block_number>  &dst,
00322                  const Vector<nonblock_number> &src) const;
00323 
00331     template <typename block_number,
00332               typename nonblock_number>
00333     void Tvmult (Vector<nonblock_number>    &dst,
00334                  const BlockVector<block_number> &src) const;
00335 
00343     template <typename nonblock_number>
00344     void Tvmult (Vector<nonblock_number>       &dst,
00345                  const Vector<nonblock_number> &src) const;
00347 
00367     template <class BlockVectorType>
00368     void precondition_Jacobi (BlockVectorType       &dst,
00369                               const BlockVectorType &src,
00370                               const number           omega = 1.) const;
00371 
00379     template <typename number2>
00380     void precondition_Jacobi (Vector<number2>       &dst,
00381                               const Vector<number2> &src,
00382                               const number           omega = 1.) const;
00384 
00430     void print_formatted (std::ostream       &out,
00431                           const unsigned int  precision   = 3,
00432                           const bool          scientific  = true,
00433                           const unsigned int  width       = 0,
00434                           const char         *zero_string = " ",
00435                           const double        denominator = 1.) const;
00437 
00443     DeclException0 (ExcBlockDimensionMismatch);
00445 
00446   private:
00456     SmartPointer<const BlockSparsityPattern,BlockSparseMatrix<number> > sparsity_pattern;
00457 };
00458 
00459 
00460 
00462 /* ------------------------- Template functions ---------------------- */
00463 
00464 
00465 
00466 template <typename number>
00467 inline
00468 BlockSparseMatrix<number> &
00469 BlockSparseMatrix<number>::operator = (const double d)
00470 {
00471   Assert (d==0, ExcScalarAssignmentOnlyForZeroValue());
00472 
00473   for (unsigned int r=0; r<this->n_block_rows(); ++r)
00474     for (unsigned int c=0; c<this->n_block_cols(); ++c)
00475       this->block(r,c) = d;
00476 
00477   return *this;
00478 }
00479 
00480 
00481 
00482 template <typename number>
00483 template <typename block_number>
00484 inline
00485 void
00486 BlockSparseMatrix<number>::vmult (BlockVector<block_number>       &dst,
00487                                   const BlockVector<block_number> &src) const
00488 {
00489   BaseClass::vmult_block_block (dst, src);
00490 }
00491 
00492 
00493 
00494 template <typename number>
00495 template <typename block_number,
00496           typename nonblock_number>
00497 inline
00498 void
00499 BlockSparseMatrix<number>::vmult (BlockVector<block_number>     &dst,
00500                                   const Vector<nonblock_number> &src) const
00501 {
00502   BaseClass::vmult_block_nonblock (dst, src);
00503 }
00504 
00505 
00506 
00507 template <typename number>
00508 template <typename block_number,
00509           typename nonblock_number>
00510 inline
00511 void
00512 BlockSparseMatrix<number>::vmult (Vector<nonblock_number>         &dst,
00513                                   const BlockVector<block_number> &src) const
00514 {
00515   BaseClass::vmult_nonblock_block (dst, src);
00516 }
00517 
00518 
00519 
00520 template <typename number>
00521 template <typename nonblock_number>
00522 inline
00523 void
00524 BlockSparseMatrix<number>::vmult (Vector<nonblock_number>       &dst,
00525                                   const Vector<nonblock_number> &src) const
00526 {
00527   BaseClass::vmult_nonblock_nonblock (dst, src);
00528 }
00529 
00530 
00531 
00532 template <typename number>
00533 template <typename block_number>
00534 inline
00535 void
00536 BlockSparseMatrix<number>::Tvmult (BlockVector<block_number>       &dst,
00537                                    const BlockVector<block_number> &src) const
00538 {
00539   BaseClass::Tvmult_block_block (dst, src);
00540 }
00541 
00542 
00543 
00544 template <typename number>
00545 template <typename block_number,
00546           typename nonblock_number>
00547 inline
00548 void
00549 BlockSparseMatrix<number>::Tvmult (BlockVector<block_number>     &dst,
00550                                    const Vector<nonblock_number> &src) const
00551 {
00552   BaseClass::Tvmult_block_nonblock (dst, src);
00553 }
00554 
00555 
00556 
00557 template <typename number>
00558 template <typename block_number,
00559           typename nonblock_number>
00560 inline
00561 void
00562 BlockSparseMatrix<number>::Tvmult (Vector<nonblock_number>         &dst,
00563                                    const BlockVector<block_number> &src) const
00564 {
00565   BaseClass::Tvmult_nonblock_block (dst, src);
00566 }
00567 
00568 
00569 
00570 template <typename number>
00571 template <typename nonblock_number>
00572 inline
00573 void
00574 BlockSparseMatrix<number>::Tvmult (Vector<nonblock_number>       &dst,
00575                                    const Vector<nonblock_number> &src) const
00576 {
00577   BaseClass::Tvmult_nonblock_nonblock (dst, src);
00578 }
00579 
00580 
00581 
00582 template <typename number>
00583 template <class BlockVectorType>
00584 inline
00585 void
00586 BlockSparseMatrix<number>::
00587 precondition_Jacobi (BlockVectorType       &dst,
00588                      const BlockVectorType &src,
00589                      const number           omega) const
00590 {
00591   Assert (this->n_block_rows() == this->n_block_cols(), ExcNotQuadratic());
00592   Assert (dst.n_blocks() == this->n_block_rows(),
00593           ExcDimensionMismatch(dst.n_blocks(), this->n_block_rows()));
00594   Assert (src.n_blocks() == this->n_block_cols(),
00595           ExcDimensionMismatch(src.n_blocks(), this->n_block_cols()));
00596 
00597                                    // do a diagonal preconditioning. uses only
00598                                    // the diagonal blocks of the matrix
00599   for (unsigned int i=0; i<this->n_block_rows(); ++i)
00600     this->block(i,i).precondition_Jacobi (dst.block(i),
00601                                           src.block(i),
00602                                           omega);
00603 }
00604 
00605 
00606 
00607 template <typename number>
00608 template <typename number2>
00609 inline
00610 void
00611 BlockSparseMatrix<number>::
00612 precondition_Jacobi (Vector<number2>       &dst,
00613                      const Vector<number2> &src,
00614                      const number           omega) const
00615 {
00616                                    // check number of blocks. the sizes of the
00617                                    // single block is checked in the function
00618                                    // we call
00619   Assert (this->n_block_cols() == 1,
00620           ExcMessage ("This function only works if the matrix has "
00621                       "a single block"));
00622   Assert (this->n_block_rows() == 1,
00623           ExcMessage ("This function only works if the matrix has "
00624                       "a single block"));
00625 
00626                                    // do a diagonal preconditioning. uses only
00627                                    // the diagonal blocks of the matrix
00628   this->block(0,0).precondition_Jacobi (dst, src, omega);
00629 }
00630 
00631 
00632 DEAL_II_NAMESPACE_CLOSE
00633 
00634 #endif    // __deal2__block_sparse_matrix_h
00635 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines

deal.II documentation generated on Thu May 17 2012 20:04:40 by doxygen 1.7.3