// $Id$
// Version: $Name$
//
-// Copyright (C) 2001, 2002, 2003, 2004, 2005 by the deal.II authors
+// Copyright (C) 2001, 2002, 2003, 2004, 2005, 2006 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
const unsigned int col,
const double prefix = 1.,
const bool transpose = false);
+
+ /**
+ * Add an entry like with enter,
+ * but use PointerMatrixAux for
+ * matrices not having functions
+ * vmult_add() and TVmult_add().
+ */
+ template <class MATRIX>
+ void enter_aux (VectorMemory<Vector<number> >& mem,
+ const MATRIX &matrix,
+ const unsigned int row,
+ const unsigned int col,
+ const double prefix = 1.,
+ const bool transpose = false);
+
/**
* Delete all entries, i.e. reset
{
public:
/**
- * Constructor initializing all
- * data fields.
+ * Constructor initializing
+ * all data fields. A
+ * PointerMatrix object is
+ * generated for
+ * <tt>matrix</tt>.
*/
template<class MATRIX>
Entry (const MATRIX& matrix,
unsigned row, unsigned int col,
double prefix, bool transpose);
+
+ /**
+ * Constructor initializing
+ * all data fields. A
+ * PointerMatrixAux object is
+ * generated for
+ * <tt>matrix</tt>.
+ */
+ template<class MATRIX>
+ Entry (const MATRIX& matrix,
+ unsigned row, unsigned int col,
+ double prefix, bool transpose,
+ VectorMemory<Vector<number> >& mem);
+
/**
* Copy constructor
* invalidating the old
const double prefix = 1.,
const bool transpose = false);
+ /**
+ * Enter a block. This calls
+ * BlockMatrixArray::enter_aux(). Remember
+ * that the diagonal blocks
+ * should actually be inverse
+ * matrices or preconditioners.
+ */
+ template <class MATRIX>
+ void enter_aux (VectorMemory<Vector<double> >& mem,
+ const MATRIX &matrix,
+ const unsigned int row,
+ const unsigned int col,
+ const double prefix = 1.,
+ const bool transpose = false);
+
/**
* Preconditioning.
*/
template <typename number>
template <class MATRIX>
inline
-BlockMatrixArray<number>::Entry::Entry (const MATRIX& matrix,
- unsigned row, unsigned int col,
- double prefix, bool transpose)
+BlockMatrixArray<number>::Entry::Entry (
+ const MATRIX& m,
+ unsigned int row,
+ unsigned int col,
+ double prefix,
+ bool transpose)
:
row (row),
col (col),
prefix (prefix),
transpose (transpose),
- matrix (new PointerMatrix<MATRIX, Vector<number> >(&matrix, typeid(*this).name()))
+ matrix (new PointerMatrix<MATRIX, Vector<number> >(&m, typeid(*this).name()))
+{}
+
+
+
+template <typename number>
+template <class MATRIX>
+inline
+BlockMatrixArray<number>::Entry::Entry (
+ const MATRIX& m,
+ unsigned row,
+ unsigned int col,
+ double prefix,
+ bool transpose,
+ VectorMemory<Vector<number> >& mem)
+ :
+ row (row),
+ col (col),
+ prefix (prefix),
+ transpose (transpose),
+ matrix (new PointerMatrixAux<MATRIX, Vector<number> >(&mem, &m, typeid(*this).name()))
{}
template <class MATRIX>
inline
void
-BlockMatrixArray<number>::enter (const MATRIX& matrix,
- unsigned row, unsigned int col,
- double prefix, bool transpose)
+BlockMatrixArray<number>::enter (
+ const MATRIX& matrix,
+ unsigned int row,
+ unsigned int col,
+ double prefix,
+ bool transpose)
{
Assert (mem != 0, ExcNotInitialized());
Assert(row<n_block_rows(), ExcIndexRange(row, 0, n_block_rows()));
}
+template <typename number>
+template <class MATRIX>
+inline
+void
+BlockMatrixArray<number>::enter_aux (
+ VectorMemory<Vector<number> >& mem,
+ const MATRIX& matrix,
+ unsigned int row,
+ unsigned int col,
+ double prefix,
+ bool transpose)
+{
+ Assert(row<n_block_rows(), ExcIndexRange(row, 0, n_block_rows()));
+ Assert(col<n_block_cols(), ExcIndexRange(col, 0, n_block_cols()));
+ entries.push_back(Entry(matrix, row, col, prefix, transpose, mem));
+}
+
+
template<typename number>
struct BlockMatrixArrayPointerMatrixLess
{
+template <typename number>
+template <class MATRIX>
+inline
+void
+BlockTrianglePrecondition<number>::enter_aux (
+ VectorMemory<Vector<double> >& mem,
+ const MATRIX& matrix,
+ unsigned int row,
+ unsigned int col,
+ double prefix,
+ bool transpose)
+{
+ BlockMatrixArray<number>::enter_aux(mem, matrix, row, col, prefix, transpose);
+}
+
+
+
#endif // DOXYGEN
#endif
// $Id$
// Version: $Name$
//
-// Copyright (C) 2002, 2003, 2004, 2005 by the deal.II authors
+// Copyright (C) 2002, 2003, 2004, 2005, 2006 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
#include <base/subscriptor.h>
#include <base/smartpointer.h>
+template<class VECTOR> class VectorMemory;
+
/*! @addtogroup Matrix2
*@{
*/
};
+/**
+ * A pointer to be used as a matrix. This class stores a pointer to a
+ * matrix and can be used as a matrix itself in iterative methods.
+ *
+ * The main purpose for the existence of this class is its base class,
+ * which only has a vector as template argument. Therefore, this
+ * interface provides an abstract base class for matrices.
+ *
+ * This class differs form PointerMatrix by its additional
+ * VectorMemory object and by the fact that it implements the
+ * functions vmult_add() and Tvmult_add() only using vmult() and
+ * Tvmult() of the MATRIX.
+ *
+ * @author Guido Kanschat 2006
+ */
+template<class MATRIX, class VECTOR>
+class PointerMatrixAux : public PointerMatrixBase<VECTOR>
+{
+ public:
+ /**
+ * Constructor. The pointer in the
+ * argument is stored in this
+ * class. As usual, the lifetime of
+ * <tt>*M</tt> must be longer than the
+ * one of the PointerMatrixAux.
+ *
+ * If <tt>M</tt> is zero, no
+ * matrix is stored.
+ */
+ PointerMatrixAux (VectorMemory<VECTOR>* mem,
+ const MATRIX* M=0);
+
+ /**
+ * Constructor. The name argument
+ * is used to identify the
+ * SmartPointer for this object.
+ */
+ PointerMatrixAux(VectorMemory<VECTOR>* mem,
+ const char* name);
+
+ /**
+ * Constructor. <tt>M</tt> points
+ * to a matrix which must live
+ * longer than the
+ * PointerMatrixAux. The name
+ * argument is used to identify
+ * the SmartPointer for this
+ * object.
+ */
+ PointerMatrixAux(VectorMemory<VECTOR>* mem,
+ const MATRIX* M,
+ const char* name);
+
+ // Use doc from base class
+ virtual void clear();
+
+ /**
+ * Return whether the object is
+ * empty.
+ */
+ bool empty () const;
+
+ /**
+ * Assign a new matrix
+ * pointer. Deletes the old pointer
+ * and releases its matrix.
+ * @see SmartPointer
+ */
+ const PointerMatrixAux& operator= (const MATRIX* M);
+
+ /**
+ * Matrix-vector product.
+ */
+ virtual void vmult (VECTOR& dst,
+ const VECTOR& src) const;
+
+ /**
+ * Tranposed matrix-vector product.
+ */
+ virtual void Tvmult (VECTOR& dst,
+ const VECTOR& src) const;
+
+ /**
+ * Matrix-vector product, adding to
+ * <tt>dst</tt>.
+ */
+ virtual void vmult_add (VECTOR& dst,
+ const VECTOR& src) const;
+
+ /**
+ * Tranposed matrix-vector product,
+ * adding to <tt>dst</tt>.
+ */
+ virtual void Tvmult_add (VECTOR& dst,
+ const VECTOR& src) const;
+
+ private:
+ /**
+ * Return the address of the
+ * matrix for comparison.
+ */
+ virtual const void* get() const;
+
+ /**
+ * Object for getting the
+ * auxiliary vector.
+ */
+ SmartPointer<VectorMemory<VECTOR> > mem;
+
+ /**
+ * The pointer to the actual matrix.
+ */
+ SmartPointer<const MATRIX> m;
+};
+
+
/**
* This function helps you creating a PointerMatrixBase object if you
* do not want to provide the full template arguments of
}
+//----------------------------------------------------------------------//
+
+
+template<class MATRIX, class VECTOR>
+PointerMatrixAux<MATRIX, VECTOR>::PointerMatrixAux (
+ VectorMemory<VECTOR>* mem,
+ const MATRIX* M)
+ : mem(mem, typeid(*this).name()),
+ m(M, typeid(*this).name())
+{}
+
+
+template<class MATRIX, class VECTOR>
+PointerMatrixAux<MATRIX, VECTOR>::PointerMatrixAux (
+ VectorMemory<VECTOR>* mem,
+ const char* name)
+ : mem(mem, name),
+ m(0, name)
+{}
+
+
+template<class MATRIX, class VECTOR>
+PointerMatrixAux<MATRIX, VECTOR>::PointerMatrixAux (
+ VectorMemory<VECTOR>* mem,
+ const MATRIX* M,
+ const char* name)
+ : mem(mem, name),
+ m(M, name)
+{}
+
+
+template<class MATRIX, class VECTOR>
+inline void
+PointerMatrixAux<MATRIX, VECTOR>::clear ()
+{
+ m = 0;
+}
+
+
+template<class MATRIX, class VECTOR>
+inline const PointerMatrixAux<MATRIX, VECTOR>&
+PointerMatrixAux<MATRIX, VECTOR>::operator= (const MATRIX* M)
+{
+ m = M;
+ return *this;
+}
+
+
+template<class MATRIX, class VECTOR>
+inline bool
+PointerMatrixAux<MATRIX, VECTOR>::empty () const
+{
+ if (m == 0)
+ return true;
+ return m->empty();
+}
+
+template<class MATRIX, class VECTOR>
+inline void
+PointerMatrixAux<MATRIX, VECTOR>::vmult (VECTOR& dst,
+ const VECTOR& src) const
+{
+ Assert (m != 0, ExcNotInitialized());
+ m->vmult (dst, src);
+}
+
+
+template<class MATRIX, class VECTOR>
+inline void
+PointerMatrixAux<MATRIX, VECTOR>::Tvmult (VECTOR& dst,
+ const VECTOR& src) const
+{
+ Assert (m != 0, ExcNotInitialized());
+ m->Tvmult (dst, src);
+}
+
+
+template<class MATRIX, class VECTOR>
+inline void
+PointerMatrixAux<MATRIX, VECTOR>::vmult_add (VECTOR& dst,
+ const VECTOR& src) const
+{
+ Assert (m != 0, ExcNotInitialized());
+ VECTOR* v = mem->alloc();
+ v->reinit(dst);
+ m->vmult (*v, src);
+ dst += *v;
+ mem->free(v);
+}
+
+
+template<class MATRIX, class VECTOR>
+inline void
+PointerMatrixAux<MATRIX, VECTOR>::Tvmult_add (VECTOR& dst,
+ const VECTOR& src) const
+{
+ Assert (m != 0, ExcNotInitialized());
+ VECTOR* v = mem->alloc();
+ v->reinit(dst);
+ m->Tvmult (*v, src);
+ dst += *v;
+ mem->free(v);
+}
+
+
+template<class MATRIX, class VECTOR>
+inline const void*
+PointerMatrixAux<MATRIX, VECTOR>::get () const
+{
+ return m;
+}
+
+
#endif