diff --git a/doc/manual/examples.cc b/doc/manual/examples.cc
new file mode 100755
index 0000000000000000000000000000000000000000..27b9cb810f19bb6c5af8de72ef945182a382539e
--- /dev/null
+++ b/doc/manual/examples.cc
@@ -0,0 +1,239 @@
+#include <config.h>
+
+#define DIMENSION 2
+#define ALUGRID 1
+//#define UGGRID 1
+#define EMBEDDED_PYTHON 1
+#undef EMBEDDED_PYTHON
+
+#define FE_VERBOSE
+
+// disable embedded python if python was not found
+#ifndef HAVE_PYTHON
+    #undef EMBEDDED_PYTHON
+#endif
+
+#if EMBEDDED_PYTHON
+    #include <Python.h>
+#endif
+
+#include <fstream>
+#include <iostream>
+#include <iostream>
+#include <vector>
+
+#include <tr1/memory>
+
+// dune includes ******************************************************
+#ifdef UGGRID
+    #include <dune/grid/uggrid.hh>
+    #include <dune/grid/io/file/dgfparser/dgfug.hh>
+#endif
+
+#ifdef ALUGRID
+    #include <dune/grid/alugrid.hh>
+    #include <dune/grid/io/file/dgfparser/dgfalu.hh>
+#endif
+
+#include <dune/common/configparser.hh>
+#include <dune/common/shared_ptr.hh>
+
+#include <dune/istl/bcrsmatrix.hh>
+
+#include <dune/grid/common/genericreferenceelements.hh>
+
+// dune-fufem includes ******************************************************
+#include <dune/fufem/boundarypatch.hh>
+#include <dune/fufem/prolongboundarypatch.hh>
+
+#if EMBEDDED_PYTHON
+    #include <dune/fufem/functions/pythonfunction.hh>
+#endif
+#include <dune/fufem/functions/functions.hh>
+#include <dune/fufem/functions/basisgridfunction.hh>
+
+#include <dune/fufem/functionspacebases/p1nodalbasis.hh>
+
+#include <dune/fufem/functiontools/namedfunctionmap.hh>
+#include <dune/fufem/functiontools/boundarydofs.hh>
+
+#include <dune/fufem/assemblers/operatorassembler.hh>
+#include <dune/fufem/assemblers/functionalassembler.hh>
+#include <dune/fufem/assemblers/transferoperatorassembler.hh>
+
+#include <dune/fufem/assemblers/localassemblers/laplaceassembler.hh>
+#include <dune/fufem/assemblers/localassemblers/massassembler.hh>
+#include <dune/fufem/assemblers/localassemblers/lumpedmassassembler.hh>
+#include <dune/fufem/assemblers/localassemblers/l2functionalassembler.hh>
+
+#if HAVE_AMIRAMESH
+    #include <dune/fufem/functiontools/amirameshbasiswriter.hh>
+#endif
+#include <dune/fufem/functiontools/vtkbasiswriter.hh>
+
+// dune-solvers includes ***************************************************
+#include <dune/solvers/norms/norm.hh>
+#include <dune/solvers/norms/sumnorm.hh>
+#include <dune/solvers/norms/energynorm.hh>
+#include <dune/solvers/norms/diagnorm.hh>
+#include <dune/solvers/norms/fullnorm.hh>
+#include <dune/solvers/norms/blocknorm.hh>
+
+// dune-tnnmg includes *******************************************************
+#include <dune/tnnmg/problem-classes/blocknonlineargsproblem.hh>
+#include <dune/tnnmg/problem-classes/nonlinearity.hh>
+#include <dune/tnnmg/nonlinearities/zerononlinearity.hh>
+#include <dune/tnnmg/nonlinearities/shiftednonlinearity.hh>
+#include <dune/tnnmg/iterationsteps/genericnonlinearjacobi.hh>
+#include <dune/tnnmg/solvers/scalartnnmg.hh>
+
+#include "laplace.hh"
+
+template <int dim, int N>
+class ConstantFunction :
+    public Dune::FunctionBase<double,double,dim,N>
+{
+    public:
+        ConstantFunction(double c):
+            c_(c)
+        {}
+
+        double eval(int comp, const Dune::FieldVector<double,dim>& x) const
+        {
+            return c_;
+        }
+
+        void evalall(const Dune::FieldVector<double,dim>& x, Dune::FieldVector<double,N>& y) const
+        {
+            y = eval(0,x);
+        }
+
+        double c_;
+};
+
+template <int dim, int N>
+class AbsFunction :
+    public Dune::FunctionBase<double,double,dim,N>
+{
+    public:
+        AbsFunction()
+        {}
+
+        double eval(int comp, const Dune::FieldVector<double,dim>& x) const
+        {
+            return fabs(x[0]);
+        }
+
+        void evalall(const Dune::FieldVector<double,dim>& x, Dune::FieldVector<double,N>& y) const
+        {
+            y = eval(0,x);
+        }
+};
+
+int main (int argc, char *argv[]) try
+{
+#ifdef EMBEDDED_PYTHON
+    EmbeddedPython::start();
+#endif
+
+    Dune::ConfigParser parset;
+
+    // parse command line once to obtain a possible parset argument
+    // parse parameter file
+    // parse command line a second time in order to allow overwriting parameters
+    parset.parseCmd(argc, argv);
+    parset.parseFile(parset.get("file", "adaptive_laplace.parset"));
+    parset.parseCmd(argc, argv);
+
+    // The grid dimension
+    const int dim = DIMENSION;
+
+    // define grid type
+#ifdef UGGRID
+    typedef Dune::UGGrid<dim> GridType;
+#endif
+#ifdef ALUGRID
+    typedef Dune::ALUSimplexGrid<dim,dim> GridType;
+#endif
+
+    typedef NamedFunctionMap<GridType::Codim<0>::Geometry::GlobalCoordinate,Dune::FieldVector<double,1> > Functions;
+    typedef Functions::Function Function;
+
+    Functions functions;
+    functions["one"] = new ConstantFunction<dim,1>(1.0);
+    functions["zero"] = new ConstantFunction<dim,1>(0.0);
+    functions["abs"] = new AbsFunction<dim,1>;
+
+    // Generate the grid
+    Dune::GridPtr<GridType> gridptr;
+    if (dim==2)
+        gridptr=Dune::GridPtr<GridType>(parset.get("grid.dgffile2d", "coarse.dgf").c_str());
+    else if (dim==3)
+        gridptr=Dune::GridPtr<GridType>(parset.get("grid.dgffile3d", "coarse.dgf").c_str());
+
+    GridType& grid = *gridptr;
+
+// grid specific settings
+#ifdef UGGRID
+    grid.setRefinementType(GridType::LOCAL);
+    grid.setClosureType(GridType::NONE);
+#endif
+#ifdef ALUGRID
+    if (parset.hasKey("grid.restore"))
+    {
+        double time;
+        grid.readGrid<Dune::xdr>(parset.get<std::string>("grid.backupfile"), time);
+        std::cout << "Grid restored. Time was " << time << std::endl;
+    }
+    else
+#endif
+    {
+        //initial refinement of grid
+        for (int i=0; i<parset.get("grid.refine", 0) or i<parset.get("laplace.refinement.minlevel",0); ++i)
+            grid.globalRefine(1);
+        std::cout << "Grid refined." << std::endl;
+    }
+
+#ifdef EMBEDDED_PYTHON
+    std::string importDictName = formatString(std::string("importAsDuneFunction%1dD"), dim);
+    if (parset.hasKey("python.inline"))
+    {
+        EmbeddedPython::run(parset.get("python.inline", ""));
+        EmbeddedPython::importFunctionsFromModule(functions, "__main__", importDictName);
+    }
+    if (parset.hasKey("python.module"))
+        EmbeddedPython::importFunctionsFromModule(functions, parset.get("python.module", "laplace"), importDictName);
+#endif
+
+    LaplaceProblem<GridType> laplaceProblem(grid, parset.sub("laplace"), functions);
+
+    bool refined;
+
+    do
+    {
+        laplaceProblem.assemble();
+        laplaceProblem.solve();
+        refined = laplaceProblem.adapt();
+    }
+    while(refined);
+
+
+#ifdef EMBEDDED_PYTHON
+    EmbeddedPython::stop();
+#endif
+
+
+    // delete created functions
+    Functions::iterator it = functions.begin();
+    Functions::iterator end = functions.end();
+    for(; it!=end; ++it)
+        delete it->second;
+
+    return 0;
+}
+
+catch (Dune::Exception e)
+{
+    std::cout << e << std::endl;
+}
+
diff --git a/doc/manual/laplace.hh b/doc/manual/laplace.hh
new file mode 100755
index 0000000000000000000000000000000000000000..008211706b0faf0beeaaecee3d08d6658f362897
--- /dev/null
+++ b/doc/manual/laplace.hh
@@ -0,0 +1,376 @@
+#ifndef TNNMG_EXAMPLE_LAPLACE_HH
+#define TNNMG_EXAMPLE_LAPLACE_HH
+
+#include <fstream>
+#include <iostream>
+#include <iostream>
+#include <vector>
+
+#include <tr1/memory>
+
+// dune includes ******************************************************
+#include <dune/common/parametertree.hh>
+#include <dune/common/shared_ptr.hh>
+
+#include <dune/grid/utility/grapedataioformattypes.hh>
+
+#include <dune/istl/bcrsmatrix.hh>
+
+#include <dune/grid/common/genericreferenceelements.hh>
+
+// dune-fufem includes ******************************************************
+#include <dune/fufem/boundarypatch.hh>
+#include <dune/fufem/prolongboundarypatch.hh>
+
+#if EMBEDDED_PYTHON
+    #include <dune/fufem/functions/pythonfunction.hh>
+#endif
+#include <dune/fufem/functions/functions.hh>
+#include <dune/fufem/functions/basisgridfunction.hh>
+
+#include <dune/fufem/functionspacebases/p1nodalbasis.hh>
+
+#include <dune/fufem/functiontools/namedfunctionmap.hh>
+#include <dune/fufem/functiontools/boundarydofs.hh>
+#include <dune/fufem/functiontools/gridfunctionadaptor.hh>
+
+#include <dune/fufem/assemblers/operatorassembler.hh>
+#include <dune/fufem/assemblers/functionalassembler.hh>
+#include <dune/fufem/assemblers/transferoperatorassembler.hh>
+
+#include <dune/fufem/assemblers/localassemblers/laplaceassembler.hh>
+#include <dune/fufem/assemblers/localassemblers/massassembler.hh>
+#include <dune/fufem/assemblers/localassemblers/lumpedmassassembler.hh>
+#include <dune/fufem/assemblers/localassemblers/l2functionalassembler.hh>
+
+#if HAVE_AMIRAMESH
+    #include <dune/fufem/functiontools/amirameshbasiswriter.hh>
+#endif
+#include <dune/fufem/functiontools/vtkbasiswriter.hh>
+
+// dune-solvers includes ***************************************************
+#include <dune/solvers/norms/norm.hh>
+#include <dune/solvers/norms/sumnorm.hh>
+#include <dune/solvers/norms/energynorm.hh>
+#include <dune/solvers/norms/diagnorm.hh>
+#include <dune/solvers/norms/fullnorm.hh>
+#include <dune/solvers/norms/blocknorm.hh>
+
+#include <dune/solvers/iterationsteps/blockgsstep.hh>
+#include <dune/solvers/iterationsteps/multigridstep.hh>
+
+#include <dune/solvers/solvers/loopsolver.hh>
+
+// dune-tnnmg includes *******************************************************
+//#include <dune/tnnmg/problem-classes/blocknonlineargsproblem.hh>
+#include <dune/tnnmg/problem-classes/blocknonlineartnnmgproblem.hh>
+#include <dune/tnnmg/problem-classes/nonlinearity.hh>
+#include <dune/tnnmg/nonlinearities/zerononlinearity.hh>
+#include <dune/tnnmg/iterationsteps/genericnonlineargs.hh>
+#include <dune/tnnmg/iterationsteps/tnnmgstep.hh>
+
+template<class T1>
+std::string formatString(const std::string format, const T1& t1)
+{
+    char buffer[1000];
+    sprintf(buffer, format.c_str(), t1);
+    return std::string(buffer);
+}
+
+template<class T1, class T2>
+std::string formatString(const std::string format, const T1& t1, const T2& t2)
+{
+    char buffer[1000];
+    sprintf(buffer, format.c_str(), t1, t2);
+    return std::string(buffer);
+}
+
+template <class GridType>
+class LaplaceProblem
+{
+    public:
+        static const int dim = GridType::dimension;
+        static const int block_size = 1;
+
+        typedef Dune::FieldMatrix<double,block_size,block_size> FineLocalMatrix;
+        typedef Dune::FieldVector<double,block_size> FineLocalVector;
+
+        typedef Dune::BCRSMatrix< FineLocalMatrix > FineMatrixType;
+        typedef Dune::BlockVector< FineLocalVector > FineVectorType;
+
+        /* in the linear Laplace problem there is no nonlinearity ;-) */
+        typedef ZeroNonlinearity<FineLocalVector,FineLocalMatrix> NonlinearityType;
+
+        /* we set the problemtypes:
+           the ConvexProblem depends on the NonlinearityType and the type of the fine grid matrix*/
+        typedef ConvexProblem<NonlinearityType, FineMatrixType> ConvexProblemType;
+        typedef BlockNonlinearTNNMGProblem<ConvexProblemType> TNNMGProblemType;
+
+        typedef typename TNNMGProblemType::Linearization::VectorType CoarseVectorType;
+        typedef typename TNNMGProblemType::Linearization::MatrixType CoarseMatrixType;
+
+        /* use a linear Gauss-Seidel as fine grid smoother */
+        typedef GenericNonlinearGS<TNNMGProblemType> FineSmoother;
+        typedef TruncatedNonsmoothNewtonMultigrid<TNNMGProblemType, FineSmoother> TNNMGType;
+
+        /* define the linear coarse grid solver type:
+           use linear Gauss-Seidel as smoother for linear multigrid */
+        typedef BlockGSStep<CoarseMatrixType,CoarseVectorType> CoarseSmoother;
+
+        typedef ScalarTNNMG::Transfer Transfer;
+
+        typedef NamedFunctionMap<typename GridType::template Codim<0>::Geometry::GlobalCoordinate,Dune::FieldVector<double,1> > FunctionMap;
+        typedef typename FunctionMap::Function Function;
+
+        typedef P1NodalBasis<typename GridType::LeafGridView> NonconformingBasis;
+        typedef ConformingBasis<NonconformingBasis> Basis;
+
+        LaplaceProblem(GridType& grid, Dune::ParameterTree& parset, const FunctionMap& functions) :
+            grid_(grid),
+            parset_(parset),
+            functions_(functions),
+            ncBasis_(grid.leafView()),
+            basis_(ncBasis_)
+        {
+            // create energy norm
+            norm_ = NormPointer(new EnergyNorm<FineMatrixType,FineVectorType>(stiffnessMatrix_));
+
+            // create a zero- nonlinearity, i.e. no nonlinearity ;-)
+            phi_ = NonlinearityPointer(new NonlinearityType());
+
+            // initialize solution vector
+            u_.resize(basis_.size());
+            u_ = 0.0;
+
+            refinementStep_ = 0;
+
+            typedef  BoundaryPatchSegmentIndexProperty<typename GridType::LevelGridView> Property;
+
+            coarseBoundaryPatch_.setup(grid.levelView(0));
+            coarseBoundaryPatch_.insertFacesByProperty(Property(1));
+        }
+
+        ~LaplaceProblem()
+        {
+            for(int l=0; l<transfer_.size(); ++l)
+                delete transfer_[l];
+
+            delete norm_;
+        }
+
+        template <class AssembleBasis>
+        void assembleRHS(const AssembleBasis& basis, const int rhsIntegrationOrder)
+        {
+            FunctionalAssembler<AssembleBasis> assembler(basis);
+
+            {
+                std::cout << "assembling right hand side" << std::endl;
+                const Function& rhsFunction = *functions_.getFunction(parset_.get("rhs.function", "one"));
+                L2FunctionalAssembler<GridType> l2Assembler(rhsFunction, rhsIntegrationOrder);
+                assembler.assemble(l2Assembler, rhs_);
+            }
+
+            {
+                std::cout << "assembling mass weight vector" << std::endl;
+                L2FunctionalAssembler<GridType> l2Assembler(*functions_.getFunction("one"), 1);
+                assembler.assemble(l2Assembler, massWeights_);
+            }
+        }
+
+        // assemble the whole problem
+        void assemble()
+        {
+            // setup boundary patch and boundary dofs
+            PatchProlongator<GridType>::prolong(coarseBoundaryPatch_, boundaryPatch_);
+
+            // setup boundary dofs
+            constructBoundaryDofs(boundaryPatch_, basis_, isBoundary_);
+
+            // assemble right hand side with given integration order
+            assembleRHS(basis_, parset_.get("rhs.integrationorder", 2));
+
+            // create a global operator assembler to be used for the matrices
+            OperatorAssembler<Basis, Basis> assembler(basis_, basis_);
+            typedef typename Basis::LocalFiniteElement FE;
+
+            {
+                // assemble stiffness matrix
+                std::cout << "assembling stiffness matrix" << std::endl;
+                LaplaceAssembler<GridType,FE,FE> laplaceAssembler;
+                assembler.assemble(laplaceAssembler, stiffnessMatrix_);
+            }
+
+            if (parset_.get("masslumping",true))
+            {
+                // assemble lumped mass matrix
+                std::cout << "assembling lumped mass matrix" << std::endl;
+                LumpedMassAssembler<GridType,FE,FE> lumpedMassAssembler;
+                assembler.assemble(lumpedMassAssembler, massMatrix_, true);
+            }
+            else
+            {
+                // assemble mass matrix
+                std::cout << "assembling mass matrix" << std::endl;
+                MassAssembler<GridType,FE,FE> massAssembler;
+                assembler.assemble(massAssembler, massMatrix_);
+            }
+
+            // multiply matrices by given coefficients and add them
+            stiffnessMatrix_ *= parset_.get("stiffnessfactor", 1.0);
+            massMatrix_ *= parset_.get("massfactor", 1.0);
+            stiffnessMatrix_ += massMatrix_;
+
+            // delete old transfer operators
+            for(int l=0; l<transfer_.size(); ++l)
+                delete transfer_[l];
+
+            // create new transfer operators
+            transfer_.resize(grid_.maxLevel());
+            for(int l=0; l<transfer_.size(); ++l)
+                transfer_[l] = new Transfer;
+
+            if (parset_.get("refinement.type", "local") == "global")
+            {
+                // assemble transfer operators on uniformly refined grid
+                // this is only for demonstration reasons
+                // the all at once aasembler below does also
+                // work for uniformly refined grids
+                for(int l=0; l<transfer_.size(); ++l)
+                    transfer_[l]->setup(grid_, l,l+1);
+            }
+            else
+            {
+                // assemble transfer operator hierarchy on adaptively refined grid
+                TransferOperatorAssembler<GridType> transferOperatorAssembler(grid_);
+                transferOperatorAssembler.assembleOperatorPointerHierarchy(transfer_);
+            }
+
+            // compute entries of solution vector for dirichlet nodes
+            // as interpolation of given dirichlet function
+            Functions::interpolate(basis_, u_, *functions_.getFunction(parset_.get("dirichletfunction", "zero")),isBoundary_);
+//            GenericGridFunctionInterpolator<Basis,FineVectorType,Function>::interpolate(basis_, u_, *functions_.getFunction(parset_.get("dirichletfunction", "zero")), isBoundary_);
+        }
+
+        void solve()
+        {
+            std::cout << "u.size(): " << u_.size() << std::endl;
+            ConvexProblemType P(1.0, stiffnessMatrix_, 0.0, massWeights_, *phi_, rhs_, u_);
+            TNNMGProblemType tnnmgProblem(parset_.sub("tnnmg"), P);
+
+            /* the "nonlinear" or rather fine grid smoother */
+            FineSmoother fine_smoother;
+            /* the linear or coarse smoothers */
+            CoarseSmoother presmoother, postsmoother;
+            /* the linear (coarse grid) iterative solver step */
+            MultigridStep<CoarseMatrixType, CoarseVectorType> mgStep;
+            {
+                mgStep.setTransferOperators(transfer_);
+                mgStep.setNumberOfLevels(mgStep.mgTransfer_.size()+1);
+
+                mgStep.setSmoother(&presmoother, &postsmoother);
+
+                /* use one simple smoother iteration as base solver */
+                CoarseSmoother* baseSolverStep = new CoarseSmoother;
+                EnergyNorm<CoarseMatrixType, CoarseVectorType>* baseEnergyNorm = new EnergyNorm<CoarseMatrixType, CoarseVectorType>(*baseSolverStep);
+                mgStep.basesolver_ = new LoopSolver<CoarseVectorType>(baseSolverStep,
+                                                                        1,
+                                                                        parset_.get("mmg.tol", 1e-12),
+                                                                        baseEnergyNorm,
+                                                                        Solver::QUIET);
+                mgStep.setMGType(1, 3, 3);
+                mgStep.setNumberOfLevels(transfer_.size()+1);
+            }
+
+            Dune::BitSetVector<1> ignoreNodes(isBoundary_);
+            for (int i=0; i<ignoreNodes.size(); ++i)
+            {
+                if (basis_.isConstrained(i))
+                    ignoreNodes[i].set();
+            }
+
+            TNNMGType tnnmgstep(mgStep, fine_smoother);
+            tnnmgstep.setProblem(P.u, tnnmgProblem);
+
+            tnnmgstep.ignoreNodes_ = &ignoreNodes;
+            tnnmgstep.setSmoothingSteps(1, 1, 0);
+
+            LoopSolver<VectorType> tnnmgsolver(&tnnmgstep,
+                                           parset_.get("tnnmg.maxiter", 100),
+                                           parset_.get("tnnmg.tol", 1e-12),
+                                           norm_,
+                                           Solver::FULL);
+
+            tnnmgsolver.preprocess();
+            tnnmgsolver.solve();
+            u_ = tnnmgstep.getSol();
+        }
+
+        void writeSolution()
+        {
+#if HAVE_AMIRAMESH
+            if (parset_.get("solutionfile_am", "") != "")
+            {
+                AmiraMeshBasisWriter<Basis> writer(basis_);
+                writer.addP1Interpolation(u_);
+                writer.write(formatString(parset_.get("solutionfile_am", ""), refinementStep_));
+            }
+#endif
+            if (parset_.get("solutionfile_vtk", "") != "")
+            {
+                VTKBasisWriter<Basis> writer(basis_);
+                writer.addP1Interpolation(u_);
+                writer.write(formatString(parset_.get("solutionfile_vtk", ""), refinementStep_));
+            }
+        }
+
+        bool adapt()
+        {
+            if (parset_.get("refinement.maxlevel",8) > grid_.maxLevel())
+            {
+                GridFunctionAdaptor<Basis> adaptor(basis_, true);
+                
+                grid_.globalRefine(1);
+                
+                ncBasis_.update(grid_.leafView());
+                basis_.update(grid_.leafView());
+
+                adaptor.adapt(u_);
+                return true;
+            }
+            return false;
+        }
+
+    private:
+
+//        const Dune::ParameterTree& parset_;
+        Dune::ParameterTree& parset_;
+        GridType& grid_;
+        const FunctionMap& functions_;
+        int refinementStep_;
+
+        NonconformingBasis ncBasis_;
+        Basis basis_;
+
+        BoundaryPatchBase< typename Basis::GridView > boundaryPatch_;
+
+        BoundaryPatchBase< typename GridType::LevelGridView > coarseBoundaryPatch_;
+
+        FineMatrixType stiffnessMatrix_;
+        FineMatrixType massMatrix_;
+
+        FineVectorType u_;
+        FineVectorType rhs_;
+        FineVectorType massWeights_;
+
+        typedef typename Dune::shared_ptr< NonlinearityType > NonlinearityPointer;
+        NonlinearityPointer phi_;
+
+        Dune::BitSetVector<1> isBoundary_;
+
+        typedef Norm<FineVectorType>* NormPointer;
+        NormPointer norm_;
+
+        std::vector<Transfer*> transfer_;
+};
+#endif
+
diff --git a/doc/manual/src/cube_centered.dgf b/doc/manual/src/cube_centered.dgf
new file mode 100755
index 0000000000000000000000000000000000000000..9f8108798c9e7de87b0c037c5c01f3f4cf71a811
--- /dev/null
+++ b/doc/manual/src/cube_centered.dgf
@@ -0,0 +1,25 @@
+DGF
+Vertex        % the verticies of the grid
+-1   -1    1       % vertex 0
+ 1   -1    1       % vertex 1
+-1   -1   -1       % vertex 4
+ 1    1    1       % vertex 3
+ 
+-1    1    1       % vertex 2
+ 1   -1   -1       % vertex 5
+-1    1   -1       % vertex 6
+ 1    1   -1       % vertex 7
+
+#
+SIMPLEX       % a simplex grid
+0 1 4 2
+2 5 6 4
+1 2 5 4
+
+1 3 4 7
+5 6 7 4
+1 5 7 4
+#
+BOUNDARYDOMAIN
+default 1     % all other boundary segments have id 1
+#
\ No newline at end of file
diff --git a/doc/manual/src/dune_tnnmg_how2.cc b/doc/manual/src/dune_tnnmg_how2.cc
new file mode 100644
index 0000000000000000000000000000000000000000..848e50d0128ba311b499ceb07f4430988bf1065a
--- /dev/null
+++ b/doc/manual/src/dune_tnnmg_how2.cc
@@ -0,0 +1,27 @@
+#ifdef HAVE_CONFIG_H
+# include "config.h"     
+#endif
+#include <iostream>
+#include"dune/common/mpihelper.hh" // An initializer of MPI
+#include"dune/common/exceptions.hh" // We use exceptions
+
+int main(int argc, char** argv)
+{
+  try{
+    //Maybe initialize Mpi
+    Dune::MPIHelper& helper = Dune::MPIHelper::instance(argc, argv);
+    std::cout << "Hello World! This is dune-tnnmg-how2." << std::endl;
+    if(Dune::MPIHelper::isFake)
+      std::cout<< "This is a sequential program." << std::endl;
+    else
+      std::cout<<"I am rank "<<helper.rank()<<" of "<<helper.size()
+        <<" processes!"<<std::endl;
+    return 0;
+  }
+  catch (Dune::Exception &e){
+    std::cerr << "Dune reported error: " << e << std::endl;
+  }
+  catch (...){
+    std::cerr << "Unknown exception thrown!" << std::endl;
+  }
+}
diff --git a/doc/manual/src/examples.cc b/doc/manual/src/examples.cc
new file mode 100644
index 0000000000000000000000000000000000000000..0646df693d06280497618fc1b78c3712c3d12a99
--- /dev/null
+++ b/doc/manual/src/examples.cc
@@ -0,0 +1,240 @@
+#include <config.h>
+
+#define DIMENSION 2
+#define ALUGRID 1
+//#define UGGRID 1
+#define EMBEDDED_PYTHON 1
+#undef EMBEDDED_PYTHON
+
+#define FE_VERBOSE
+
+// disable embedded python if python was not found
+#ifndef HAVE_PYTHON
+    #undef EMBEDDED_PYTHON
+#endif
+
+#if EMBEDDED_PYTHON
+    #include <Python.h>
+#endif
+
+#include <fstream>
+#include <iostream>
+#include <iostream>
+#include <vector>
+
+#include <tr1/memory>
+
+// dune includes ******************************************************
+#ifdef UGGRID
+    #include <dune/grid/uggrid.hh>
+    #include <dune/grid/io/file/dgfparser/dgfug.hh>
+#endif
+
+#ifdef ALUGRID
+    #include <dune/grid/alugrid.hh>
+    #include <dune/grid/io/file/dgfparser/dgfalu.hh>
+#endif
+
+#include <dune/common/parametertree.hh>
+#include <dune/common/parametertreeparser.hh>
+#include <dune/common/shared_ptr.hh>
+
+#include <dune/istl/bcrsmatrix.hh>
+
+#include <dune/grid/common/genericreferenceelements.hh>
+
+// dune-fufem includes ******************************************************
+#include <dune/fufem/boundarypatch.hh>
+#include <dune/fufem/prolongboundarypatch.hh>
+
+#if EMBEDDED_PYTHON
+    #include <dune/fufem/functions/pythonfunction.hh>
+#endif
+#include <dune/fufem/functions/basisgridfunction.hh>
+
+#include <dune/fufem/functionspacebases/p1nodalbasis.hh>
+
+#include <dune/fufem/functiontools/namedfunctionmap.hh>
+#include <dune/fufem/functiontools/boundarydofs.hh>
+
+#include <dune/fufem/assemblers/operatorassembler.hh>
+#include <dune/fufem/assemblers/functionalassembler.hh>
+#include <dune/fufem/assemblers/transferoperatorassembler.hh>
+
+#include <dune/fufem/assemblers/localassemblers/laplaceassembler.hh>
+#include <dune/fufem/assemblers/localassemblers/massassembler.hh>
+#include <dune/fufem/assemblers/localassemblers/lumpedmassassembler.hh>
+#include <dune/fufem/assemblers/localassemblers/l2functionalassembler.hh>
+
+#if HAVE_AMIRAMESH
+    #include <dune/fufem/functiontools/amirameshbasiswriter.hh>
+#endif
+#include <dune/fufem/functiontools/vtkbasiswriter.hh>
+
+// dune-solvers includes ***************************************************
+#include <dune/solvers/norms/norm.hh>
+#include <dune/solvers/norms/sumnorm.hh>
+#include <dune/solvers/norms/energynorm.hh>
+#include <dune/solvers/norms/diagnorm.hh>
+#include <dune/solvers/norms/fullnorm.hh>
+#include <dune/solvers/norms/blocknorm.hh>
+
+// dune-tnnmg includes *******************************************************
+#include <dune/tnnmg/problem-classes/blocknonlineargsproblem.hh>
+#include <dune/tnnmg/problem-classes/nonlinearity.hh>
+#include <dune/tnnmg/nonlinearities/zerononlinearity.hh>
+#include <dune/tnnmg/nonlinearities/shiftednonlinearity.hh>
+#include <dune/tnnmg/iterationsteps/genericnonlinearjacobi.hh>
+#include <dune/tnnmg/solvers/scalartnnmg.hh>
+
+#include "laplace.hh"
+
+template <int dim, int N>
+class ConstantFunction :
+    public Dune::FunctionBase<double,double,dim,N>
+{
+    public:
+        ConstantFunction(double c):
+            c_(c)
+        {}
+
+        double eval(int comp, const Dune::FieldVector<double,dim>& x) const
+        {
+            return c_;
+        }
+
+        void evalall(const Dune::FieldVector<double,dim>& x, Dune::FieldVector<double,N>& y) const
+        {
+            y = eval(0,x);
+        }
+
+        double c_;
+};
+
+template <int dim, int N>
+class AbsFunction :
+    public Dune::FunctionBase<double,double,dim,N>
+{
+    public:
+        AbsFunction()
+        {}
+
+        double eval(int comp, const Dune::FieldVector<double,dim>& x) const
+        {
+            return fabs(x[0]);
+        }
+
+        void evalall(const Dune::FieldVector<double,dim>& x, Dune::FieldVector<double,N>& y) const
+        {
+            y = eval(0,x);
+        }
+};
+
+int main (int argc, char *argv[]) try
+{
+#ifdef EMBEDDED_PYTHON
+    EmbeddedPython::start();
+#endif
+
+    Dune::ParameterTree parset;
+
+    // parse command line once to obtain a possible parset argument
+    // parse parameter file
+    // parse command line a second time in order to allow overwriting parameters
+    Dune::ParameterTreeParser::readOptions(argc, argv,parset);
+    Dune::ParameterTreeParser::readINITree(parset.get("file", "adaptive_laplace.parset"),parset);
+    Dune::ParameterTreeParser::readOptions(argc, argv,parset);
+
+    // The grid dimension
+    const int dim = DIMENSION;
+
+    // define grid type
+#ifdef UGGRID
+    typedef Dune::UGGrid<dim> GridType;
+#endif
+#ifdef ALUGRID
+    typedef Dune::ALUSimplexGrid<dim,dim> GridType;
+#endif
+
+    typedef NamedFunctionMap<GridType::Codim<0>::Geometry::GlobalCoordinate,Dune::FieldVector<double,1> > Functions;
+    typedef Functions::Function Function;
+
+    Functions functions;
+    functions["one"] = new ConstantFunction<dim,1>(1.0);
+    functions["zero"] = new ConstantFunction<dim,1>(0.0);
+    functions["abs"] = new AbsFunction<dim,1>;
+
+    // Generate the grid
+    Dune::GridPtr<GridType> gridptr;
+    if (dim==2)
+        gridptr=Dune::GridPtr<GridType>(parset.get("grid.dgffile2d", "coarse.dgf").c_str());
+    else if (dim==3)
+        gridptr=Dune::GridPtr<GridType>(parset.get("grid.dgffile3d", "coarse.dgf").c_str());
+
+    GridType& grid = *gridptr;
+
+// grid specific settings
+#ifdef UGGRID
+    grid.setRefinementType(GridType::LOCAL);
+    grid.setClosureType(GridType::NONE);
+#endif
+#ifdef ALUGRID
+    if (parset.hasKey("grid.restore"))
+    {
+        double time;
+        grid.readGrid<Dune::xdr>(parset.get<std::string>("grid.backupfile"), time);
+        std::cout << "Grid restored. Time was " << time << std::endl;
+    }
+    else
+#endif
+    {
+        //initial refinement of grid
+        for (int i=0; i<parset.get("grid.refine", 0) or i<parset.get("laplace.refinement.minlevel",0); ++i)
+            grid.globalRefine(1);
+        std::cout << "Grid refined." << std::endl;
+    }
+
+#ifdef EMBEDDED_PYTHON
+    std::string importDictName = formatString(std::string("importAsDuneFunction%1dD"), dim);
+    if (parset.hasKey("python.inline"))
+    {
+        EmbeddedPython::run(parset.get("python.inline", ""));
+        EmbeddedPython::importFunctionsFromModule(functions, "__main__", importDictName);
+    }
+    if (parset.hasKey("python.module"))
+        EmbeddedPython::importFunctionsFromModule(functions, parset.get("python.module", "laplace"), importDictName);
+#endif
+
+    LaplaceProblem<GridType> laplaceProblem(grid, parset.sub("laplace"), functions);
+
+    bool refined;
+
+    do
+    {
+        laplaceProblem.assemble();
+        laplaceProblem.solve();
+        refined = laplaceProblem.adapt();
+        laplaceProblem.writeSolution();
+    }
+    while(refined);
+
+
+#ifdef EMBEDDED_PYTHON
+    EmbeddedPython::stop();
+#endif
+
+
+    // delete created functions
+    Functions::iterator it = functions.begin();
+    Functions::iterator end = functions.end();
+    for(; it!=end; ++it)
+        delete it->second;
+
+    return 0;
+}
+
+catch (Dune::Exception e)
+{
+    std::cout << e << std::endl;
+}
+
diff --git a/doc/manual/src/laplace.hh b/doc/manual/src/laplace.hh
new file mode 100644
index 0000000000000000000000000000000000000000..d08ba536a490b4a9dbe5fc382e576ea58a4c22da
--- /dev/null
+++ b/doc/manual/src/laplace.hh
@@ -0,0 +1,378 @@
+#ifndef TNNMG_EXAMPLE_LAPLACE_HH
+#define TNNMG_EXAMPLE_LAPLACE_HH
+
+#include <fstream>
+#include <iostream>
+#include <iostream>
+#include <vector>
+
+#include <tr1/memory>
+
+// dune includes ******************************************************
+#include <dune/common/parametertree.hh>
+#include <dune/common/shared_ptr.hh>
+
+#include <dune/grid/utility/grapedataioformattypes.hh>
+
+#include <dune/istl/bcrsmatrix.hh>
+
+#include <dune/grid/common/genericreferenceelements.hh>
+
+// dune-fufem includes ******************************************************
+#include <dune/fufem/boundarypatch.hh>
+#include <dune/fufem/prolongboundarypatch.hh>
+
+#if EMBEDDED_PYTHON
+    #include <dune/fufem/functions/pythonfunction.hh>
+#endif
+#include <dune/fufem/functions/functions.hh>
+#include <dune/fufem/functions/basisgridfunction.hh>
+
+#include <dune/fufem/functionspacebases/p1nodalbasis.hh>
+
+#include <dune/fufem/functiontools/namedfunctionmap.hh>
+#include <dune/fufem/functiontools/boundarydofs.hh>
+#include <dune/fufem/functiontools/gridfunctionadaptor.hh>
+
+#include <dune/fufem/assemblers/operatorassembler.hh>
+#include <dune/fufem/assemblers/functionalassembler.hh>
+#include <dune/fufem/assemblers/transferoperatorassembler.hh>
+
+#include <dune/fufem/assemblers/localassemblers/laplaceassembler.hh>
+#include <dune/fufem/assemblers/localassemblers/massassembler.hh>
+#include <dune/fufem/assemblers/localassemblers/lumpedmassassembler.hh>
+#include <dune/fufem/assemblers/localassemblers/l2functionalassembler.hh>
+
+#if HAVE_AMIRAMESH
+    #include <dune/fufem/functiontools/amirameshbasiswriter.hh>
+#endif
+#include <dune/fufem/functiontools/vtkbasiswriter.hh>
+
+// dune-solvers includes ***************************************************
+#include <dune/solvers/norms/norm.hh>
+#include <dune/solvers/norms/sumnorm.hh>
+#include <dune/solvers/norms/energynorm.hh>
+#include <dune/solvers/norms/diagnorm.hh>
+#include <dune/solvers/norms/fullnorm.hh>
+#include <dune/solvers/norms/blocknorm.hh>
+
+#include <dune/solvers/iterationsteps/blockgsstep.hh>
+#include <dune/solvers/iterationsteps/multigridstep.hh>
+
+#include <dune/solvers/solvers/loopsolver.hh>
+
+// dune-tnnmg includes *******************************************************
+//#include <dune/tnnmg/problem-classes/blocknonlineargsproblem.hh>
+#include <dune/tnnmg/problem-classes/blocknonlineartnnmgproblem.hh>
+#include <dune/tnnmg/problem-classes/nonlinearity.hh>
+#include <dune/tnnmg/nonlinearities/zerononlinearity.hh>
+#include <dune/tnnmg/iterationsteps/genericnonlineargs.hh>
+#include <dune/tnnmg/iterationsteps/tnnmgstep.hh>
+
+template<class T1>
+std::string formatString(const std::string format, const T1& t1)
+{
+    char buffer[1000];
+    sprintf(buffer, format.c_str(), t1);
+    return std::string(buffer);
+}
+
+template<class T1, class T2>
+std::string formatString(const std::string format, const T1& t1, const T2& t2)
+{
+    char buffer[1000];
+    sprintf(buffer, format.c_str(), t1, t2);
+    return std::string(buffer);
+}
+
+template <class GridType>
+class LaplaceProblem
+{
+    public:
+        static const int dim = GridType::dimension;
+        static const int block_size = 1;
+
+        typedef Dune::FieldMatrix<double,block_size,block_size> FineLocalMatrix;
+        typedef Dune::FieldVector<double,block_size> FineLocalVector;
+
+        typedef Dune::BCRSMatrix< FineLocalMatrix > FineMatrixType;
+        typedef Dune::BlockVector< FineLocalVector > FineVectorType;
+
+        /* in the linear Laplace problem there is no nonlinearity ;-) */
+        typedef ZeroNonlinearity<FineLocalVector,FineLocalMatrix> NonlinearityType;
+
+        /* we set the problemtypes:
+           the ConvexProblem depends on the NonlinearityType and the type of the fine grid matrix*/
+        typedef ConvexProblem<NonlinearityType, FineMatrixType> ConvexProblemType;
+        typedef BlockNonlinearTNNMGProblem<ConvexProblemType> TNNMGProblemType;
+
+        typedef typename TNNMGProblemType::Linearization::VectorType CoarseVectorType;
+        typedef typename TNNMGProblemType::Linearization::MatrixType CoarseMatrixType;
+
+        /* use a linear Gauss-Seidel as fine grid smoother */
+        typedef GenericNonlinearGS<TNNMGProblemType> FineSmoother;
+        typedef TruncatedNonsmoothNewtonMultigrid<TNNMGProblemType, FineSmoother> TNNMGType;
+
+        /* define the linear coarse grid solver type:
+           use linear Gauss-Seidel as smoother for linear multigrid */
+        typedef BlockGSStep<CoarseMatrixType,CoarseVectorType> CoarseSmoother;
+
+        typedef ScalarTNNMG::Transfer Transfer;
+
+        typedef NamedFunctionMap<typename GridType::template Codim<0>::Geometry::GlobalCoordinate,Dune::FieldVector<double,1> > FunctionMap;
+        typedef typename FunctionMap::Function Function;
+
+        typedef P1NodalBasis<typename GridType::LeafGridView> NonconformingBasis;
+        typedef ConformingBasis<NonconformingBasis> Basis;
+
+        LaplaceProblem(GridType& grid, Dune::ParameterTree& parset, const FunctionMap& functions) :
+            grid_(grid),
+            parset_(parset),
+            functions_(functions),
+            ncBasis_(grid.leafView()),
+            basis_(ncBasis_),
+            coarseBoundaryPatch_(grid.levelView(0),true)
+        {
+            // create energy norm
+            norm_ = NormPointer(new EnergyNorm<FineMatrixType,FineVectorType>(stiffnessMatrix_));
+
+            // create a zero- nonlinearity, i.e. no nonlinearity ;-)
+            phi_ = NonlinearityPointer(new NonlinearityType());
+
+            // initialize solution vector
+            u_.resize(basis_.size());
+            u_ = 0.0;
+
+            refinementStep_ = 0;
+
+//            typedef  BoundaryPatchSegmentIndexProperty<typename GridType::LevelGridView> Property;
+
+//            coarseBoundaryPatch_.setup(grid.levelView(0));
+//            coarseBoundaryPatch_.insertFacesByProperty(Property(1));
+        }
+
+        ~LaplaceProblem()
+        {
+            for(int l=0; l<transfer_.size(); ++l)
+                delete transfer_[l];
+
+            delete norm_;
+        }
+
+        template <class AssembleBasis>
+        void assembleRHS(const AssembleBasis& basis, const int rhsIntegrationOrder)
+        {
+            FunctionalAssembler<AssembleBasis> assembler(basis);
+
+            {
+                std::cout << "assembling right hand side" << std::endl;
+                const Function& rhsFunction = *functions_.getFunction(parset_.get("rhs.function", "one"));
+                L2FunctionalAssembler<GridType> l2Assembler(rhsFunction, rhsIntegrationOrder);
+                assembler.assemble(l2Assembler, rhs_);
+            }
+
+            {
+                std::cout << "assembling mass weight vector" << std::endl;
+                L2FunctionalAssembler<GridType> l2Assembler(*functions_.getFunction("one"), 1);
+                assembler.assemble(l2Assembler, massWeights_);
+            }
+        }
+
+        // assemble the whole problem
+        void assemble()
+        {
+            // setup boundary patch and boundary dofs
+            PatchProlongator<GridType>::prolong(coarseBoundaryPatch_, boundaryPatch_);
+
+            // setup boundary dofs
+            constructBoundaryDofs(boundaryPatch_, basis_, isBoundary_);
+            std::cout << "no of boundary nodes: " << isBoundary_.count() << std::endl;
+
+            // assemble right hand side with given integration order
+            assembleRHS(basis_, parset_.get("rhs.integrationorder", 2));
+
+            // create a global operator assembler to be used for the matrices
+            OperatorAssembler<Basis, Basis> assembler(basis_, basis_);
+            typedef typename Basis::LocalFiniteElement FE;
+
+            {
+                // assemble stiffness matrix
+                std::cout << "assembling stiffness matrix" << std::endl;
+                LaplaceAssembler<GridType,FE,FE> laplaceAssembler;
+                assembler.assemble(laplaceAssembler, stiffnessMatrix_);
+            }
+
+            if (parset_.get("masslumping",true))
+            {
+                // assemble lumped mass matrix
+                std::cout << "assembling lumped mass matrix" << std::endl;
+                LumpedMassAssembler<GridType,FE,FE> lumpedMassAssembler;
+                assembler.assemble(lumpedMassAssembler, massMatrix_, true);
+            }
+            else
+            {
+                // assemble mass matrix
+                std::cout << "assembling mass matrix" << std::endl;
+                MassAssembler<GridType,FE,FE> massAssembler;
+                assembler.assemble(massAssembler, massMatrix_);
+            }
+
+            // multiply matrices by given coefficients and add them
+            stiffnessMatrix_ *= parset_.get("stiffnessfactor", 1.0);
+            massMatrix_ *= parset_.get("massfactor", 1.0);
+            stiffnessMatrix_ += massMatrix_;
+
+            // delete old transfer operators
+            for(int l=0; l<transfer_.size(); ++l)
+                delete transfer_[l];
+
+            // create new transfer operators
+            transfer_.resize(grid_.maxLevel());
+            for(int l=0; l<transfer_.size(); ++l)
+                transfer_[l] = new Transfer;
+
+            if (parset_.get("refinement.type", "local") == "global")
+            {
+                // assemble transfer operators on uniformly refined grid
+                // this is only for demonstration reasons
+                // the all at once aasembler below does also
+                // work for uniformly refined grids
+                for(int l=0; l<transfer_.size(); ++l)
+                    transfer_[l]->setup(grid_, l,l+1);
+            }
+            else
+            {
+                // assemble transfer operator hierarchy on adaptively refined grid
+                TransferOperatorAssembler<GridType> transferOperatorAssembler(grid_);
+                transferOperatorAssembler.assembleOperatorPointerHierarchy(transfer_);
+            }
+
+            // compute entries of solution vector for dirichlet nodes
+            // as interpolation of given dirichlet function
+            Functions::interpolate(basis_, u_, *functions_.getFunction(parset_.get("dirichletfunction", "zero")),isBoundary_);
+//            GenericGridFunctionInterpolator<Basis,FineVectorType,Function>::interpolate(basis_, u_, *functions_.getFunction(parset_.get("dirichletfunction", "zero")), isBoundary_);
+        }
+
+        void solve()
+        {
+            std::cout << "u.size(): " << u_.size() << std::endl;
+            ConvexProblemType P(1.0, stiffnessMatrix_, 0.0, massWeights_, *phi_, rhs_, u_);
+            TNNMGProblemType tnnmgProblem(parset_.sub("tnnmg"), P);
+
+            /* the "nonlinear" or rather fine grid smoother */
+            FineSmoother fine_smoother;
+            /* the linear or coarse smoothers */
+            CoarseSmoother presmoother, postsmoother;
+            /* the linear (coarse grid) iterative solver step */
+            MultigridStep<CoarseMatrixType, CoarseVectorType> coarseSolverStep;
+            {
+                coarseSolverStep.setTransferOperators(transfer_);
+                coarseSolverStep.setNumberOfLevels(coarseSolverStep.mgTransfer_.size()+1);
+
+                coarseSolverStep.setSmoother(&presmoother, &postsmoother);
+
+                /* use one simple smoother iteration as base solver */
+                CoarseSmoother* baseSolverStep = new CoarseSmoother;
+                EnergyNorm<CoarseMatrixType, CoarseVectorType>* baseEnergyNorm = new EnergyNorm<CoarseMatrixType, CoarseVectorType>(*baseSolverStep);
+                coarseSolverStep.basesolver_ = new LoopSolver<CoarseVectorType>(baseSolverStep,
+                                                                        1,
+                                                                        parset_.get("mmg.tol", 1e-12),
+                                                                        baseEnergyNorm,
+                                                                        Solver::QUIET);
+                coarseSolverStep.setMGType(1, 3, 3);
+                coarseSolverStep.setNumberOfLevels(transfer_.size()+1);
+            }
+
+            Dune::BitSetVector<1> ignoreNodes(isBoundary_);
+            for (int i=0; i<ignoreNodes.size(); ++i)
+            {
+                if (basis_.isConstrained(i))
+                    ignoreNodes[i].set();
+            }
+
+            TNNMGType tnnmgstep(coarseSolverStep, fine_smoother);
+            tnnmgstep.setProblem(P.u, tnnmgProblem);
+
+            tnnmgstep.ignoreNodes_ = &ignoreNodes;
+            tnnmgstep.setSmoothingSteps(1, 1, 0);
+
+            LoopSolver<VectorType> tnnmgsolver(&tnnmgstep,
+                                           parset_.get("tnnmg.maxiter", 100),
+                                           parset_.get("tnnmg.tol", 1e-12),
+                                           norm_,
+                                           Solver::FULL);
+
+            tnnmgsolver.preprocess();
+            tnnmgsolver.solve();
+            u_ = tnnmgstep.getSol();
+        }
+
+        void writeSolution()
+        {
+#if HAVE_AMIRAMESH
+            if (parset_.get("solutionfile_am", "") != "")
+            {
+                AmiraMeshBasisWriter<Basis> writer(basis_);
+                writer.addP1Interpolation(u_);
+                writer.write(formatString(parset_.get("solutionfile_am", ""), refinementStep_));
+            }
+#endif
+            if (parset_.get("solutionfile_vtk", "") != "")
+            {
+                VTKBasisWriter<Basis> writer(basis_);
+                writer.addP1Interpolation(u_);
+                writer.write(formatString(parset_.get("solutionfile_vtk", ""), refinementStep_));
+            }
+        }
+
+        bool adapt()
+        {
+            if (parset_.get("refinement.maxlevel",8) > grid_.maxLevel())
+            {
+                GridFunctionAdaptor<Basis> adaptor(basis_, true);
+                
+                grid_.globalRefine(1);
+                
+                ncBasis_.update(grid_.leafView());
+                basis_.update(grid_.leafView());
+
+                adaptor.adapt(u_);
+                return true;
+            }
+            return false;
+        }
+
+    private:
+
+//        const Dune::ParameterTree& parset_;
+        Dune::ParameterTree& parset_;
+        GridType& grid_;
+        const FunctionMap& functions_;
+        int refinementStep_;
+
+        NonconformingBasis ncBasis_;
+        Basis basis_;
+
+        BoundaryPatchBase< typename Basis::GridView > boundaryPatch_;
+
+        BoundaryPatchBase< typename GridType::LevelGridView > coarseBoundaryPatch_;
+
+        FineMatrixType stiffnessMatrix_;
+        FineMatrixType massMatrix_;
+
+        FineVectorType u_;
+        FineVectorType rhs_;
+        FineVectorType massWeights_;
+
+        typedef typename Dune::shared_ptr< NonlinearityType > NonlinearityPointer;
+        NonlinearityPointer phi_;
+
+        Dune::BitSetVector<1> isBoundary_;
+
+        typedef Norm<FineVectorType>* NormPointer;
+        NormPointer norm_;
+
+        std::vector<Transfer*> transfer_;
+};
+#endif
+
diff --git a/doc/manual/src/laplace.parset b/doc/manual/src/laplace.parset
new file mode 100755
index 0000000000000000000000000000000000000000..47383ffadf0c3fd01946c950792f7407b4670866
--- /dev/null
+++ b/doc/manual/src/laplace.parset
@@ -0,0 +1,116 @@
+# parameter
+#
+[python]
+# module = "laplace"
+inline = "
+
+from math import *
+
+def dist2(x, y):
+    d = 0
+    for i in range(0, len(x)):
+        d += (x[i]-y[i])**2
+    return d
+
+class Disc:
+    def __init__(self, x0, r, cin=1, cout=0):
+        self.x0 = x0
+        self.r = r
+        self.cin = cin
+        self.cout = cout
+
+    def __call__(self, x, i=0):
+        if (dist2(x,self.x0) < self.r**2):
+            return self.cin
+        return self.cout
+
+def SinSin(x, i=0):
+    r = 1
+    for xx in x:
+       r *= sin(2*pi*xx)
+#         r *= sin(xx)
+    return r
+
+def Diri(x, i=0):
+    return abs(x[0])
+#     if x[0]==0:
+#         return 0
+#     return x[0]/abs(x[0])
+
+importAsDuneFunction2D={}
+importAsDuneFunction2D['SinSin'] = SinSin
+# importAsDuneFunction2D['Peak'] = Disc([0.0,0.0], .5, 1, 0)
+importAsDuneFunction2D['Peak'] = Disc([0.0,0.0], .2, 10, 0)
+importAsDuneFunction2D['Diri'] = Diri
+"
+
+[grid]
+# dgffile2d = ../grids/square_centered_symmetric.dgf
+dgffile2d = square_centered.dgf
+# dgffile2d = ../grids/reentrant_corner.dgf
+# dgffile2d = ../grids/circle.dgf
+# dgffile2d = ../grids/reentrant_corner_circle.dgf
+dgffile3d = cube_centered.dgf
+
+# finefile = output/grid.am
+refine = 0
+
+# backupfile = adaptive_laplace.backup
+
+
+[laplace]
+masslumping = false
+stiffnessfactor = 1.0
+massfactor = 0.0
+# solutionfile_am = adaptive_laplace_%05d.am
+solutionfile_vtk = adaptive_laplace_%05d
+dirichletfunction = zero
+# dirichletfunction = abs
+# dirichletfunction = Diri
+
+
+[laplace.rhs]
+function = one
+# function = Disc
+# function = Peak
+# function = SinSin
+# function = zero
+integrationorder = 4
+
+
+[laplace.refinement]
+type = local
+# type = global
+minlevel = 0
+maxlevel = 4
+tol = 5e-2
+extension = refinedp1
+averagemarking.sigma = 0.9
+# averagemarking.sigma = 0.64
+# averagemarking.sigma = 2
+
+
+[laplace.tnnmg]
+maxiter = 500
+tol = 1e-12
+
+# 3 nonlinear per and posts mothing steps
+finepre = 3
+finepost = 3
+
+# V(3,3)-cycle with 6 smoothin steps on coarsest level
+linearfinesmoothing = 1
+pre = 3
+post = 3
+coarse = 6
+
+# accept fraction >=0.5 of exact correction in nonlinear Gauß-Seidel
+gs_acceptance = 0.5
+gs_tol = 1e+10
+
+# accept fraction >=0.9 of exact correction in linesearch
+linesearch_acceptance = 0.9
+linesearch_tol = 1e-10
+
+verbosity = 20
+indent = '   '
diff --git a/doc/manual/src/square_centered.dgf b/doc/manual/src/square_centered.dgf
new file mode 100755
index 0000000000000000000000000000000000000000..7e04d3fef6f5297f08096867308ee711ff69e5df
--- /dev/null
+++ b/doc/manual/src/square_centered.dgf
@@ -0,0 +1,14 @@
+DGF
+Vertex        % the verticies of the grid
+-1   -1       % vertex 0
+ 1   -1       % vertex 1
+-1    1       % vertex 2
+ 1    1       % vertex 3
+#
+SIMPLEX       % a simplex grid
+0 1 2         % triangle 0
+1 3 2         % triangle 1
+#
+BOUNDARYDOMAIN
+default 1     % all other boundary segments have id 1
+#
\ No newline at end of file