Skip to content
Snippets Groups Projects
Commit ad5e4167 authored by Oliver Sander's avatar Oliver Sander
Browse files

Removed the obsolete class MGTransfer. From now on, please use the

equivalent MultiGridTransfer instead.

[[Imported from SVN: r1992]]
parent aaab6860
Branches
Tags
No related merge requests found
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#include <dune/fem/dofmanager.hh>
#include <dune/fem/lagrangebase.hh>
#include <dune/common/functionspace.hh>
template<class DiscFuncType>
template<class FunctionSpaceType>
void Dune::MGTransfer<DiscFuncType>::setup(const FunctionSpaceType& cS, const FunctionSpaceType& fS)
{
coarseLevel = cS.level();
fineLevel = fS.level();
typedef typename FunctionSpaceType::GridType GridType;
int rows = fS.size();
int cols = cS.size();
const GridType& grid = fS.getGrid();
//
matrix_.resize(rows, cols, GridType::dimension*10);
matrix_.clear();
typedef typename GridType::template codim<0>::LevelIterator ElementIterator;
typedef typename FunctionSpaceType::BaseFunctionSetType BaseFunctionSetType;
ElementIterator cIt = grid.template lbegin<0>(coarseLevel);
ElementIterator cEndIt = grid.template lend<0>(coarseLevel);
for (; cIt != cEndIt; ++cIt) {
const BaseFunctionSetType& coarseBaseSet = fS.getBaseFunctionSet( *cIt );
const int numCoarseBaseFct = coarseBaseSet.getNumberOfBaseFunctions();
typedef typename GridType::template codim<0>::Entity EntityType;
typedef typename EntityType::HierarchicIterator HierarchicIterator;
HierarchicIterator fIt = cIt->hbegin(fineLevel);
HierarchicIterator fEndIt = cIt->hend(fineLevel);
for (; fIt != fEndIt; ++fIt) {
if (fIt->level()==cIt->level())
continue;
const BaseFunctionSetType& fineBaseSet = fS.getBaseFunctionSet( *fIt );
const int numFineBaseFct = fineBaseSet.getNumberOfBaseFunctions();
for (int i=0; i<numCoarseBaseFct; i++) {
int globalCoarse = fS.mapToGlobal(*cIt, i);
for (int j=0; j<numFineBaseFct; j++) {
int globalFine = fS.mapToGlobal(*fIt, j);
// Evaluate coarse grid base function at the location of the fine grid dof
// first determine local fine grid dof position
FieldVector<double, GridType::dimension> local = cIt->geometry().local(fIt->geometry()[j]);
// Evaluate coarse grid base function
typename FunctionSpaceType::Range value;
FieldVector<int, 0> diffVariable;
coarseBaseSet.evaluate(i, diffVariable, local, value);
matrix_.set(globalFine, globalCoarse, value[0]);
}
}
}
}
}
template<class DiscFuncType>
void Dune::MGTransfer<DiscFuncType>::prolong(const DiscFuncType& f, DiscFuncType& t) const
{
assert(t.getFunctionSpace().size() == matrix_.rows());
assert(f.getFunctionSpace().size() == matrix_.cols());
typedef typename DiscFuncType::DofIteratorType Iterator;
typedef typename DiscFuncType::ConstDofIteratorType ConstIterator;
typedef typename SparseRowMatrix<double>::ColumnIterator ColumnIterator;
Iterator tIt = t.dbegin();
ConstIterator fIt = f.dbegin();
for(int row=0; row<matrix_.rows(); row++) {
*tIt = 0.0;
/** \todo Remove casts */
ColumnIterator cIt = const_cast<SparseRowMatrix<double>&>(matrix_).template rbegin(row);
ColumnIterator cEndIt = const_cast<SparseRowMatrix<double>&>(matrix_).template rend(row);
for(; cIt!=cEndIt; ++cIt) {
*tIt += *cIt * fIt[cIt.col()];
}
++tIt;
}
}
template<class DiscFuncType>
void Dune::MGTransfer<DiscFuncType>::restrict (const DiscFuncType & f, DiscFuncType& t) const
{
assert(f.getFunctionSpace().size() == matrix_.rows());
assert(t.getFunctionSpace().size() == matrix_.cols());
typedef typename DiscFuncType::DofIteratorType IteratorType;
typedef typename DiscFuncType::ConstDofIteratorType ConstIteratorType;
typedef typename SparseRowMatrix<double>::ColumnIterator ColumnIterator;
t.clear();
IteratorType tIt = t.dbegin();
ConstIteratorType fIt = f.dbegin();
for (int row=0; row<matrix_.rows(); row++) {
/** \todo Remove casts */
ColumnIterator cIt = const_cast<SparseRowMatrix<double>&>(matrix_).template rbegin(row);
ColumnIterator cEndIt = const_cast<SparseRowMatrix<double>&>(matrix_).template rend(row);
for(; cIt!=cEndIt; ++cIt)
tIt[cIt.col()] += fIt[row] * (*cIt);
}
}
template<class DiscFuncType>
Dune::SparseRowMatrix<double> Dune::MGTransfer<DiscFuncType>::
galerkinRestrict(const Dune::SparseRowMatrix<double>& fineMat) const
{
typedef typename SparseRowMatrix<double>::ColumnIterator ColumnIterator;
SparseRowMatrix<double> result(matrix_.cols() ,matrix_.cols() , fineMat.numNonZeros());
for (int row=0; row<fineMat.rows(); row++) {
ColumnIterator cIt = const_cast<SparseRowMatrix<double>&>(fineMat).template rbegin(row);
ColumnIterator cEndIt = const_cast<SparseRowMatrix<double>&>(fineMat).template rend(row);
for(; cIt!=cEndIt; ++cIt) {
double mvalue = *cIt;
ColumnIterator tciIt = const_cast<SparseRowMatrix<double>&>(matrix_).template rbegin(row);
ColumnIterator tciEndIt = const_cast<SparseRowMatrix<double>&>(matrix_).template rend(row);
for (; tciIt!=tciEndIt; ++tciIt) {
double fac = mvalue* (*tciIt);
ColumnIterator tcjIt = const_cast<SparseRowMatrix<double>&>(matrix_).template rbegin(cIt.col());
ColumnIterator tcjEndIt = const_cast<SparseRowMatrix<double>&>(matrix_).template rend(cIt.col());
for (; tcjIt!=tcjEndIt; ++tcjIt) {
//result.add( tciIt.col(), tcjIt.col(), fac* (*tcjIt));
result.add( tcjIt.col(), tciIt.col(), fac* (*tcjIt));
}
}
}
}
return result;
}
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#ifndef DUNE_MG_TRANSFER_HH
#define DUNE_MG_TRANSFER_HH
#include <dune/fem/feop/spmatrix.hh>
namespace Dune {
/** \brief Standard multigrid restriction and prolongation operator
*
* Currently only works for first-order Lagrangian elements!
*/
template<class DiscFuncType>
class MGTransfer {
//typedef typename DiscFuncType::FunctionSpaceType FunctionSpaceType;
public:
/** \brief Sets up the operator between levels cL and fL
*
* \param fS The function space hierarchy between levels of which we're mapping
*/
template <class FunctionSpaceType>
void setup(const FunctionSpaceType& cS, const FunctionSpaceType& fS);
/** \brief Restrict level fL of f and store the result in level cL of t
*/
void restrict (const DiscFuncType & f, DiscFuncType &t) const;
/** \brief Prolong level cL of f and store the result in level fL of t
*/
void prolong(const DiscFuncType& f, DiscFuncType &t) const;
/** \brief Galerkin assemble a coarse stiffness matrix
*/
SparseRowMatrix<double> galerkinRestrict(const SparseRowMatrix<double>& fineMat) const;
const SparseRowMatrix<double>& getMatrix() const {return matrix_;}
protected:
int coarseLevel;
int fineLevel;
/** \todo Should we extract the value type from DiscFuncType? */
SparseRowMatrix<double> matrix_;
};
}
#include "mgtransfer.cc"
#endif
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment