moved to new module dune-solvers

[[Imported from SVN: r2346]]

dune-solvers/iterationsteps/blockgsstep.cc

+template<class OperatorType, class DiscFuncType, class BitVectorType>
+inline
+DiscFuncType BlockGSStep<OperatorType, DiscFuncType, BitVectorType>::getSol()
+{
+    return *(this->x_);
+}
+
+template<class OperatorType, class DiscFuncType, class BitVectorType>
+inline 
+void BlockGSStep<OperatorType, DiscFuncType, BitVectorType>::
+residual(int index, VectorBlock& r) const
+{
+    const OperatorType& mat = *this->mat_;
+
+    typedef typename OperatorType::row_type RowType;
+    const RowType& row = mat[index];
+
+    typedef typename RowType::ConstIterator ColumnIterator;
+
+    r = (*this->rhs_)[index];
+
+    /* The following loop subtracts
+     * \f[ sum_i = \sum_j A_{ij}w_j \f]
+     */
+    ColumnIterator cIt    = row.template begin();
+    ColumnIterator cEndIt = row.template end();
+    
+    for (; cIt!=cEndIt; ++cIt) {
+        // r_i -= A_ij x_j
+        cIt->mmv((*this->x_)[cIt.index()], r);
+    }
+
+}
+
+template<class OperatorType, class DiscFuncType, class BitVectorType>
+inline
+void BlockGSStep<OperatorType, DiscFuncType, BitVectorType>::iterate()
+{
+    const OperatorType& mat = *this->mat_;
+
+    for (int i=0; i<this->x_->size(); i++) {
+
+        /** \todo Handle more general boundary conditions */
+        if ((*this->ignoreNodes_)[i][0]) {
+            continue;
+        }
+
+        VectorBlock r;
+        residual(i, r);
+
+        // Compute x_i += A_{i,i}^{-1} r[i] 
+        VectorBlock v;
+        VectorBlock& x = (*this->x_)[i];
+
+        // No obstacle --> solve linear problem
+        mat[i][i].solve(v, r);
+        
+        // Add correction
+        x += v;
+
+    }
+
+}

dune-solvers/iterationsteps/blockgsstep.hh

+#ifndef DUNE_BLOCK_GAUSS_SEIDEL_STEP_HH
+#define DUNE_BLOCK_GAUSS_SEIDEL_STEP_HH
+
+#include <dune/common/bitsetvector.hh>
+
+#include "lineariterationstep.hh"
+
+template<class OperatorType, 
+         class DiscFuncType,
+         class BitVectorType = Dune::BitSetVector<DiscFuncType::block_type::dimension> >
+         class BlockGSStep : public LinearIterationStep<OperatorType, DiscFuncType, BitVectorType>
+    {
+
+        typedef typename DiscFuncType::block_type VectorBlock;
+
+        enum {BlockSize = VectorBlock::dimension};
+
+    public:
+        
+        //! Default constructor.  Doesn't init anything
+        BlockGSStep() {}
+
+        //! Constructor with a linear problem
+        BlockGSStep(const OperatorType& mat, DiscFuncType& x, DiscFuncType& rhs)
+            : LinearIterationStep<OperatorType,DiscFuncType>(mat, x, rhs)
+        {}
+
+        //! Perform one iteration
+        virtual void iterate();
+
+        virtual DiscFuncType getSol();
+
+        virtual void residual(int index, VectorBlock& r) const;
+
+    };
+
+
+#include "blockgsstep.cc"
+
+#endif