diff --git a/AMDiS/Radiosity/3d/GL/glu.h b/AMDiS/Radiosity/3d/GL/glu.h
deleted file mode 100755
index 8aa962d658c5fafc6b21f199c02f89f9dc68709a..0000000000000000000000000000000000000000
--- a/AMDiS/Radiosity/3d/GL/glu.h
+++ /dev/null
@@ -1,533 +0,0 @@
-/* $Id: glu.h,v 1.1 2005/07/27 13:28:46 weichmann Exp $ */
-
-/*
- * Mesa 3-D graphics library
- * Version:  3.3
- *
- * Copyright (C) 1995-1999  Brian Paul
- *
- * This library is free software; you can redistribute it and/or
- * modify it under the terms of the GNU Library General Public
- * License as published by the Free Software Foundation; either
- * version 2 of the License, or (at your option) any later version.
- *
- * This library is distributed in the hope that it will be useful,
- * but WITHOUT ANY WARRANTY; without even the implied warranty of
- * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
- * Library General Public License for more details.
- *
- * You should have received a copy of the GNU Library General Public
- * License along with this library; if not, write to the Free
- * Software Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
- */
-
-
-#ifndef __glu_h__
-#define __glu_h__
-
-
-#if defined(USE_MGL_NAMESPACE)
-#include "glu_mangle.h"
-#endif
-
-#include "GL/gl.h"
-
-
-#ifdef __cplusplus
-extern "C" {
-#endif
-
-
-	/* to facilitate clean DLL building ... */
-#if !defined(OPENSTEP) && (defined(__WIN32__) || defined(__CYGWIN32__))
-#	if defined(_MSC_VER) && defined(BUILD_GLU32) /* tag specify we're building mesa as a DLL */
-#		define GLUAPI __declspec(dllexport)
-#	elif defined(_MSC_VER) && defined(_DLL) /* tag specifying we're building for DLL runtime support */
-#		define GLUAPI __declspec(dllimport)
-#	else /* for use with static link lib build of Win32 edition only */
-#		define GLUAPI extern
-#	endif /* _STATIC_MESA support */
-#else
-#	define GLUAPI extern
-#endif /* WIN32 / CYGWIN32 bracket */
-
-#ifdef macintosh
-	#pragma enumsalwaysint on
-	#if PRAGMA_IMPORT_SUPPORTED
-	#pragma import on
-	#endif
-#endif
-
-#ifndef GLUAPI
-#define GLUAPI
-#endif
-
-#ifndef GLAPIENTRY
-#define GLAPIENTRY
-#endif
-
-#ifndef GLCALLBACK
-#define GLCALLBACK
-#endif
-
-
-#define GLU_VERSION_1_1			1
-#define GLU_VERSION_1_2			1
-
-
-#define GLU_TRUE			1
-#define GLU_FALSE			0
-
-
-/* Normal vectors */
-#define GLU_SMOOTH			100000
-#define GLU_FLAT			100001
-#define GLU_NONE			100002
-
-/* Quadric draw styles */
-#define GLU_POINT			100010
-#define GLU_LINE			100011
-#define GLU_FILL			100012
-#define GLU_SILHOUETTE			100013
-
-/* Quadric orientation */
-#define GLU_OUTSIDE			100020
-#define GLU_INSIDE			100021
-
-/* Tessellator */
-#define GLU_TESS_BEGIN			100100
-#define GLU_TESS_VERTEX			100101
-#define GLU_TESS_END			100102
-#define GLU_TESS_ERROR			100103
-#define GLU_TESS_EDGE_FLAG		100104
-#define GLU_TESS_COMBINE		100105
-
-#define GLU_TESS_BEGIN_DATA		100106
-#define GLU_TESS_VERTEX_DATA		100107
-#define GLU_TESS_END_DATA		100108
-#define GLU_TESS_ERROR_DATA		100109
-#define GLU_TESS_EDGE_FLAG_DATA		100110
-#define GLU_TESS_COMBINE_DATA		100111
-
-/* Winding rules */
-#define GLU_TESS_WINDING_ODD		100130
-#define GLU_TESS_WINDING_NONZERO	100131
-#define GLU_TESS_WINDING_POSITIVE	100132
-#define GLU_TESS_WINDING_NEGATIVE	100133
-#define GLU_TESS_WINDING_ABS_GEQ_TWO	100134
-
-/* Tessellation properties */
-#define GLU_TESS_WINDING_RULE		100140
-#define GLU_TESS_BOUNDARY_ONLY		100141
-#define GLU_TESS_TOLERANCE		100142
-
-/* Tessellation errors */
-#define GLU_TESS_ERROR1			100151  /* Missing gluBeginPolygon */
-#define GLU_TESS_ERROR2			100152  /* Missing gluBeginContour */
-#define GLU_TESS_ERROR3			100153  /* Missing gluEndPolygon */
-#define GLU_TESS_ERROR4			100154  /* Missing gluEndContour */
-#define GLU_TESS_ERROR5			100155  /* */
-#define GLU_TESS_ERROR6			100156  /* */
-#define GLU_TESS_ERROR7			100157  /* */
-#define GLU_TESS_ERROR8			100158  /* */
-
-/* NURBS */
-#define GLU_AUTO_LOAD_MATRIX		100200
-#define GLU_CULLING			100201
-#define GLU_PARAMETRIC_TOLERANCE	100202
-#define GLU_SAMPLING_TOLERANCE		100203
-#define GLU_DISPLAY_MODE		100204
-#define GLU_SAMPLING_METHOD		100205
-#define GLU_U_STEP			100206
-#define GLU_V_STEP			100207
-
-#define GLU_PATH_LENGTH			100215
-#define GLU_PARAMETRIC_ERROR		100216
-#define GLU_DOMAIN_DISTANCE		100217
-
-#define GLU_MAP1_TRIM_2			100210
-#define GLU_MAP1_TRIM_3			100211
-
-#define GLU_OUTLINE_POLYGON		100240
-#define GLU_OUTLINE_PATCH		100241
-
-#define GLU_NURBS_ERROR1	100251   /* spline order un-supported */
-#define GLU_NURBS_ERROR2	100252   /* too few knots */
-#define GLU_NURBS_ERROR3	100253   /* valid knot range is empty */
-#define GLU_NURBS_ERROR4	100254   /* decreasing knot sequence */
-#define GLU_NURBS_ERROR5	100255   /* knot multiplicity > spline order */
-#define GLU_NURBS_ERROR6	100256   /* endcurve() must follow bgncurve() */
-#define GLU_NURBS_ERROR7	100257   /* bgncurve() must precede endcurve() */
-#define GLU_NURBS_ERROR8	100258   /* ctrlarray or knot vector is NULL */
-#define GLU_NURBS_ERROR9 	100259   /* can't draw pwlcurves */
-#define GLU_NURBS_ERROR10	100260   /* missing gluNurbsCurve() */
-#define GLU_NURBS_ERROR11	100261   /* missing gluNurbsSurface() */
-#define GLU_NURBS_ERROR12	100262   /* endtrim() must precede endsurface() */
-#define GLU_NURBS_ERROR13	100263   /* bgnsurface() must precede endsurface() */
-#define GLU_NURBS_ERROR14	100264   /* curve of improper type passed as trim curve */
-#define GLU_NURBS_ERROR15	100265   /* bgnsurface() must precede bgntrim() */
-#define GLU_NURBS_ERROR16	100266   /* endtrim() must follow bgntrim() */
-#define GLU_NURBS_ERROR17	100267   /* bgntrim() must precede endtrim()*/
-#define GLU_NURBS_ERROR18	100268   /* invalid or missing trim curve*/
-#define GLU_NURBS_ERROR19	100269   /* bgntrim() must precede pwlcurve() */
-#define GLU_NURBS_ERROR20	100270   /* pwlcurve referenced twice*/
-#define GLU_NURBS_ERROR21	100271   /* pwlcurve and nurbscurve mixed */
-#define GLU_NURBS_ERROR22	100272   /* improper usage of trim data type */
-#define GLU_NURBS_ERROR23	100273   /* nurbscurve referenced twice */
-#define GLU_NURBS_ERROR24	100274   /* nurbscurve and pwlcurve mixed */
-#define GLU_NURBS_ERROR25	100275   /* nurbssurface referenced twice */
-#define GLU_NURBS_ERROR26	100276   /* invalid property */
-#define GLU_NURBS_ERROR27	100277   /* endsurface() must follow bgnsurface() */
-#define GLU_NURBS_ERROR28	100278   /* intersecting or misoriented trim curves */
-#define GLU_NURBS_ERROR29	100279   /* intersecting trim curves */
-#define GLU_NURBS_ERROR30	100280   /* UNUSED */
-#define GLU_NURBS_ERROR31	100281   /* unconnected trim curves */
-#define GLU_NURBS_ERROR32	100282   /* unknown knot error */
-#define GLU_NURBS_ERROR33	100283   /* negative vertex count encountered */
-#define GLU_NURBS_ERROR34	100284   /* negative byte-stride */
-#define GLU_NURBS_ERROR35	100285   /* unknown type descriptor */
-#define GLU_NURBS_ERROR36	100286   /* null control point reference */
-#define GLU_NURBS_ERROR37	100287   /* duplicate point on pwlcurve */
-
-/* GLU 1.3 and later */
-#define GLU_NURBS_MODE ?
-
-
-/* Errors */
-#define GLU_INVALID_ENUM		100900
-#define GLU_INVALID_VALUE		100901
-#define GLU_OUT_OF_MEMORY		100902
-#define GLU_INCOMPATIBLE_GL_VERSION	100903
-
-/* GLU 1.1 and later */
-#define GLU_VERSION			100800
-#define GLU_EXTENSIONS			100801
-
-
-
-/*** GLU 1.0 tessellation - obsolete! ***/
-
-/* Contour types */
-#define GLU_CW				100120
-#define GLU_CCW				100121
-#define GLU_INTERIOR			100122
-#define GLU_EXTERIOR			100123
-#define GLU_UNKNOWN			100124
-
-/* Tessellator */
-#define GLU_BEGIN			GLU_TESS_BEGIN
-#define GLU_VERTEX			GLU_TESS_VERTEX
-#define GLU_END				GLU_TESS_END
-#define GLU_ERROR			GLU_TESS_ERROR
-#define GLU_EDGE_FLAG			GLU_TESS_EDGE_FLAG
-
-
-/*
- * These are the GLU 1.1 typedefs.  GLU 1.3 has different ones!
- */
-#if defined(__BEOS__)
-    /* The BeOS does something funky and makes these typedefs in one
-     * of its system headers.
-     */
-#else
-    typedef struct GLUquadric GLUquadricObj;
-    typedef struct GLUnurbs GLUnurbsObj;
-
-    /* FIXME: We need to implement the other 1.3 typedefs - GH */
-    typedef struct GLUtesselator GLUtesselator;
-    typedef GLUtesselator GLUtriangulatorObj;
-#endif
-
-
-
-#if defined(__BEOS__) || defined(__QUICKDRAW__)
-#pragma export on
-#endif
-
-
-/*
- *
- * Miscellaneous functions
- *
- */
-
-GLUAPI void GLAPIENTRY gluLookAt( GLdouble eyex, GLdouble eyey, GLdouble eyez,
-                                  GLdouble centerx, GLdouble centery,
-                                  GLdouble centerz,
-                                  GLdouble upx, GLdouble upy, GLdouble upz );
-
-
-GLUAPI void GLAPIENTRY gluOrtho2D( GLdouble left, GLdouble right,
-                                   GLdouble bottom, GLdouble top );
-
-
-GLUAPI void GLAPIENTRY gluPerspective( GLdouble fovy, GLdouble aspect,
-                                       GLdouble zNear, GLdouble zFar );
-
-
-GLUAPI void GLAPIENTRY gluPickMatrix( GLdouble x, GLdouble y,
-                                      GLdouble width, GLdouble height,
-                                      const GLint viewport[4] );
-
-GLUAPI GLint GLAPIENTRY gluProject( GLdouble objx, GLdouble objy, GLdouble objz,
-                                    const GLdouble modelMatrix[16],
-                                    const GLdouble projMatrix[16],
-                                    const GLint viewport[4],
-                                    GLdouble *winx, GLdouble *winy,
-                                    GLdouble *winz );
-
-GLUAPI GLint GLAPIENTRY gluUnProject( GLdouble winx, GLdouble winy,
-                                      GLdouble winz,
-                                      const GLdouble modelMatrix[16],
-                                      const GLdouble projMatrix[16],
-                                      const GLint viewport[4],
-                                      GLdouble *objx, GLdouble *objy,
-                                      GLdouble *objz );
-
-GLUAPI const GLubyte* GLAPIENTRY gluErrorString( GLenum errorCode );
-
-
-
-/*
- *
- * Mipmapping and image scaling
- *
- */
-
-GLUAPI GLint GLAPIENTRY gluScaleImage( GLenum format,
-                                       GLint widthin, GLint heightin,
-                                       GLenum typein, const void *datain,
-                                       GLint widthout, GLint heightout,
-                                       GLenum typeout, void *dataout );
-
-GLUAPI GLint GLAPIENTRY gluBuild1DMipmaps( GLenum target, GLint components,
-                                           GLint width, GLenum format,
-                                           GLenum type, const void *data );
-
-GLUAPI GLint GLAPIENTRY gluBuild2DMipmaps( GLenum target, GLint components,
-                                           GLint width, GLint height,
-                                           GLenum format,
-                                           GLenum type, const void *data );
-
-
-
-/*
- *
- * Quadrics
- *
- */
-
-GLUAPI GLUquadricObj* GLAPIENTRY gluNewQuadric( void );
-
-GLUAPI void GLAPIENTRY gluDeleteQuadric( GLUquadricObj *state );
-
-GLUAPI void GLAPIENTRY gluQuadricDrawStyle( GLUquadricObj *quadObject,
-                                            GLenum drawStyle );
-
-GLUAPI void GLAPIENTRY gluQuadricOrientation( GLUquadricObj *quadObject,
-                                              GLenum orientation );
-
-GLUAPI void GLAPIENTRY gluQuadricNormals( GLUquadricObj *quadObject,
-                                          GLenum normals );
-
-GLUAPI void GLAPIENTRY gluQuadricTexture( GLUquadricObj *quadObject,
-                                          GLboolean textureCoords );
-
-GLUAPI void GLAPIENTRY gluQuadricCallback( GLUquadricObj *qobj,
-                                           GLenum which,
-                                           void (GLCALLBACK *fn)() );
-
-GLUAPI void GLAPIENTRY gluCylinder( GLUquadricObj *qobj,
-                                    GLdouble baseRadius,
-                                    GLdouble topRadius,
-                                    GLdouble height,
-                                    GLint slices, GLint stacks );
-
-GLUAPI void GLAPIENTRY gluSphere( GLUquadricObj *qobj,
-                                  GLdouble radius, GLint slices,
-                                  GLint stacks );
-
-GLUAPI void GLAPIENTRY gluDisk( GLUquadricObj *qobj,
-                                GLdouble innerRadius, GLdouble outerRadius,
-                                GLint slices, GLint loops );
-
-GLUAPI void GLAPIENTRY gluPartialDisk( GLUquadricObj *qobj, GLdouble innerRadius,
-                                       GLdouble outerRadius, GLint slices,
-                                       GLint loops, GLdouble startAngle,
-                                       GLdouble sweepAngle );
-
-
-
-/*
- *
- * Nurbs
- *
- */
-
-GLUAPI GLUnurbsObj* GLAPIENTRY gluNewNurbsRenderer( void );
-
-GLUAPI void GLAPIENTRY gluDeleteNurbsRenderer( GLUnurbsObj *nobj );
-
-GLUAPI void GLAPIENTRY gluLoadSamplingMatrices( GLUnurbsObj *nobj,
-                                                const GLfloat modelMatrix[16],
-                                                const GLfloat projMatrix[16],
-                                                const GLint viewport[4] );
-
-GLUAPI void GLAPIENTRY gluNurbsProperty( GLUnurbsObj *nobj, GLenum property,
-                                         GLfloat value );
-
-GLUAPI void GLAPIENTRY gluGetNurbsProperty( GLUnurbsObj *nobj, GLenum property,
-                                            GLfloat *value );
-
-GLUAPI void GLAPIENTRY gluBeginCurve( GLUnurbsObj *nobj );
-
-GLUAPI void GLAPIENTRY gluEndCurve( GLUnurbsObj * nobj );
-
-GLUAPI void GLAPIENTRY gluNurbsCurve( GLUnurbsObj *nobj, GLint nknots,
-                                      GLfloat *knot, GLint stride,
-                                      GLfloat *ctlarray, GLint order,
-                                      GLenum type );
-
-GLUAPI void GLAPIENTRY gluBeginSurface( GLUnurbsObj *nobj );
-
-GLUAPI void GLAPIENTRY gluEndSurface( GLUnurbsObj * nobj );
-
-GLUAPI void GLAPIENTRY gluNurbsSurface( GLUnurbsObj *nobj,
-                                        GLint sknot_count, GLfloat *sknot,
-                                        GLint tknot_count, GLfloat *tknot,
-                                        GLint s_stride, GLint t_stride,
-                                        GLfloat *ctlarray,
-                                        GLint sorder, GLint torder,
-                                        GLenum type );
-
-GLUAPI void GLAPIENTRY gluBeginTrim( GLUnurbsObj *nobj );
-
-GLUAPI void GLAPIENTRY gluEndTrim( GLUnurbsObj *nobj );
-
-GLUAPI void GLAPIENTRY gluPwlCurve( GLUnurbsObj *nobj, GLint count,
-                                    GLfloat *array, GLint stride,
-                                    GLenum type );
-
-GLUAPI void GLAPIENTRY gluNurbsCallback( GLUnurbsObj *nobj, GLenum which,
-                                         void (GLCALLBACK *fn)() );
-
-
-
-/*
- *
- * Polygon tessellation
- *
- */
-
-GLUAPI GLUtesselator* GLAPIENTRY gluNewTess( void );
-
-GLUAPI void GLAPIENTRY gluDeleteTess( GLUtesselator *tobj );
-
-GLUAPI void GLAPIENTRY gluTessBeginPolygon( GLUtesselator *tobj,
-					    void *polygon_data );
-
-GLUAPI void GLAPIENTRY gluTessBeginContour( GLUtesselator *tobj );
-
-GLUAPI void GLAPIENTRY gluTessVertex( GLUtesselator *tobj, GLdouble coords[3],
-				      void *vertex_data );
-
-GLUAPI void GLAPIENTRY gluTessEndContour( GLUtesselator *tobj );
-
-GLUAPI void GLAPIENTRY gluTessEndPolygon( GLUtesselator *tobj );
-
-GLUAPI void GLAPIENTRY gluTessProperty( GLUtesselator *tobj, GLenum which,
-					GLdouble value );
-
-GLUAPI void GLAPIENTRY gluTessNormal( GLUtesselator *tobj, GLdouble x,
-				      GLdouble y, GLdouble z );
-
-GLUAPI void GLAPIENTRY gluTessCallback( GLUtesselator *tobj, GLenum which,
-					void (GLCALLBACK *fn)() );
-
-GLUAPI void GLAPIENTRY gluGetTessProperty( GLUtesselator *tobj, GLenum which,
-					   GLdouble *value );
-
-/*
- *
- * Obsolete 1.0 tessellation functions
- *
- */
-
-GLUAPI void GLAPIENTRY gluBeginPolygon( GLUtesselator *tobj );
-
-GLUAPI void GLAPIENTRY gluNextContour( GLUtesselator *tobj, GLenum type );
-
-GLUAPI void GLAPIENTRY gluEndPolygon( GLUtesselator *tobj );
-
-
-
-/*
- *
- * New functions in GLU 1.1
- *
- */
-
-GLUAPI const GLubyte* GLAPIENTRY gluGetString( GLenum name );
-
-
-
-/*
- *
- * GLU 1.3 functions
- *
- */
-
-GLUAPI GLboolean GLAPIENTRY
-gluCheckExtension(const char *extName, const GLubyte *extString);
-
-
-GLUAPI GLint GLAPIENTRY
-gluBuild3DMipmaps( GLenum target, GLint internalFormat, GLsizei width,
-                   GLsizei height, GLsizei depth, GLenum format,
-                   GLenum type, const void *data );
-
-GLUAPI GLint GLAPIENTRY
-gluBuild1DMipmapLevels( GLenum target, GLint internalFormat, GLsizei width,
-                        GLenum format, GLenum type, GLint level, GLint base,
-                        GLint max, const void *data );
-
-GLUAPI GLint GLAPIENTRY
-gluBuild2DMipmapLevels( GLenum target, GLint internalFormat, GLsizei width,
-                        GLsizei height, GLenum format, GLenum type,
-                        GLint level, GLint base, GLint max,
-                        const void *data );
-
-GLUAPI GLint GLAPIENTRY
-gluBuild3DMipmapLevels( GLenum target, GLint internalFormat, GLsizei width,
-                        GLsizei height, GLsizei depth, GLenum format,
-                        GLenum type, GLint level, GLint base, GLint max,
-                        const void *data );
-
-GLUAPI GLint GLAPIENTRY
-gluUnProject4( GLdouble winx, GLdouble winy, GLdouble winz, GLdouble clipw,
-               const GLdouble modelMatrix[16], const GLdouble projMatrix[16],
-               const GLint viewport[4], GLclampd zNear, GLclampd zFar,
-               GLdouble *objx, GLdouble *objy, GLdouble *objz,
-               GLdouble *objw );
-
-
-
-#if defined(__BEOS__) || defined(__QUICKDRAW__)
-#pragma export off
-#endif
-
-
-#ifdef macintosh
-	#pragma enumsalwaysint reset
-	#if PRAGMA_IMPORT_SUPPORTED
-	#pragma import off
-	#endif
-#endif
-
-
-#ifdef __cplusplus
-}
-#endif
-
-
-#endif /* __glu_h__ */
diff --git a/AMDiS/levelset/Makefile b/AMDiS/levelset/Makefile
deleted file mode 100644
index f95338fef65849e81459289f24390184a81433bc..0000000000000000000000000000000000000000
--- a/AMDiS/levelset/Makefile
+++ /dev/null
@@ -1,46 +0,0 @@
-
-SHELL = /bin/sh
-CXX = g++
-DEFS = -DPACKAGE=\"AMDiS\" -DVERSION=\"0.1\" -DHAVE_DLFCN_H=1 
-DEFAULT_INCLUDES =  -I.
-INCLUDES = -I/solhome/vey/sourcen/AMDiS//src
-CPPFLAGS = 
-top_builddir = ../
-LIBS = 
-LIBTOOL = $(SHELL) $(top_builddir)/libtool
-
-LDADD = /solhome/vey/sourcen/AMDiS//lib/libAMDiS_debug.la
-#LDADD = /solhome/vey/sourcen/AMDiS//lib/libAMDiS.la
-PCXXFLAGS = -g -O0
-#PCXXFLAGS = -O2
-PLDFLAGS = -static
-#PLDFLAGS = 
-
-CXXFLAGS = -ftemplate-depth-30 $(PCXXFLAGS)
-LDFLAGS = $(PLDFLAGS)
-
-all : 
-	make $(PROGRAMS)
-
-clean: 
-	-rm -rf *.o
-	-rm -rf $(PROGRAMS)
-
-.cc.o: src/$*.cc
-	$(CXX) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(AM_CPPFLAGS) $(CPPFLAGS) $(AM_CXXFLAGS) $(CXXFLAGS) -c -o $*.o src/$*.cc
-
-
-# ============================================================================
-# ===== user programs ========================================================
-# ============================================================================
-
-PROGRAMS = levelset
-
-VPATH = ./src/:./
-
-# ===== ellipt ==============================================================
-
-LEVELSET_OFILES = levelset.o
-
-levelset: $(LEVELSET_OFILES)
-	$(LIBTOOL) --mode=link $(CXX) $(CXXFLAGS) $(PCXXFLAGS) -o levelset levelset.o $(LDADD)
diff --git a/AMDiS/levelset/init/levelset.dat.2d b/AMDiS/levelset/init/levelset.dat.2d
deleted file mode 100644
index b6c8d1c4b9b94cbed5e3a2bc8424860a3d382b26..0000000000000000000000000000000000000000
--- a/AMDiS/levelset/init/levelset.dat.2d
+++ /dev/null
@@ -1,63 +0,0 @@
-dimension of world:   2
-
-levelsetMesh->macro file name:            ./macro/macro.stand.2d
-levelsetMesh->global refinements:         0
-
-levelset->space->polynomial degree:       2
-
-
-%-------------------------- problem parameters -----------------------------
-levelset->epsilon:                      1.e-8
-levelset->alpha:                        1.0
-levelset->problem number:               0
-
-
-
-
-levelset->space->mesh:                    levelsetMesh
-
-levelset->space->solver:                  2 % 1: BICGSTAB 2: CG 3: GMRES 4: ODIR 5: ORES
-levelset->space->solver->max iteration:   1000
-levelset->space->solver->restart:         10  %  only used for GMRES
-levelset->space->solver->tolerance:       1.e-8
-levelset->space->solver->info:            8
-levelset->space->solver->precon:          1   % 0: no precon 1: diag precon
-
-levelset->estimator C0:                  1.0
-levelset->estimator C1:                  1.0
-levelset->estimator C2:                  1.0
-levelset->estimator C3:                  1.0
-
-levelset->theta:                         1.0
-
-levelset->adapt->tolerance:              0.01
-
-levelset->adapt->timestep:               0.01
-levelset->adapt->startTime:              0.0
-levelset->adapt->endTime:                0.05
-levelset->adapt->rel initial error:      0.5
-levelset->adapt->rel space error:        0.5
-levelset->adapt->rel time error:         0.5
-levelset->adapt->strategy:               0   % 0=explicit, 1=implicit
-levelset->adapt->max iteration:          1
-levelset->adapt->info:                   8
-
-levelset->initial->adapt->strategy:      0    % 0=none, 1=GR, 2=MS, 3=ES, 4=GERS
-levelset->initial->adapt->max iteration: 10
-levelset->initial->adapt->info:          8
-
-levelset->space->dim:                    2
-levelset->space->adapt->strategy:        3    % 0=none, 1=GR, 2=MS, 3=ES, 4=GERS
-levelset->space->adapt->ESTheta:        0.9
-levelset->space->adapt->ESThetaC:      0.05
-levelset->space->adapt->max iteration:   2
-levelset->space->adapt->coarsen allowed: 1   % 0|1
-levelset->space->adapt->info:            8
-
-WAIT: 0
-
-
-
-
-
-
diff --git a/AMDiS/levelset/macro/macro.stand.1d b/AMDiS/levelset/macro/macro.stand.1d
deleted file mode 100644
index 31eb45ca4222eec6566dc69fd74a8d2791798424..0000000000000000000000000000000000000000
--- a/AMDiS/levelset/macro/macro.stand.1d
+++ /dev/null
@@ -1,15 +0,0 @@
-DIM: 1
-DIM_OF_WORLD: 1
-
-number of elements: 1  
-number of vertices: 2  
-
-element vertices:
-0 1 
-
-element boundaries:
-1 1 
-
-vertex coordinates:
- 0.0  
- 1.0
diff --git a/AMDiS/levelset/macro/macro.stand.2d b/AMDiS/levelset/macro/macro.stand.2d
deleted file mode 100644
index f81cc45794cc9986247ae5953d71dbd6e6c86d33..0000000000000000000000000000000000000000
--- a/AMDiS/levelset/macro/macro.stand.2d
+++ /dev/null
@@ -1,30 +0,0 @@
-DIM: 2
-DIM_OF_WORLD: 2
-
-number of elements: 4  
-number of vertices: 5  
-
-element vertices:
-0 1 4 
-1 2 4  
-2 3 4 
-3 0 4 
-
-element boundaries:
-0 0 1 
-0 0 1 
-0 0 1 
-0 0 1 
-
-vertex coordinates:
- 0.0 0.0  
- 1.0 0.0
- 1.0 1.0
- 0.0 1.0
- 0.5 0.5
-
-element neighbours:
-1 3 -1 
-2 0 -1 
-3 1 -1 
-0 2 -1 
diff --git a/AMDiS/levelset/macro/macro.stand.3d b/AMDiS/levelset/macro/macro.stand.3d
deleted file mode 100644
index 2bde72fd05fdaf113903d2f185a727eb9ed2a7da..0000000000000000000000000000000000000000
--- a/AMDiS/levelset/macro/macro.stand.3d
+++ /dev/null
@@ -1,39 +0,0 @@
-DIM:          3 
-DIM_OF_WORLD: 3
-
-number of vertices: 8
-number of elements: 6
-
-vertex coordinates:
-  0.0  0.0  0.0
-  1.0  0.0  0.0
-  0.0  0.0  1.0
-  1.0  0.0  1.0
-  1.0  1.0  0.0
-  1.0  1.0  1.0
-  0.0  1.0  0.0
-  0.0  1.0  1.0
-
-element vertices:
-  0    5    4    1
-  0    5    3    1
-  0    5    3    2
-  0    5    4    6
-  0    5    7    6
-  0    5    7    2
-
-element boundaries:
-  1    1    0    0
-  1    1    0    0 
-  1    1    0    0
-  1    1    0    0
-  1    1    0    0
-  1    1    0    0
-
-element neighbours:
- -1   -1    1    3
- -1   -1    0    2
- -1   -1    5    1
- -1   -1    4    0
- -1   -1    3    5
- -1   -1    2    4
diff --git a/AMDiS/levelset/macro/macro_big.stand.1d b/AMDiS/levelset/macro/macro_big.stand.1d
deleted file mode 100644
index a1bf307b5696d1a6c9fd91f4f764d44778595745..0000000000000000000000000000000000000000
--- a/AMDiS/levelset/macro/macro_big.stand.1d
+++ /dev/null
@@ -1,15 +0,0 @@
-DIM: 1
-DIM_OF_WORLD: 1
-
-number of vertices: 2
-number of elements: 1
-
-vertex coordinates:
--1.0
- 1.0  
-
-element vertices:
-0 1
-
-element boundaries:
-1 1
\ No newline at end of file
diff --git a/AMDiS/levelset/macro/macro_big.stand.2d b/AMDiS/levelset/macro/macro_big.stand.2d
deleted file mode 100644
index 5aee8edf9ec2ee0c631e776d4b926eee4a86e7a1..0000000000000000000000000000000000000000
--- a/AMDiS/levelset/macro/macro_big.stand.2d
+++ /dev/null
@@ -1,24 +0,0 @@
-DIM: 2
-DIM_OF_WORLD: 2
-
-number of vertices: 5
-number of elements: 4
-
-vertex coordinates:
--1.0 -1.0  
- 1.0 -1.0
- 1.0  1.0
--1.0  1.0
- 0.0  0.0
-
-element vertices:
-0 1 4
-1 2 4 
-2 3 4
-3 0 4
-
-element boundaries:
-0 0 1
-0 0 1
-0 0 1
-0 0 1
diff --git a/AMDiS/levelset/macro/macro_big.stand.3d b/AMDiS/levelset/macro/macro_big.stand.3d
deleted file mode 100644
index 4deeef781ba49d903a3da21dd9f320a38793bd2d..0000000000000000000000000000000000000000
--- a/AMDiS/levelset/macro/macro_big.stand.3d
+++ /dev/null
@@ -1,40 +0,0 @@
-DIM:          3 
-DIM_OF_WORLD: 3
-
-number of vertices: 8
-number of elements: 6
-
-vertex coordinates:
- -1.0 -1.0 -1.0
-  1.0 -1.0 -1.0
- -1.0 -1.0  1.0
-  1.0 -1.0  1.0
-  1.0  1.0 -1.0
-  1.0  1.0  1.0
- -1.0  1.0 -1.0
- -1.0  1.0  1.0
-
-element vertices:
-  0    5    4    1
-  0    5    3    1
-  0    5    3    2
-  0    5    4    6
-  0    5    7    6
-  0    5    7    2
-
-element boundaries:
-  1    1    0    0
-  1    1    0    0 
-  1    1    0    0
-  1    1    0    0
-  1    1    0    0
-  1    1    0    0
-
-element neighbours:
- -1   -1    1    3
- -1   -1    0    2
- -1   -1    5    1
- -1   -1    4    0
- -1   -1    3    5
- -1   -1    2    4
-
diff --git a/AMDiS/levelset/src/levelset.cc b/AMDiS/levelset/src/levelset.cc
deleted file mode 100644
index 035255e5a14f20fab409f0ab724490e42efc681b..0000000000000000000000000000000000000000
--- a/AMDiS/levelset/src/levelset.cc
+++ /dev/null
@@ -1,669 +0,0 @@
-#include "AMDiS.h"
-
-using namespace std;
-
-
-
-
-static double epsNorm1(const WorldVector x, const double eps)
-{
-  static double result = 0.0;
-  static int k;
-  
-  result = 0.0;
-  for(k=0; k<x.getSize(); k++){
-    result +=  x[k]*x[k];
-  }
-  result += eps*eps;
-  result = sqrt( result );
-  TEST_EXIT( result > 1.e-10)( "eps_norm = %e \n", result);
-  result = 1.0/result;  
-
-  return(result);
-}
-
-
-static double epsNorm(const WorldVector x, const double eps)
-{
-  static double result = 0.0;
-  static int  k;
-  
-  result = 0.0;
-
-  for(k=0; k<x.getSize(); k++){
-    result +=  x[k]*x[k];
-  }
-  result += eps*eps;
-  result = sqrt( result );
-
-  return(result);
-}
-
-static void epsNormalize(const WorldVector x, const double eps, WorldVector y)
-{
-  static int k;
-  
-  for(k=0; k<x.getSize(); k++)
-    y[k] = x[k] * epsNorm1(x,eps);
-  
-  return;
-}
-
-
-
-class Circle : public AbstractFunction<double, double>
-{
-public:
-  Circle(double radius_) : radius(radius_) {};
-
-  const double& f(const double& alpha) const
-  {
-    return radius;
-  };
-
-protected:
-  double radius;
-};
-
-
-class PerturbedCircle : public AbstractFunction<double, double>
-{
-public:
-  PerturbedCircle(double radius_, double strength_, int period_) : 
-    radius(radius_), strength(strength_), period(period_) {};
-
-  const double& f(const double &alpha) const
-  {
-    static double result; 
-    result = 1.0 + strength * cos( period * alpha );
-    result *= radius;
-    return result;
-  };
-
-protected:
-  double radius;
-  double strength;
-  int period;
-};
-
-
-/****************************************************************************/
-/*       signed distance function for levelset 0 circle with radius         */
-/****************************************************************************/
-
-
-static double signedDistStar(const WorldVector& x, const WorldVector& midPoint, 
-			     AbstractFunction<double, double> *radius, double eps)
-{
-  WorldVector  x_trans;
-  double norm_x;
-  double result;
-  double alpha;
-  
-  int k;
-  for(k=0; k < x.getSize(); k++) {
-    x_trans[k] = x[k] - midPoint[k];
-  }
-
-  //x_trans = (WorldVector&) x - (WorldVector&) mid_point;
-  
-  norm_x = epsNorm(x_trans, eps);
-  TEST_EXIT( norm_x > 1.e-15)( "eps_norm = %e \n", norm_x);
-  
-  if( x_trans[0] >= 0.0 )
-    alpha = acos( x_trans[0] / norm_x );
-  else 
-    alpha = 2.0*M_PI - acos( x_trans[0] / norm_x );
-  
-  result = ( norm_x - (*radius)(alpha)); 
-    
-  return(result);
-};
-
-
-
-/****************************************************************************/
-/* anisotropy function and derivative                                        */
-/****************************************************************************/
-
-class Gamma : public AbstractFunction<double, WorldVector>
-{
-public:
-  Gamma(double alpha_) : alpha(alpha_) {};
-
-  const double& f(const WorldVector& x) const
-  {
-    static double result;
-    int  k;
-
-    result = alpha * x[0] * x[0] + 1.e-10;
-    for( k=1; k<x.getSize(); k++)
-      result += x[k]*x[k];
-
-    result = sqrt(result);
-    return result;
-  };
-
-protected:
-  double alpha;
-};
-
-
-
-
-class Gamma_z : public AbstractFunction<WorldVector, WorldVector>
-{
-public:
-  Gamma_z(double alpha_) : alpha(alpha_) {};
-
-  const WorldVector& f(const WorldVector& x) const
-  {
-    static WorldVector result;
-    int  k;
-    double norm_1 = 0.0;
-    static Gamma gamma(alpha);
-    double norm = gamma(x);
-
-    TEST_EXIT( norm  > 1.e-15)( "gamma(z) = %e \n", norm);
-    norm_1 = 1.0/norm;
-
-    result[0] =  alpha * x[0] * norm_1;
-
-    for( k=1; k<x.getSize(); k++)
-      result[k] = x[k] * norm_1;
-
-    return result;
-  };
-
-protected:
-  double alpha;
-};
-
-
-/****************************************************************************/
-/*     Abstract functions used by operators                                 */
-/****************************************************************************/
-
-/*      function for second order operator A     */
-
-class ProjectionGradNorm : public AbstractFunction<WorldMatrix, WorldVector>
-{
-public:
-  ProjectionGradNorm( double epsilon_) : epsilon(epsilon_) {};
-
-  const WorldMatrix& f(const WorldVector& x) const
-  {
-    static WorldMatrix result;
-    static WorldVector gradNorm;
-    static double      norm;
-    int k,l;
-
-    norm  = epsNorm(x, epsilon);
-    epsNormalize(x, epsilon,  gradNorm);
- 
-    for(k = 0; k<x.getSize(); k++)
-      for(l = 0; l<= k; l++){
-	if( k == l )  
-	  result[k][l] = 1.0 - gradNorm[k]*gradNorm[l];
-	else
-	  result[k][l] =  -gradNorm[k]*gradNorm[l];
-	result[k][l] *= norm;
-	result[l][k] = result[k][l];
-      }
-
-    return result;
-  };
-
-protected:
-  double alpha;
-  double epsilon;
-};
-
-
-/*      function for second order operator B     */
-
-class GammaNorm : public AbstractFunction<double, WorldVector>
-{
-public:
-  GammaNorm( double epsilon_, double alpha_) : epsilon(epsilon_), alpha(alpha_) {};
-
-  const double& f(const WorldVector& x) const
-  {
-    static double result;
-    static WorldVector gradNorm;
-    static double      norm1;
-    static Gamma gamma(alpha);
-    int k,l;
-
-    norm1  = epsNorm1(x, epsilon);
-    epsNormalize(x, epsilon,  gradNorm);
-
-    result = gamma(gradNorm) * norm1;
-    return result;
-  };
-
-protected:
-  double epsilon;
-  double alpha;
-};
-
-
-
-
-
-/*        function for rhs F    */
-class GammaGradNorm : public AbstractFunction<WorldVector, WorldVector>
-{
-public:
-  GammaGradNorm(double alpha_, double epsilon_) : alpha(alpha_), epsilon(epsilon_) {};
-
-  const WorldVector& f(const WorldVector& x) const
-  {
-    static WorldVector result;
-    static WorldVector gradNorm;
-    static Gamma gamma(alpha);
-
-    epsNormalize(x, epsilon,  gradNorm);
-    result = gamma( gradNorm ) *  gradNorm;
-
-    return result;
-  };
-
-protected:
-  double alpha;
-  double epsilon;
-};
-
-
-/*        function for rhs G   */
-class Gamma_zGradNorm : public AbstractFunction<WorldVector, WorldVector>
-{
-public:
-  Gamma_zGradNorm(double alpha_, double epsilon_) : alpha(alpha_), epsilon(epsilon_) {};
-
-  const WorldVector& f(const WorldVector& x) const
-  {
-    static WorldVector result;
-    static WorldVector gradNorm;
-    static Gamma_z gamma_z(alpha);
-
-    epsNormalize(x, epsilon,  gradNorm);
-    result = gamma_z( gradNorm );
-
-    return result;
-  };
-
-protected:
-  double alpha;
-  double epsilon;
-};
-
-
-/**********************************************/
-/* create stationary problem                  */
-/**********************************************/
-
-class LevelsetSpace : public Problem<double>
-{
-public:
-  LevelsetSpace() 
-    : Problem<double>("levelset->space")
-  {
-    GET_PARAMETER(1, name + "->lambda", "%f", &lambda);
-
-    A   = NEW DOFMatrix(feSpace, "A");
-    B   = NEW DOFMatrix(feSpace, "B");
-    M   = NEW DOFMatrix(feSpace, "M");
-    M_lump   = NEW DOFMatrix(feSpace, "M_lump");
-    M_1 = NEW DOFMatrix(feSpace, "M_1");
-
-    f   = NEW DOFVector<double>(feSpace, "f");
-    g   = NEW DOFVector<double>(feSpace, "g");
-    rhs = NEW DOFVector<double>(feSpace, "rhs");
-  };
-
-  ~LevelsetSpace() 
-  {
-    DELETE A;
-    DELETE B;
-    DELETE M;
-    DELETE M_lump;
-    DELETE M_1;
-    DELETE f;
-    DELETE g;
-  };
-
-  DOFMatrix *getMatrixA() { return A; };
-  DOFMatrix *getMatrixB() { return B; };
-  DOFMatrix *getMatrixM() { return M; };
-  DOFVector<double> *getVectorF() { return f; };
-  DOFVector<double> *getVectorG() { return g; };
-
-  static int buildAfterCoarsenFunction(ElInfo *elInfo) 
-  {
-    const char *bound;
-
-    if(traversePtr->useGetBound)
-      bound = traversePtr->feSpace->getBasisFcts()->getBound(elInfo, NULL);
-    else
-      bound = NULL;
-  
-    traversePtr->A->assemble(1.0, elInfo, bound);
-    traversePtr->B->assemble(1.0, elInfo, bound);
-    traversePtr->M->assemble(1.0, elInfo, bound);
-    traversePtr->f->assemble(1.0, elInfo, bound);
-    traversePtr->g->assemble(1.0, elInfo, bound);
-
-    return 0;
-  };
-
-  void buildAfterCoarsen(Flag flag)
-  {
-    FUNCNAME("LevelsetSpace::buildAfterCoarsen");
-
-    MSG("%d DOFs for %s\n", feSpace->getAdmin()->getUsedSize(), feSpace->getName().c_str());
-
-    mesh->dofCompress();
-
-    Flag assembleFlag = 
-      flag | A->getAssembleFlag() | B->getAssembleFlag() | M->getAssembleFlag() | 
-      f->getAssembleFlag() | g->getAssembleFlag();
-
-    if(useGetBound)
-      assembleFlag |= Mesh::FILL_BOUND;
-
-    A->clear();
-    B->clear();
-    M->clear();
-    f->set(0.0);
-    g->set(0.0);
-   
-  
-    traversePtr = this;
-
-    mesh->traverse(-1, assembleFlag, &buildAfterCoarsenFunction);
-    //dirichletBound(gFct, fh, uh, NULL);           /*  boundary values         */
-  };
-  
-  void solve() {};
-
-  double estimate() { return 0.0; };
-
-private:
-
-#if 0
-
-void matrixLump(DOF_MATRIX *matrix, DOF_MATRIX *matrix_diag)
-{
-  int k,i,icol;
-  MATRIX_ROW *row;
-
-  FUNCNAME("matrix_lump");
-  clear_dof_matrix(matrix_diag);
-  FOR_ALL_DOFS(gd->fe_space->admin,    
-	       matrix_diag->matrix_row[dof] = get_matrix_row(gd->fe_space);
-	       matrix_diag->matrix_row[dof]->entry[0]=0.0;
-	       row = matrix->matrix_row[dof];
-	       while( row ) {
-		 for(i=0; i<ROW_LENGTH; i++){
-		   icol = row->col[i];
-		   if( ENTRY_USED(icol) )
-		     matrix_diag->matrix_row[dof]->entry[0] +=
-		       matrix->matrix_row[dof]->entry[i];
-		   else{
-		     if (icol == NO_MORE_ENTRIES)
-		       break;
-		   }
-		 }
-		 row = row->next;
-	       }
-	       matrix_diag->matrix_row[dof]->col[0] = dof;
-	       matrix_diag->matrix_row[dof]->col[1] = NO_MORE_ENTRIES;
-	       for( k = 2; k < ROW_LENGTH; k++)
-	       matrix_diag->matrix_row[dof]->col[k] = UNUSED_ENTRY;
-
-	       matrix_diag->matrix_row[dof]->next = nil;
-	       )
- return;
-}
-
-#endif
-
-/****************************************************************************/
-/*invertDiagMatrix(): inverts a positive definit diagonal matrix          */
-/****************************************************************************/
-
-void invertDiagMatrix(DOFMatrix *matrix, DOFMatrix *matrix_1)
-{
-  int k;
-  FUNCNAME("invert_matrix_lump");
-  matrix_1->clear();
-
-  DOFMatrix::Iterator matrixIt(matrix, USED_DOFS);
-  double entry = 0.0;
-  DegreeOfFreedom dof = 0;
-
-  for(matrixIt.reset(); !matrixIt.end(); ++matrixIt) {
-    dof = matrixIt.getPosition();
-    if( (*matrix)[dof][0].entry > 1.e-10 ){
-      entry = 1.0 / (*matrix)[dof][0].entry;
-      matrix_1->addSparseDofEntry(1.0, dof, dof, entry);
-    }
-    else 
-      ERROR_EXIT("matrix[%d][%d] = %f", dof, dof, (*matrix)[dof][0].entry );
-
-  }
-   
-
- return;
-}
-
-
-
-
-
-private:
-  DOFMatrix *A, *B, *M, *M_lump, *M_1;
-  DOFVector<double> *rhs, *f, *g;
-  
-  double *timestepPtr;
-
-  double lambda;
-
-  static LevelsetSpace *traversePtr;
-
-  friend class Levelset;
-};
-
-LevelsetSpace *LevelsetSpace::traversePtr = NULL;
-
-/**********************************************/
-/* create instationary problem                */
-/**********************************************/
-
-class Levelset : public ProblemInstat<double>
-{
-public:
-  Levelset(LevelsetSpace *levelsetSpace)
-    : ProblemInstat<double>("levelset", (Problem<double>*) levelsetSpace)
-  {
-    
-    timestep = adaptInstat->getTimestepPtr();
-    levelsetSpace->timestepPtr = timestep;
-
-    int problemNumber = -1;
-    GET_PARAMETER(1, name + "->problem number","%d", &problemNumber);
-    epsilon = 0.1;
-    GET_PARAMETER(1, name + "->epsilon","%f", &epsilon);
-    alpha = 1.0;
-    GET_PARAMETER(1, name + "->alpha","%f", &alpha);
-
-
-    switch(problemNumber) {
-    case 0:
-      initialFunction = new InitialFunction0(epsilon);
-      break;
-
-    case 1:
-      initialFunction = new InitialFunction1(epsilon);
-      break;
-
-    case -1:
-      ERROR_EXIT("no problem number specified\n");
-      break;
-    default: ERROR_EXIT("non existing problem number %d\n", problemNumber);
-    }
-
-  };
-
-
-  // === timestep ===================================================
-  void initTimestep() 
-  { 
-    uhOld->copy(*(problemStat->getUh())); 
-  };
-
-  void closeTimestep() 
-  {
-    FUNCNAME("closeTimestep");  
-  };
-
-  // === initial functions ==================================
-  void solveInitial() 
-  {
-    //problemStat->getMesh()->dofCompress();
-    //problemStat->getBoundaryFunction()->setTime(0.0);
-    //problemStat->getUh()->interpol(initialFunction);
-  };
-
-  void setEvaluationTimes() {};
- 
-  double *getTimestepPtr() { return timestep; };
-
-
-  double getAlpha(){ return alpha; };
-  double getEpsilon(){ return epsilon; };
-
-
-private:
-  AbstractFunction<double, WorldVector> *initialFunction;
-
-  double *timestep;
-  double epsilon;
-  double alpha;
-
-private:
-  class InitialFunction0 : public AbstractFunction<double, WorldVector> 
-  {
-  public:
-    InitialFunction0(double epsilon_) : epsilon(epsilon_) {}
-
-    const double& f(const WorldVector& x) const
-    {
-      static WorldVector midPoint(DEFAULT_VALUE, 2.0);
-      static double result;
-      static Circle circle1(1.0);
-      result = signedDistStar( x, midPoint, &circle1, epsilon ); 
-      
-    };  
-
-  private:
-    double epsilon;
-
-  };
-
-  class InitialFunction1 : public AbstractFunction<double, WorldVector> 
-  {
-  public:
-    InitialFunction1(double epsilon_) : epsilon(epsilon_) {}
-
-    const double& f(const WorldVector& x) const
-    {
-      static WorldVector midPoint(DEFAULT_VALUE, 2.0);
-      static double result;
-      static PerturbedCircle perturbedCircle1(1.0, 0.3, 6);
-      result = signedDistStar( x, midPoint, &perturbedCircle1, epsilon ); 
-    };  
-
-  private:
-    double epsilon;
-  };
- 
-};
-
-
-
- 
-int main(int argc, char** argv)
-{
-  TEST_EXIT(argc == 2)("usage: levelset initfile\n");
-
-  Parameters::init(false, argv[1]);
-
-  LevelsetSpace *levelsetSpace = new LevelsetSpace;
-  
-  Levelset *levelset = new Levelset(levelsetSpace);
-  double epsilon = levelset->getEpsilon();
-  double alpha   = levelset->getAlpha();
-
-
-  // create  mass matrix M operator
-  Operator *operatorM = NEW Operator(Operator::MATRIX_OPERATOR | 
-				     Operator::VECTOR_OPERATOR,
-				     levelsetSpace->getFESpace());
-
-  operatorM->addZeroOrderTerm(NEW SimpleZeroOrderTerm);
-
-  operatorM->setUhOld(levelset->getUhOld());
-
-  levelsetSpace->addMatrixOperator(operatorM);
-  levelsetSpace->addVectorOperator(operatorM);
-
-  // create  mass matrix B operator
-  Operator *operatorB = NEW Operator(Operator::MATRIX_OPERATOR,
-				     levelsetSpace->getFESpace());
-
-  operatorB->addSecondOrderTerm(NEW FctGradient_SOT(levelset->getUhOld(),
-						    new GammaNorm(alpha, epsilon),
-						    100));
-
-  // create  mass matrix A operator
-  Operator *operatorA = NEW Operator(Operator::MATRIX_OPERATOR,
-				     levelsetSpace->getFESpace());
-
-  operatorA->addSecondOrderTerm(NEW MatrixGradient_SOT(levelset->getUhOld(),
-						    new ProjectionGradNorm(epsilon),
-						    100));
-
-  // create rhs F 
-  Operator *operatorF = NEW Operator(Operator::VECTOR_OPERATOR,
-				     levelsetSpace->getFESpace());
-
-  operatorF->addFirstOrderTerm(NEW  VectorGradient_FOT(levelset->getUhOld(),
-						       new Gamma_zGradNorm(alpha,epsilon),
-						       100, false));
-
-
-  // create rhs G 
-  Operator *operatorG = NEW Operator(Operator::VECTOR_OPERATOR,
-				     levelsetSpace->getFESpace());
-
-  operatorG->addFirstOrderTerm(NEW  VectorGradient_FOT(levelset->getUhOld(),
-						       new GammaGradNorm(alpha,epsilon),
-						       100, false));
-
-
-  levelsetSpace->getMatrixA()->addOperator(operatorA);
-  levelsetSpace->getMatrixB()->addOperator(operatorB);
-  levelsetSpace->getMatrixM()->addOperator(operatorM);
-  levelsetSpace->getVectorF()->addOperator(operatorF);
-  levelsetSpace->getVectorG()->addOperator(operatorG);
-
-
-  levelsetSpace->addVectorOperator(operatorM);
-
-
-  levelset->adaptMethodInstat();
-
-  return 0;
-}
diff --git a/AMDiS/libtool b/AMDiS/libtool
index 6a7f3c767c5190ca326df932f5aaf9521edace61..95d2fe8903c9787ee9c06e44082aec7d871eaa1b 100755
--- a/AMDiS/libtool
+++ b/AMDiS/libtool
@@ -30,10 +30,10 @@
 # the same distribution terms that you use for the rest of that program.
 
 # A sed program that does not truncate output.
-SED="/bin/sed"
+SED="/usr/bin/sed"
 
 # Sed that helps us avoid accidentally triggering echo(1) options like -n.
-Xsed="/bin/sed -e 1s/^X//"
+Xsed="/usr/bin/sed -e 1s/^X//"
 
 # The HP-UX ksh and POSIX shell print the target directory to stdout
 # if CDPATH is set.
@@ -44,7 +44,7 @@ available_tags=" CXX F77"
 
 # ### BEGIN LIBTOOL CONFIG
 
-# Libtool was configured on host NWRW15:
+# Libtool was configured on host deimos101:
 
 # Shell to use when invoking shell scripts.
 SHELL="/bin/sh"
@@ -66,12 +66,12 @@ fast_install=yes
 
 # The host system.
 host_alias=
-host=i686-pc-linux-gnu
+host=x86_64-unknown-linux-gnu
 host_os=linux-gnu
 
 # The build system.
 build_alias=
-build=i686-pc-linux-gnu
+build=x86_64-unknown-linux-gnu
 build_os=linux-gnu
 
 # An echo program that does not interpret backslashes.
@@ -82,13 +82,13 @@ AR="ar"
 AR_FLAGS="cru"
 
 # A C compiler.
-LTCC="gcc"
+LTCC="/licsoft/libraries/openmpi/1.2.6/64bit/bin/mpicc"
 
 # LTCC compiler flags.
 LTCFLAGS="-g -O2"
 
 # A language-specific compiler.
-CC="gcc"
+CC="/licsoft/libraries/openmpi/1.2.6/64bit/bin/mpicc"
 
 # Is the compiler the GNU C compiler?
 with_gcc=yes
@@ -97,7 +97,7 @@ with_gcc=yes
 EGREP="grep -E"
 
 # The linker used to build libraries.
-LD="/usr/bin/ld"
+LD="/usr/x86_64-suse-linux/bin/ld -m elf_x86_64"
 
 # Whether we need hard or soft links.
 LN_S="ln -s"
@@ -171,7 +171,7 @@ dlopen_self=unknown
 dlopen_self_static=unknown
 
 # Compiler flag to prevent dynamic linking.
-link_static_flag="-static"
+link_static_flag=""
 
 # Compiler flag to turn off builtin functions.
 no_builtin_flag=" -fno-builtin"
@@ -325,10 +325,10 @@ variables_saved_for_relink="PATH LD_LIBRARY_PATH LD_RUN_PATH GCC_EXEC_PREFIX COM
 link_all_deplibs=unknown
 
 # Compile-time system search path for libraries
-sys_lib_search_path_spec=" /u/witkowski/local/lib/i386-redhat-linux/4.1.2/ /u/witkowski/local/lib/ /u/witkowski/local/intel/mkl/10.0.1.014/lib/32/i386-redhat-linux/4.1.2/ /u/witkowski/local/intel/mkl/10.0.1.014/lib/32/ /usr/lib/gcc/i386-redhat-linux/4.1.2/ /usr/lib/gcc/i386-redhat-linux/4.1.2/ /usr/lib/gcc/i386-redhat-linux/4.1.2/../../../../i386-redhat-linux/lib/i386-redhat-linux/4.1.2/ /usr/lib/gcc/i386-redhat-linux/4.1.2/../../../../i386-redhat-linux/lib/ /usr/lib/gcc/i386-redhat-linux/4.1.2/../../../i386-redhat-linux/4.1.2/ /usr/lib/gcc/i386-redhat-linux/4.1.2/../../../ /lib/i386-redhat-linux/4.1.2/ /lib/ /usr/lib/i386-redhat-linux/4.1.2/ /usr/lib/"
+sys_lib_search_path_spec=" /fastfs/wir/local/lib/x86_64-suse-linux/4.1.2/ /fastfs/wir/local/lib/../lib64/ /usr/lib64/gcc/x86_64-suse-linux/4.1.2/ /usr/lib/gcc/x86_64-suse-linux/4.1.2/ /usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../../../x86_64-suse-linux/lib/x86_64-suse-linux/4.1.2/ /usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../../../x86_64-suse-linux/lib/../lib64/ /usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../../x86_64-suse-linux/4.1.2/ /usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../../../lib64/ /lib/x86_64-suse-linux/4.1.2/ /lib/../lib64/ /usr/lib/x86_64-suse-linux/4.1.2/ /usr/lib/../lib64/ /fastfs/wir/local/lib/ /usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../../../x86_64-suse-linux/lib/ /usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../../ /lib/ /usr/lib/"
 
 # Run-time system search path for libraries
-sys_lib_dlsearch_path_spec="/lib /usr/lib /usr/lib/mysql /usr/lib/octave-2.9.9 /usr/lib/qt-3.3/lib /usr/lib/qt4/lib /usr/lib/xulrunner-1.9.2 "
+sys_lib_dlsearch_path_spec="/lib /usr/lib /usr/X11R6/lib64/Xaw3d /usr/X11R6/lib64 /usr/X11R6/lib/Xaw3d /usr/X11R6/lib /usr/x86_64-suse-linux/lib /usr/local/lib64 /usr/local/lib /opt/kde3/lib64 /opt/kde3/lib /opt/gnome/lib64 /opt/gnome/lib /lib64 /lib /usr/lib64 /usr/lib /opt/cluster/intel/cce/9.1.042/lib /opt/cluster/intel/fce/9.1.036/lib /opt/cluster/Pathscale3.0/lib/2.9.99 /opt/cluster/Pathscale3.0/lib/2.9.99/32 /work/licsoft/compilers/pgi/linux86-64/6.2/lib /work/licsoft/compilers/pgi/linux86-64/6.2/libso "
 
 # Fix the shell variable $srcfile for the compiler.
 fix_srcfile_path=""
@@ -6760,7 +6760,7 @@ build_old_libs=`case $build_libtool_libs in yes) $echo no;; *) $echo yes;; esac`
 # End:
 # ### BEGIN LIBTOOL TAG CONFIG: CXX
 
-# Libtool was configured on host NWRW15:
+# Libtool was configured on host deimos101:
 
 # Shell to use when invoking shell scripts.
 SHELL="/bin/sh"
@@ -6782,12 +6782,12 @@ fast_install=yes
 
 # The host system.
 host_alias=
-host=i686-pc-linux-gnu
+host=x86_64-unknown-linux-gnu
 host_os=linux-gnu
 
 # The build system.
 build_alias=
-build=i686-pc-linux-gnu
+build=x86_64-unknown-linux-gnu
 build_os=linux-gnu
 
 # An echo program that does not interpret backslashes.
@@ -6798,13 +6798,13 @@ AR="ar"
 AR_FLAGS="cru"
 
 # A C compiler.
-LTCC="gcc"
+LTCC="/licsoft/libraries/openmpi/1.2.6/64bit/bin/mpicc"
 
 # LTCC compiler flags.
 LTCFLAGS="-g -O2"
 
 # A language-specific compiler.
-CC="g++"
+CC="/licsoft/libraries/openmpi/1.2.6/64bit/bin/mpicxx"
 
 # Is the compiler the GNU C compiler?
 with_gcc=yes
@@ -6813,7 +6813,7 @@ with_gcc=yes
 EGREP="grep -E"
 
 # The linker used to build libraries.
-LD="/usr/bin/ld"
+LD="/usr/x86_64-suse-linux/bin/ld -m elf_x86_64"
 
 # Whether we need hard or soft links.
 LN_S="ln -s"
@@ -6887,7 +6887,7 @@ dlopen_self=unknown
 dlopen_self_static=unknown
 
 # Compiler flag to prevent dynamic linking.
-link_static_flag="-static"
+link_static_flag=""
 
 # Compiler flag to turn off builtin functions.
 no_builtin_flag=" -fno-builtin"
@@ -6942,11 +6942,11 @@ striplib="strip --strip-unneeded"
 
 # Dependencies to place before the objects being linked to create a
 # shared library.
-predep_objects="/usr/lib/gcc/i386-redhat-linux/4.1.2/../../../crti.o /usr/lib/gcc/i386-redhat-linux/4.1.2/crtbeginS.o"
+predep_objects="/usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../../../lib64/crti.o /usr/lib64/gcc/x86_64-suse-linux/4.1.2/crtbeginS.o"
 
 # Dependencies to place after the objects being linked to create a
 # shared library.
-postdep_objects="/usr/lib/gcc/i386-redhat-linux/4.1.2/crtendS.o /usr/lib/gcc/i386-redhat-linux/4.1.2/../../../crtn.o"
+postdep_objects="/usr/lib64/gcc/x86_64-suse-linux/4.1.2/crtendS.o /usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../../../lib64/crtn.o"
 
 # Dependencies to place before the objects being linked to create a
 # shared library.
@@ -6954,11 +6954,11 @@ predeps=""
 
 # Dependencies to place after the objects being linked to create a
 # shared library.
-postdeps="-lstdc++ -lm -lgcc_s -lc -lgcc_s"
+postdeps="-lmpi_cxx -lmpi -lopen-rte -lopen-pal -libverbs -lrt -lnuma -ldl -lnsl -lutil -ldl -lstdc++ -lm -lgcc_s -lpthread -lc -lgcc_s"
 
 # The library search path used internally by the compiler when linking
 # a shared library.
-compiler_lib_search_path="-L/u/witkowski/local/lib -L/u/witkowski/local/intel/mkl/10.0.1.014/lib/32 -L/usr/lib/gcc/i386-redhat-linux/4.1.2 -L/usr/lib/gcc/i386-redhat-linux/4.1.2 -L/usr/lib/gcc/i386-redhat-linux/4.1.2/../../.."
+compiler_lib_search_path="-L/usr/lib64 -L/licsoft/libraries/openmpi/1.2.6/64bit/lib -L/usr/lib64/gcc/x86_64-suse-linux/4.1.2 -L/usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../../../lib64 -L/lib/../lib64 -L/usr/lib/../lib64 -L/fastfs/wir/local/lib -L/usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../../../x86_64-suse-linux/lib -L/usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../.."
 
 # Method to check whether dependent libraries are shared objects.
 deplibs_check_method="pass_all"
@@ -7038,10 +7038,10 @@ variables_saved_for_relink="PATH LD_LIBRARY_PATH LD_RUN_PATH GCC_EXEC_PREFIX COM
 link_all_deplibs=unknown
 
 # Compile-time system search path for libraries
-sys_lib_search_path_spec=" /u/witkowski/local/lib/i386-redhat-linux/4.1.2/ /u/witkowski/local/lib/ /u/witkowski/local/intel/mkl/10.0.1.014/lib/32/i386-redhat-linux/4.1.2/ /u/witkowski/local/intel/mkl/10.0.1.014/lib/32/ /usr/lib/gcc/i386-redhat-linux/4.1.2/ /usr/lib/gcc/i386-redhat-linux/4.1.2/ /usr/lib/gcc/i386-redhat-linux/4.1.2/../../../../i386-redhat-linux/lib/i386-redhat-linux/4.1.2/ /usr/lib/gcc/i386-redhat-linux/4.1.2/../../../../i386-redhat-linux/lib/ /usr/lib/gcc/i386-redhat-linux/4.1.2/../../../i386-redhat-linux/4.1.2/ /usr/lib/gcc/i386-redhat-linux/4.1.2/../../../ /lib/i386-redhat-linux/4.1.2/ /lib/ /usr/lib/i386-redhat-linux/4.1.2/ /usr/lib/"
+sys_lib_search_path_spec=" /fastfs/wir/local/lib/x86_64-suse-linux/4.1.2/ /fastfs/wir/local/lib/../lib64/ /usr/lib64/gcc/x86_64-suse-linux/4.1.2/ /usr/lib/gcc/x86_64-suse-linux/4.1.2/ /usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../../../x86_64-suse-linux/lib/x86_64-suse-linux/4.1.2/ /usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../../../x86_64-suse-linux/lib/../lib64/ /usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../../x86_64-suse-linux/4.1.2/ /usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../../../lib64/ /lib/x86_64-suse-linux/4.1.2/ /lib/../lib64/ /usr/lib/x86_64-suse-linux/4.1.2/ /usr/lib/../lib64/ /fastfs/wir/local/lib/ /usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../../../x86_64-suse-linux/lib/ /usr/lib64/gcc/x86_64-suse-linux/4.1.2/../../../ /lib/ /usr/lib/"
 
 # Run-time system search path for libraries
-sys_lib_dlsearch_path_spec="/lib /usr/lib /usr/lib/mysql /usr/lib/octave-2.9.9 /usr/lib/qt-3.3/lib /usr/lib/qt4/lib /usr/lib/xulrunner-1.9.2 "
+sys_lib_dlsearch_path_spec="/lib /usr/lib /usr/X11R6/lib64/Xaw3d /usr/X11R6/lib64 /usr/X11R6/lib/Xaw3d /usr/X11R6/lib /usr/x86_64-suse-linux/lib /usr/local/lib64 /usr/local/lib /opt/kde3/lib64 /opt/kde3/lib /opt/gnome/lib64 /opt/gnome/lib /lib64 /lib /usr/lib64 /usr/lib /opt/cluster/intel/cce/9.1.042/lib /opt/cluster/intel/fce/9.1.036/lib /opt/cluster/Pathscale3.0/lib/2.9.99 /opt/cluster/Pathscale3.0/lib/2.9.99/32 /work/licsoft/compilers/pgi/linux86-64/6.2/lib /work/licsoft/compilers/pgi/linux86-64/6.2/libso "
 
 # Fix the shell variable $srcfile for the compiler.
 fix_srcfile_path=""
@@ -7065,7 +7065,7 @@ include_expsyms=""
 
 # ### BEGIN LIBTOOL TAG CONFIG: F77
 
-# Libtool was configured on host NWRW15:
+# Libtool was configured on host deimos101:
 
 # Shell to use when invoking shell scripts.
 SHELL="/bin/sh"
@@ -7087,12 +7087,12 @@ fast_install=yes
 
 # The host system.
 host_alias=
-host=i686-pc-linux-gnu
+host=x86_64-unknown-linux-gnu
 host_os=linux-gnu
 
 # The build system.
 build_alias=
-build=i686-pc-linux-gnu
+build=x86_64-unknown-linux-gnu
 build_os=linux-gnu
 
 # An echo program that does not interpret backslashes.
@@ -7103,7 +7103,7 @@ AR="ar"
 AR_FLAGS="cru"
 
 # A C compiler.
-LTCC="gcc"
+LTCC="/licsoft/libraries/openmpi/1.2.6/64bit/bin/mpicc"
 
 # LTCC compiler flags.
 LTCFLAGS="-g -O2"
@@ -7112,13 +7112,13 @@ LTCFLAGS="-g -O2"
 CC="g77"
 
 # Is the compiler the GNU C compiler?
-with_gcc=yes
+with_gcc=
 
 # An ERE matcher.
 EGREP="grep -E"
 
 # The linker used to build libraries.
-LD="/usr/bin/ld"
+LD="/usr/x86_64-suse-linux/bin/ld -m elf_x86_64"
 
 # Whether we need hard or soft links.
 LN_S="ln -s"
@@ -7346,10 +7346,10 @@ variables_saved_for_relink="PATH LD_LIBRARY_PATH LD_RUN_PATH GCC_EXEC_PREFIX COM
 link_all_deplibs=unknown
 
 # Compile-time system search path for libraries
-sys_lib_search_path_spec=" /u/witkowski/local/lib/i386-redhat-linux/3.4.6/ /u/witkowski/local/lib/ /u/witkowski/local/intel/mkl/10.0.1.014/lib/32/i386-redhat-linux/3.4.6/ /u/witkowski/local/intel/mkl/10.0.1.014/lib/32/ /usr/lib/gcc/i386-redhat-linux/3.4.6/ /usr/lib/gcc/i386-redhat-linux/3.4.6/ /usr/lib/gcc/i386-redhat-linux/3.4.6/../../../../i386-redhat-linux/lib/i386-redhat-linux/3.4.6/ /usr/lib/gcc/i386-redhat-linux/3.4.6/../../../../i386-redhat-linux/lib/ /usr/lib/gcc/i386-redhat-linux/3.4.6/../../../i386-redhat-linux/3.4.6/ /usr/lib/gcc/i386-redhat-linux/3.4.6/../../../ /lib/i386-redhat-linux/3.4.6/ /lib/ /usr/lib/i386-redhat-linux/3.4.6/ /usr/lib/"
+sys_lib_search_path_spec=" /fastfs/wir/local/lib/x86_64-suse-linux/3.3.5/ /fastfs/wir/local/lib/ /usr/lib64/gcc-lib/x86_64-suse-linux/3.3.5/ /usr/lib/gcc/x86_64-suse-linux/3.3.5/ /usr/lib64/gcc-lib/x86_64-suse-linux/3.3.5/../../../../x86_64-suse-linux/lib/x86_64-suse-linux/3.3.5/ /usr/lib64/gcc-lib/x86_64-suse-linux/3.3.5/../../../../x86_64-suse-linux/lib/ /usr/lib64/gcc-lib/x86_64-suse-linux/3.3.5/../../../x86_64-suse-linux/3.3.5/ /usr/lib64/gcc-lib/x86_64-suse-linux/3.3.5/../../../ /lib/x86_64-suse-linux/3.3.5/ /lib/ /usr/lib/x86_64-suse-linux/3.3.5/ /usr/lib/"
 
 # Run-time system search path for libraries
-sys_lib_dlsearch_path_spec="/lib /usr/lib /usr/lib/mysql /usr/lib/octave-2.9.9 /usr/lib/qt-3.3/lib /usr/lib/qt4/lib /usr/lib/xulrunner-1.9.2 "
+sys_lib_dlsearch_path_spec="/lib /usr/lib /usr/X11R6/lib64/Xaw3d /usr/X11R6/lib64 /usr/X11R6/lib/Xaw3d /usr/X11R6/lib /usr/x86_64-suse-linux/lib /usr/local/lib64 /usr/local/lib /opt/kde3/lib64 /opt/kde3/lib /opt/gnome/lib64 /opt/gnome/lib /lib64 /lib /usr/lib64 /usr/lib /opt/cluster/intel/cce/9.1.042/lib /opt/cluster/intel/fce/9.1.036/lib /opt/cluster/Pathscale3.0/lib/2.9.99 /opt/cluster/Pathscale3.0/lib/2.9.99/32 /work/licsoft/compilers/pgi/linux86-64/6.2/lib /work/licsoft/compilers/pgi/linux86-64/6.2/libso "
 
 # Fix the shell variable $srcfile for the compiler.
 fix_srcfile_path=""
diff --git a/AMDiS/src/CoarseningManager.cc b/AMDiS/src/CoarseningManager.cc
index 1b31a92ab62d798e76e9871048a27178699c8d32..5bae9d5626c7c6fdcfd03cc2c49c8775420097f7 100644
--- a/AMDiS/src/CoarseningManager.cc
+++ b/AMDiS/src/CoarseningManager.cc
@@ -33,7 +33,7 @@ namespace AMDiS {
 
       if (el->getChild(0)) {
 	// interior node of the tree
-	signed char mark = max(el->getChild(0)->getMark(), el->getChild(1)->getMark());
+	int mark = max(el->getChild(0)->getMark(), el->getChild(1)->getMark());
 	el->setMark(std::min(mark + 1, 0));
       } else {
 	// leaf node of the tree                                 
diff --git a/AMDiS/src/Debug.cc b/AMDiS/src/Debug.cc
index f73cd02d343c4cf4ea601d0a4d07f182c1cc187c..7c4513bb7bbd8f5c9b91c0bcac00698eb8e687c2 100644
--- a/AMDiS/src/Debug.cc
+++ b/AMDiS/src/Debug.cc
@@ -62,11 +62,16 @@ namespace AMDiS {
     }
 
 
-    void writeElementIndexMesh(Mesh *mesh, std::string filename)
+    void writeElementIndexMesh(Mesh *mesh, std::string filename, int level)
     {
+      FUNCNAME("debug::writeElementIndexMesh()");
+
       std::map<int, double> vec;    
       TraverseStack stack;
-      ElInfo *elInfo = stack.traverseFirst(mesh, -1, Mesh::CALL_LEAF_EL);
+      ElInfo *elInfo = 
+	level == -1 ? 
+	stack.traverseFirst(mesh, -1, Mesh::CALL_LEAF_EL) :
+	stack.traverseFirst(mesh, level, Mesh::CALL_EL_LEVEL);
       
       while (elInfo) {		  
 	int index = elInfo->getElement()->getIndex();
@@ -74,7 +79,7 @@ namespace AMDiS {
 	elInfo = stack.traverseNext(elInfo);
       }
 
-      ElementFileWriter::writeFile(vec, mesh, filename);
+      ElementFileWriter::writeFile(vec, mesh, filename, level);
     }
 
 
diff --git a/AMDiS/src/Debug.h b/AMDiS/src/Debug.h
index bdc7cfd4f6bc64bd59d28cea6800b7974f3942e1..a965f658dbb3cac321b193c3227c97f3468a671b 100644
--- a/AMDiS/src/Debug.h
+++ b/AMDiS/src/Debug.h
@@ -71,8 +71,10 @@ namespace AMDiS {
      *
      * \param[in]  feSpace   The FE space to be used.
      * \param[in]  filename  Name of the file.
+     * \param[in]  level     If level is -1, all leaf elements will be put to the
+     *                       output file, otherwise the elements with the given level.
      */
-    void writeElementIndexMesh(Mesh *mesh, std::string filename);
+    void writeElementIndexMesh(Mesh *mesh, std::string filename, int level = -1);
 
     void highlightElementIndexMesh(Mesh *mesh, int idx, std::string filename);
 
diff --git a/AMDiS/src/ElInfo2d.cc b/AMDiS/src/ElInfo2d.cc
index 7839cfc43327c7db395b9ba5ac3e1477a2b9ad8c..af832bc1664686b95237e4430419d4150bc3429c 100644
--- a/AMDiS/src/ElInfo2d.cc
+++ b/AMDiS/src/ElInfo2d.cc
@@ -411,6 +411,11 @@ namespace AMDiS {
 	
 	Element *nb = elInfoOld->neighbour[2];
 	if (nb) {
+// 	  MSG("OLDINFO HAS EL %d mit PTR %p\n", elInfoOld->getElement()->getIndex(),
+// 	      elInfoOld->getElement());
+
+// 	  MSG("1 ALS NACHBAR HABEN WIR %d mit PTR %p\n", nb->getIndex(), nb);
+
 	  TEST_EXIT_DBG(elInfoOld->oppVertex[2] == 2)
 	    ("Fill child %d of element %d (mel %d): Invalid neighbour %d!\n",
 	     ichild,
@@ -424,6 +429,8 @@ namespace AMDiS {
 	    ("Missing second child in element %d!\n", nb->getIndex());
 	 
 	  nb = nb->getSecondChild();
+	  // 	  MSG("2 ALS NACHBAR HABEN WIR %d mit PTR %p (is LEAF: %d)\n", 
+	  //	  nb->getIndex(), nb, nb->isLeaf());
 
 	  if (nb->getFirstChild()) {
 	    oppVertex[0] = 2;
@@ -443,6 +450,7 @@ namespace AMDiS {
 	    }	   
  
 	    nb = nb->getFirstChild();
+	    //	    MSG("3 ALS NACHBAR HABEN WIR %d\n", nb->getIndex());
 	  } else {
 	    oppVertex[0] = 1;
 
diff --git a/AMDiS/src/Element.cc b/AMDiS/src/Element.cc
index cb2bf142b6822abba27d60dd28115c998bf33296..7ee7089973fb239888394c696504c02210c87e62 100644
--- a/AMDiS/src/Element.cc
+++ b/AMDiS/src/Element.cc
@@ -21,6 +21,7 @@ namespace AMDiS {
     newCoord = NULL;
     elementData = NULL;
     mark = 0;
+    dofValid = true;
 
     if (mesh)
       createNewDofPtrs();
@@ -125,6 +126,7 @@ namespace AMDiS {
       child[1]->deleteElementDOFs();
   }
 
+
   Element* Element::cloneWithDOFs()
   {
     Element *el;
@@ -192,6 +194,7 @@ namespace AMDiS {
     return el;
   }
 
+
   /****************************************************************************/
   /*  ATTENTION:                                                              */
   /*  new_dof_fct() destroys new_dof !!!!!!!!!!                               */
@@ -324,9 +327,7 @@ namespace AMDiS {
 
       for (int j = 0; j < dim; j++) {
 	if (dof[i] == pdof[j]) {
-	  /****************************************************************************/
-	  /* i is a common vertex                                                     */
-	  /****************************************************************************/
+	  // i is a common vertex 
 	  ov += i;
 	  nv++;
 	  break;
@@ -337,10 +338,8 @@ namespace AMDiS {
     
     if (nv != mesh->getDim()) 
       return(-1);
-    /****************************************************************************/
-    /*  the opposite vertex is 3(6) - (sum of indices of common vertices) in    */
-    /*  2d(3d)                                                                  */
-    /****************************************************************************/
+
+    // The opposite vertex is 3(6) - (sum of indices of common vertices) in 2D(3D)
 
     switch(mesh->getDim()) {
     case 1:
@@ -384,39 +383,42 @@ namespace AMDiS {
     int nodes = mesh->getNumberOfNodes();
     int j = 0;
     SerUtil::serialize(out, nodes);
+    SerUtil::serialize(out, dofValid);
    
-    for (int pos = 0; pos <= dim; pos++) {
-      GeoIndex position = INDEX_OF_DIM(pos, dim);
-      int ndof = 0;
-
-      for (int i = 0; i < mesh->getNumberOfDOFAdmin(); i++)
-	ndof += mesh->getDofAdmin(i).getNumberOfDofs(position);	      
-
-      if (ndof > 0) {
-	for (int i = 0; i < mesh->getGeo(position); i++) {
-	  if (dof[j] != NULL) {
-	    // Create index to check if the dofs were already written.
-	    std::pair<DegreeOfFreedom, int> idx = std::make_pair(dof[j][0], pos);
-
-	    if (Mesh::serializedDOFs[idx] == NULL) {
-	      Mesh::serializedDOFs[idx] = dof[j];
-	      SerUtil::serialize(out, ndof);
-	      SerUtil::serialize(out, pos);
-	      out.write(reinterpret_cast<const char*>(dof[j]), 
-			ndof * sizeof(DegreeOfFreedom));
+    if (dofValid) {
+      for (int pos = 0; pos <= dim; pos++) {
+	GeoIndex position = INDEX_OF_DIM(pos, dim);
+	int ndof = 0;
+	
+	for (int i = 0; i < mesh->getNumberOfDOFAdmin(); i++)
+	  ndof += mesh->getDofAdmin(i).getNumberOfDofs(position);	      
+	
+	if (ndof > 0) {
+	  for (int i = 0; i < mesh->getGeo(position); i++) {
+	    if (dof[j] != NULL) {
+	      // Create index to check if the dofs were already written.
+	      std::pair<DegreeOfFreedom, int> idx = std::make_pair(dof[j][0], pos);
+	      
+	      if (Mesh::serializedDOFs[idx] == NULL) {
+		Mesh::serializedDOFs[idx] = dof[j];
+		SerUtil::serialize(out, ndof);
+		SerUtil::serialize(out, pos);
+		out.write(reinterpret_cast<const char*>(dof[j]), 
+			  ndof * sizeof(DegreeOfFreedom));
+	      } else {
+		int minusOne = -1;
+		SerUtil::serialize(out, minusOne);
+		SerUtil::serialize(out, pos);
+		out.write(reinterpret_cast<const char*>(&(dof[j][0])), 
+			  sizeof(DegreeOfFreedom));
+	      }
 	    } else {
-	      int minusOne = -1;
-	      SerUtil::serialize(out, minusOne);
+	      int zero = 0;
+	      SerUtil::serialize(out, zero);
 	      SerUtil::serialize(out, pos);
-	      out.write(reinterpret_cast<const char*>(&(dof[j][0])), 
-			sizeof(DegreeOfFreedom));
 	    }
-	  } else {
-	    int zero = 0;
-	    SerUtil::serialize(out, zero);
-	    SerUtil::serialize(out, pos);
+	    j++;
 	  }
-	  j++;
 	}
       }
     }
@@ -478,46 +480,48 @@ namespace AMDiS {
     // read dofs
     int nodes;
     SerUtil::deserialize(in, nodes);
+    SerUtil::deserialize(in, dofValid);
 
-    TEST_EXIT_DBG(nodes > 0)("Should not happen!\n");
-
+    TEST_EXIT_DBG(nodes > 0)("Should not happen!\n");      
     dof = new DegreeOfFreedom*[nodes]; 
-
-    for (int i = 0; i < nodes; i++) {
-      int nDofs, pos;
-      SerUtil::deserialize(in, nDofs);
-      SerUtil::deserialize(in, pos);
-
-      if (nDofs) {
-	if (nDofs != -1) {
-	  dof[i] = new DegreeOfFreedom[nDofs];
-	  in.read(reinterpret_cast<char*>(dof[i]), nDofs * sizeof(DegreeOfFreedom));
-
-	  // Create index to check if the dofs were alread read from file.
-	  std::pair<DegreeOfFreedom, int> idx = std::make_pair(dof[i][0], pos);
-
-	  if (Mesh::serializedDOFs[idx] != NULL) {
-	    DegreeOfFreedom *dofPtr = Mesh::serializedDOFs[idx];
-	    delete [] dof[i];
-	    dof[i] = dofPtr;
+     
+    if (dofValid) {
+      for (int i = 0; i < nodes; i++) {
+	int nDofs, pos;
+	SerUtil::deserialize(in, nDofs);
+	SerUtil::deserialize(in, pos);
+	
+	if (nDofs) {
+	  if (nDofs != -1) {
+	    dof[i] = new DegreeOfFreedom[nDofs];
+	    in.read(reinterpret_cast<char*>(dof[i]), nDofs * sizeof(DegreeOfFreedom));
+	    
+	    // Create index to check if the dofs were alread read from file.
+	    std::pair<DegreeOfFreedom, int> idx = std::make_pair(dof[i][0], pos);
+	    
+	    if (Mesh::serializedDOFs[idx] != NULL) {
+	      DegreeOfFreedom *dofPtr = Mesh::serializedDOFs[idx];
+	      delete [] dof[i];
+	      dof[i] = dofPtr;
+	    } else {
+	      Mesh::serializedDOFs[idx] = dof[i];
+	    }
+	    
 	  } else {
-	    Mesh::serializedDOFs[idx] = dof[i];
+	    DegreeOfFreedom index;
+	    SerUtil::deserialize(in, index);
+	    
+	    std::pair<DegreeOfFreedom, int> idx = std::make_pair(index, pos);
+	    TEST_EXIT(Mesh::serializedDOFs.find(idx) !=  Mesh::serializedDOFs.end())
+	      ("This should never happen!\n");
+	    dof[i] = Mesh::serializedDOFs[idx];
 	  }
-
 	} else {
-	  DegreeOfFreedom index;
-	  SerUtil::deserialize(in, index);
-      
-	  std::pair<DegreeOfFreedom, int> idx = std::make_pair(index, pos);
-	  TEST_EXIT(Mesh::serializedDOFs.find(idx) !=  Mesh::serializedDOFs.end())
-	    ("This should never happen!\n");
-	  dof[i] = Mesh::serializedDOFs[idx];
+	  dof[i] = NULL;
 	}
-      } else {
-	dof[i] = NULL;
       }
     }
-   
+      
     // read index
     SerUtil::deserialize(in, index);
 
@@ -556,20 +560,6 @@ namespace AMDiS {
   }
 
 
-  int Element::calcMemoryUsage()
-  {
-    int result = 0;
-
-    result += sizeof(Element);
-    result += mesh->getNumberOfNodes() * sizeof(DegreeOfFreedom*);
-
-    if (child[0])
-      result += child[0]->calcMemoryUsage() + child[1]->calcMemoryUsage();    
-
-    return result;
-  }
-
-
   void writeDotFile(Element *el, std::string filename, int maxLevels)
   {
     std::vector<int> graphNodes;
diff --git a/AMDiS/src/Element.h b/AMDiS/src/Element.h
index 21f185ced2be188e4ffe1d8043dbf738b7f45441..05f9e57da98bf0b425f0951440175d9635c5acad 100644
--- a/AMDiS/src/Element.h
+++ b/AMDiS/src/Element.h
@@ -120,20 +120,26 @@ namespace AMDiS {
     }
 
     /// Returns \ref dof[i][j] which is the j-th DOF of the i-th node of Element.
-    const DegreeOfFreedom getDof(int i, int j) const 
+    inline DegreeOfFreedom getDof(int i, int j) const 
     { 
+      TEST_EXIT_DBG(dofValid)("Should not happen!\n");
+
       return dof[i][j];
     }
 
     /// Returns \ref dof[i] which is a pointer to the DOFs of the i-th node.
-    const DegreeOfFreedom* getDof(int i) const 
+    inline const DegreeOfFreedom* getDof(int i) const 
     {
+      TEST_EXIT_DBG(dofValid)("Should not happen!\n");
+
       return dof[i];
     }
 
     /// Returns a pointer to the DOFs of this Element
-    const DegreeOfFreedom** getDof() const 
+    inline const DegreeOfFreedom** getDof() const 
     {
+      TEST_EXIT_DBG(dofValid)("Should not happen!\n");
+
       return const_cast<const DegreeOfFreedom**>(dof);
     }
 
@@ -143,6 +149,11 @@ namespace AMDiS {
       return mesh; 
     }
 
+    inline void setDofValid(bool b)
+    {
+      dofValid = b;
+    }
+
     /** \brief
      * Returns \ref elementData's error estimation, if Element is a leaf element
      * and has leaf data. 
@@ -207,7 +218,7 @@ namespace AMDiS {
     virtual int getGeo(GeoIndex i) const = 0;
 
     /// Returns Element's \ref mark
-    inline const signed char getMark() const 
+    inline const int getMark() const 
     { 
       return mark;
     }
@@ -328,7 +339,7 @@ namespace AMDiS {
     }
 
     /// Sets Element's \ref mark
-    inline void setMark(signed char m) 
+    inline void setMark(int m) 
     {
       mark = m;
     }
@@ -509,9 +520,6 @@ namespace AMDiS {
     /// Deserialize an element from a file.
     void deserialize(std::istream &in);
 
-    ///
-    int calcMemoryUsage();
-
     /// Sets Element's \ref dof pointer.
     void createNewDofPtrs();
  
@@ -565,7 +573,7 @@ namespace AMDiS {
      * element, this element is refined mark times. if mark is negative for
      * a leaf element, this element is coarsened -mark times.
      */
-    signed char mark;
+    int mark;
  
     /** \brief
      * If the element has a boundary edge on a curved boundary, this is a pointer
@@ -582,6 +590,14 @@ namespace AMDiS {
     /// Pointer to Element's leaf data
     ElementData* elementData;
 
+    /** \brief
+     * Is true, if the DOF pointers are valid. In sequential computations this
+     * is always the case. In parallel computations all macro elements are stored
+     * in memory though not all of them are part of the mesh (because they are owned
+     * by another rank). In this case, there are no valid DOFs on the element. 
+     */
+    bool dofValid;
+
     /** \brief
      * This map is used for deletion of all DOFs of all elements of a mesh. Once
      * a DOF-vector (all DOFS at a node, edge, etc.) is deleted, its address is
diff --git a/AMDiS/src/MacroElement.cc b/AMDiS/src/MacroElement.cc
index 49cef62516f875620ee02e387d10bbeea31aff6f..fb23bc693c22f341554efc4be2f212097bf6b200 100644
--- a/AMDiS/src/MacroElement.cc
+++ b/AMDiS/src/MacroElement.cc
@@ -25,12 +25,14 @@ namespace AMDiS {
     projection.set(NULL);
   }
 
+
   MacroElement::~MacroElement()
   {
     if (element)
       delete element;    
   }
 
+
   MacroElement& MacroElement::operator=(const MacroElement &el)
   {
     FUNCNAME("MacroElement::operator=()");
@@ -48,6 +50,7 @@ namespace AMDiS {
     return *this;
   }
 
+
   void MacroElement::serialize(std::ostream &out)
   {
     // write element-tree
@@ -90,6 +93,7 @@ namespace AMDiS {
     SerUtil::serialize(out, elType);
   }
 
+
   void MacroElement::deserialize(std::istream &in)
   {
     FUNCNAME("MacroElement::deserialize()");
@@ -105,7 +109,7 @@ namespace AMDiS {
     } else {
       if (typeName == "Line") 
 	element = new Line(NULL);
-      if (typeName == "Triangle") 
+      if (typeName == "Triangle")
 	element = new Triangle(NULL);
       if (typeName == "Tetrahedron") 
 	element = new Tetrahedron(NULL);
@@ -169,18 +173,4 @@ namespace AMDiS {
     SerUtil::deserialize(in, elType);
   }
 
-  int MacroElement::calcMemoryUsage() 
-  {
-    int result = 0;
-
-    result += sizeof(MacroElement);
-    result += coord.getSize() * sizeof(WorldVector<double>);
-    result += projection.getSize() * sizeof(Projection*);
-    result += neighbour.getSize() * sizeof(MacroElement*);
-    result += oppVertex.getSize() * sizeof(char);
-    result += element->calcMemoryUsage();
-
-    return result;
-  }
-
 }
diff --git a/AMDiS/src/MacroElement.h b/AMDiS/src/MacroElement.h
index c3f6175b3d1fd0c7f51751b0215e1a65743eb320..284f66938c496a366af32aa0d9e2758745644cc2 100644
--- a/AMDiS/src/MacroElement.h
+++ b/AMDiS/src/MacroElement.h
@@ -80,7 +80,7 @@ namespace AMDiS {
     }
 
     /// Returns the i-th opp-vertex of this MacroElement \ref oppVertex[i]
-    inline char getOppVertex(int i) const 
+    inline int getOppVertex(int i) const 
     {
       return oppVertex[i];
     }
@@ -151,7 +151,7 @@ namespace AMDiS {
     }
 
     /// Sets the i-th opp vertex to c
-    inline void  setOppVertex(int i, char c) 
+    inline void setOppVertex(int i, int c)
     {
       oppVertex[i] = c;
     }
@@ -183,9 +183,6 @@ namespace AMDiS {
       deserializedNeighbourIndices = indices;
     }
 
-    ///
-    int calcMemoryUsage();
-
   protected:
     /// Element of this MacroElement.
     Element *element;
@@ -203,7 +200,7 @@ namespace AMDiS {
     FixVec<MacroElement*, NEIGH> neighbour;
 
     /// opp vertices of this MacroElement
-    FixVec<char, NEIGH> oppVertex;
+    FixVec<int, NEIGH> oppVertex;
 
     /// index of this MacroElement
     int index;  
diff --git a/AMDiS/src/Mesh.cc b/AMDiS/src/Mesh.cc
index d530bb076876700996de846c9f09fecee3b207f6..ee3592199d1017386a7c6b53d8e269f32da00518 100644
--- a/AMDiS/src/Mesh.cc
+++ b/AMDiS/src/Mesh.cc
@@ -305,8 +305,9 @@ namespace AMDiS {
       // member of the list of macro elements to be deleted, is is inserted to the new
       // macro element list.
       if (macros.find(*elIter) == macros.end())
-	newMacroElements.push_back(*elIter);      
+	newMacroElements.push_back(*elIter);     
     }
+
     // And replace the macro element list with the new one.
     macroElements.clear();
     macroElements = newMacroElements;
@@ -337,6 +338,11 @@ namespace AMDiS {
 
 	(*macroIt)->getElement()->setElementData(elementDataPrototype->clone()); 
       }
+
+      // We will delete at least some element DOFs in the next but will keep the
+      // element object still in memory. Therefore the DOF pointer must be set to be
+      // invalide.
+      (*macroIt)->getElement()->setDofValid(false);
     }     
 
 
@@ -359,7 +365,7 @@ namespace AMDiS {
 
       if (deleteDof)
 	freeDof(const_cast<DegreeOfFreedom*>(dofsIt->first), 
-		dofsPosIndex[dofsIt->first]);
+		dofsPosIndex[dofsIt->first]);      
     }
 
 
@@ -1101,9 +1107,6 @@ namespace AMDiS {
       elIndexVecIndex[macroElements[i]->getIndex()] = i;
     }
 
-    SerUtil::deserialize(in, elementIndex);
-    SerUtil::deserialize(in, initialized);
-
     // set neighbour pointer in macro elements
     int nNeighbour = getGeo(NEIGH);
     for (int i = 0; i < static_cast<int>(macroElements.size()); i++) {
@@ -1133,7 +1136,11 @@ namespace AMDiS {
 
     serializedDOFs.clear();
 
-    
+
+    SerUtil::deserialize(in, elementIndex);
+    SerUtil::deserialize(in, initialized);
+
+   
     /// === Read periodic assoications. ===
 
     int mapSize = 0;
@@ -1155,6 +1162,8 @@ namespace AMDiS {
   {
     FUNCNAME("Mesh::initialize()");
 
+    TEST_EXIT(admin.size() > 0)("No DOF admin defined!\n");
+
     std::string macroFilename("");
     std::string valueFilename("");
     std::string periodicFilename("");
@@ -1186,14 +1195,8 @@ namespace AMDiS {
 	MacroReader::readMacro(macroFilename, this, periodicFilename, check);
 #endif
 
-      // In the sequentiell case, if there is no value file which should 
-      // be written, we can delete the information of the macro file. In
-      // parallel computations, the macro file infos are always needed for
-      // mesh repartitioning.
-#ifndef HAVE_PARALLEL_DOMAIN_AMDIS
       if (!valueFilename.length())
        	clearMacroFileInfo();
-#endif
     }
 
     initialized = true;
@@ -1432,9 +1435,6 @@ namespace AMDiS {
       result += admin[i]->getUsedSize() * sizeof(DegreeOfFreedom);
     }
     
-    for (int i = 0; i < static_cast<int>(macroElements.size()); i++)
-      result += macroElements[i]->calcMemoryUsage();    
-    
     return result;
   }
 
diff --git a/AMDiS/src/RefinementManager2d.cc b/AMDiS/src/RefinementManager2d.cc
index 02bd80a2f458a5f06a36e5aa1a0136acb34475ce..70d49dd711306e6cad914afb4ff110b38002c08b 100644
--- a/AMDiS/src/RefinementManager2d.cc
+++ b/AMDiS/src/RefinementManager2d.cc
@@ -239,17 +239,18 @@ namespace AMDiS {
   {
     FUNCNAME("RefinementManager2d::bisectTriangle()");
  
-    TEST_EXIT_DBG(mesh)("no mesh!\n");
- 
+    TEST_EXIT_DBG(mesh)("No mesh!\n");
+
     Triangle *child[2];
     child[0] = dynamic_cast<Triangle*>(mesh->createNewElement(el));
     child[1] = dynamic_cast<Triangle*>(mesh->createNewElement(el));
 
-    signed char newMark = max(0, el->getMark() - 1);
+    int newMark = max(0, el->getMark() - 1);
     child[0]->setMark(newMark);
     child[1]->setMark(newMark);
     el->setMark(0);
 
+
     // === Transfer hidden data from parent to children. ===
 
     // call of subclass-method
@@ -260,11 +261,13 @@ namespace AMDiS {
     if (newMark > 0) 
       doMoreRecursiveRefine = true;
 
+
     // === Vertex 2 is the newest vertex. ===
 
     child[0]->setDof(2, newDOFs[0]);
     child[1]->setDof(2, newDOFs[0]);
 
+
     // === The other vertices are handed on from the parent. ===
 
     for (int i_child = 0; i_child < 2; i_child++) {
@@ -272,6 +275,7 @@ namespace AMDiS {
       child[i_child]->setDof(1 - i_child, const_cast<int*>(el->getDof(i_child)));
     }
 
+
     // === There is one more leaf element, two hierachical elements and one ===
     // === more edge.                                                       ===
 
diff --git a/AMDiS/src/io/ElementFileWriter.cc b/AMDiS/src/io/ElementFileWriter.cc
index e50e8608373a62bf5fa27e6de82c07d2163f5481..d24054d2180d1380ad7eff7d3612599c4b25924a 100644
--- a/AMDiS/src/io/ElementFileWriter.cc
+++ b/AMDiS/src/io/ElementFileWriter.cc
@@ -10,10 +10,8 @@ namespace AMDiS {
 				       Mesh *mesh_,
 				       std::map<int, double> &mapvec)
     : name(name_),
-      tecplotExt(".plt"),
       amdisMeshDatExt(".elem.mesh"),
       vtkExt(".vtu"),
-      writeTecPlotFormat(0),
       writeAMDiSFormat(0),
       writeVtkFormat(0),
       appendIndex(0),
@@ -26,8 +24,6 @@ namespace AMDiS {
   {
     if (name != "") {
       GET_PARAMETER(0, name + "->output->filename", &filename);
-      GET_PARAMETER(0, name + "->output->TecPlot format", "%d", &writeTecPlotFormat);
-      GET_PARAMETER(0, name + "->output->TecPlot ext", &tecplotExt);
       GET_PARAMETER(0, name + "->output->AMDiS format", "%d", &writeAMDiSFormat);
       GET_PARAMETER(0, name + "->output->AMDiS mesh-dat ext", &amdisMeshDatExt);
       GET_PARAMETER(0, name + "->output->ParaView format", "%d", &writeVtkFormat);
@@ -38,6 +34,7 @@ namespace AMDiS {
     }
   }
 
+
   void ElementFileWriter::writeFiles(AdaptInfo *adaptInfo, bool force,
 				     int level, Flag traverseFlag, 
 				     bool (*writeElem)(ElInfo*))
@@ -66,18 +63,10 @@ namespace AMDiS {
       fn += timeStr;
     }
 
-
-    if (writeTecPlotFormat) {
-      TEST_EXIT(mesh)("no mesh\n");
-
-      writeTecPlotValues(const_cast<char*>((fn + tecplotExt).c_str()));
-      MSG("TecPlot file written to %s\n", (fn + tecplotExt).c_str());
-    }
-
     if (writeAMDiSFormat) {
       TEST_EXIT(mesh)("no mesh\n");
     
-      writeMeshDatValues(const_cast<char*>( (fn + amdisMeshDatExt).c_str()),
+      writeMeshDatValues(const_cast<char*>((fn + amdisMeshDatExt).c_str()),
 			 adaptInfo ? adaptInfo->getTime() : 0.0);
       MSG("MeshDat file written to %s\n", (fn + amdisMeshDatExt).c_str());
     }
@@ -85,7 +74,7 @@ namespace AMDiS {
     if (writeVtkFormat) {
       TEST_EXIT(mesh)("no mesh\n");
 
-      writeVtkValues(const_cast<char*>( (fn + vtkExt).c_str()));
+      writeVtkValues(const_cast<char*>((fn + vtkExt).c_str()));
       MSG("VTK file written to %s\n", (fn + vtkExt).c_str());		     
     }
   }
@@ -93,93 +82,14 @@ namespace AMDiS {
 
   void ElementFileWriter::writeFile(std::map<int, double> &vec,
 				    Mesh *mesh,
-				    std::string filename)
+				    std::string filename,
+				    int level)
   {
     ElementFileWriter efw("", mesh, vec);
-    efw.writeVtkValues(filename);
+    efw.writeVtkValues(filename, level);
   }
   
 
-  void ElementFileWriter::writeTecPlotValues(std::string filename)
-  {
-    FUNCNAME("ElementFileWriter::writeTecPlotValues()");
-    std::ofstream fout(filename.c_str());
-
-    TEST_EXIT(fout)("Could not open file %s !\n", filename.c_str());
-    fout.setf(std::ios::scientific,std::ios::floatfield);
-
-    int dim = mesh->getDim();
-    double val;
-
-
-    // === Write header. ===
-    fout << "TITLE = \"" << name.c_str() << "\"\n";
-    fout << "VARIABLES = ";
-    switch (dim) {
-    case 2: fout << "\"x\",\"y\"";
-      break;
-    case 3: fout << "\"x\",\"y\",\"z\"";
-      break;
-    default: ERROR_EXIT("illegal dimension !\n");
-      break;
-    }
-    fout << ",\"" << name.c_str() << "\"\n";
-    fout << "ZONE T=\"" << name.c_str() << "\""
-	 << ", N=" << 3*mesh->getNumberOfLeaves() 
-	 << ", E=" << mesh->getNumberOfLeaves() 
-	 << ", F=FEPOINT, ";
-    switch (dim) {
-    case 2: fout << "ET=TRIANGLE\n\n";
-      break;
-    case 3: fout << "ET=TETRAHEDRON\n\n";
-      break;
-    default: ERROR_EXIT("illegal dimension !\n");
-      break;
-    }
-
-
-    // === Write vertex coordinates and values (for each element !). ===
-    TraverseStack stack;
-
-    ElInfo *elInfo = stack.traverseFirst(mesh,
-					 -1, 
-					 Mesh::CALL_LEAF_EL | 
-					 Mesh::FILL_COORDS);
-
-    while (elInfo) {
-      // Get element value.
-      val = vec[elInfo->getElement()->getIndex()];
-    
-      // Write coordinates of all element vertices and element value.
-      for (int i = 0; i <= dim; ++i) {
-      	for (int j = 0; j < dim; ++j)
-	  fout << elInfo->getCoord(i)[j] << " ";
-	
-	fout << val << "\n";
-      }
-    
-      elInfo = stack.traverseNext(elInfo);
-    }  // end of: mesh traverse
-
-  
-    // === Write elements. ===
-    int numLeaves = mesh->getNumberOfLeaves();
-    int vertCntr = 0;
-    fout << "\n";
-
-    for (int i = 0; i<numLeaves; ++i) {
-      for (int j = 0; j<=dim; ++j) {
-	++vertCntr;
-	fout << vertCntr << " ";
-      }
-      fout << "\n";
-    }
-
-
-    fout.close();
-  }
-
-
   void ElementFileWriter::writeMeshDatValues(std::string filename, double time)
   {
     FUNCNAME("ElementFileWriter::writeMeshDatValues()");
@@ -280,7 +190,7 @@ namespace AMDiS {
   }
 
 
-  void ElementFileWriter::writeVtkValues(std::string filename)
+  void ElementFileWriter::writeVtkValues(std::string filename, int level)
   {
     FUNCNAME("ElementFileWriter::writeVtkValues()");
     std::ofstream fout(filename.c_str());
@@ -291,6 +201,18 @@ namespace AMDiS {
     int dimOfWorld = Global::getGeo(WORLD);
     int vertices = mesh->getGeo(VERTEX);
     int nElements = mesh->getNumberOfLeaves();
+    Flag traverseFlag = level == -1 ? Mesh::CALL_LEAF_EL : Mesh::CALL_EL_LEVEL;
+
+    if (level != -1) {
+      nElements = 0;
+
+      TraverseStack stack;
+      ElInfo *elInfo = stack.traverseFirst(mesh, level, Mesh::CALL_EL_LEVEL);      
+      while (elInfo) {
+	nElements++;
+	elInfo = stack.traverseNext(elInfo);
+      }
+    }
 
     fout << "<?xml version=\"1.0\"?>" << std::endl;
     fout << "<VTKFile type=\"UnstructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\">" << std::endl;
@@ -304,7 +226,7 @@ namespace AMDiS {
     // === Write vertex coordinates (every vertex for every element). ===
     TraverseStack stack;
     ElInfo *elInfo = 
-      stack.traverseFirst(mesh, -1, Mesh::CALL_LEAF_EL | Mesh::FILL_COORDS);
+      stack.traverseFirst(mesh, level, traverseFlag | Mesh::FILL_COORDS);
 
     while (elInfo) {
       // Write coordinates of all element vertices.
@@ -366,7 +288,7 @@ namespace AMDiS {
 
     fout.setf(std::ios::scientific,std::ios::floatfield);
 
-    elInfo = stack.traverseFirst(mesh, -1, Mesh::CALL_LEAF_EL | Mesh::FILL_COORDS);
+    elInfo = stack.traverseFirst(mesh, level, traverseFlag | Mesh::FILL_COORDS);
     int vc = 0;
     while (elInfo) {
       // Get element value.
diff --git a/AMDiS/src/io/ElementFileWriter.h b/AMDiS/src/io/ElementFileWriter.h
index 2b3a4164a3b3c06429fc19faeb3b72260b7ad024..0d8731ebd13f80c7fe1f8517741c9f25a466157a 100644
--- a/AMDiS/src/io/ElementFileWriter.h
+++ b/AMDiS/src/io/ElementFileWriter.h
@@ -15,7 +15,7 @@ namespace AMDiS {
    * It can be necessary to visualize data defined on elements, e.g. residual
    * error that is defined on elements and not on DOFs. This class takes as
    * input a mesh and a map, which defines for each element index a (double)
-   * value, and outputs a TecPlot/AMDiS/VTK mesh.
+   * value, and outputs a AMDiS/VTK mesh.
    */
   class ElementFileWriter : public FileWriterInterface
   {
@@ -34,17 +34,15 @@ namespace AMDiS {
     /// Simple writing procedure for one element map.
     static void writeFile(std::map<int, double> &vec,
 			  Mesh *mesh,
-			  std::string filename);
+			  std::string filename,
+			  int level = -1);
 
   protected:
-    /// Writes element data in tecplot format.
-    void writeTecPlotValues(std::string filename);
-
     /// Writes element data in AMDiS format (1 file !).
     void writeMeshDatValues(std::string filename, double time);
 
     /// Writes element data in VTK format.
-    void writeVtkValues(std::string filename);
+    void writeVtkValues(std::string filename, int level = -1);
 
   protected:
     /// Name.
@@ -53,37 +51,22 @@ namespace AMDiS {
     /// Used filename prefix.
     std::string filename;
 
-    /// TecPlot file extension.
-    std::string tecplotExt;
-
     /// AMDiS mesh-data-file extension.
     std::string amdisMeshDatExt;
 
     /// VTK file extension.
     std::string vtkExt;
 
-    /** \brief
-     * 0: Don't write TecPlot files.
-     * 1: Write TecPlot files. 
-     */
-    int writeTecPlotFormat;
-
-    /** \brief
-     * 0: Don't write AMDiS files.
-     * 1: Write AMDiS files. 
-     */
+    /// 0: Don't write AMDiS files.
+    /// 1: Write AMDiS files. 
     int writeAMDiSFormat;
 
-    /** \brief
-     * 0: Don't write VTK files.
-     * 1: Write VTK files. 
-     */
+    /// 0: Don't write VTK files.
+    /// 1: Write VTK files. 
     int writeVtkFormat;
 
-    /** \brief
-     * 0: Don't append time index to filename prefix.
-     * 1: Append time index to filename prefix.
-     */
+    /// 0: Don't append time index to filename prefix.
+    /// 1: Append time index to filename prefix.
     int appendIndex;
 
     /// Total length of appended time index.
diff --git a/AMDiS/src/parallel/ElementObjectData.cc b/AMDiS/src/parallel/ElementObjectData.cc
index e1b2a9e6c4fdb5d6a346cf3b4bff6316f6948303..a6a5d89b505400113c7cea4ea851127eb9e32b96 100644
--- a/AMDiS/src/parallel/ElementObjectData.cc
+++ b/AMDiS/src/parallel/ElementObjectData.cc
@@ -84,12 +84,14 @@ namespace AMDiS {
     }
 
 
+
     SerUtil::serialize(out, vertexOwner);
     SerUtil::serialize(out, edgeOwner);
     SerUtil::serialize(out, faceOwner);
 
 
     nSize = vertexInRank.size();
+    SerUtil::serialize(out, nSize);
     for (std::map<DegreeOfFreedom, std::map<int, ElementObjectData> >::iterator it = vertexInRank.begin();
 	 it != vertexInRank.end(); ++it) {
       SerUtil::serialize(out, it->first);
@@ -97,6 +99,7 @@ namespace AMDiS {
     }
 
     nSize = edgeInRank.size();
+    SerUtil::serialize(out, nSize);
     for (std::map<DofEdge, std::map<int, ElementObjectData> >::iterator it = edgeInRank.begin();
 	 it != edgeInRank.end(); ++it) {
       SerUtil::serialize(out, it->first);
@@ -104,6 +107,7 @@ namespace AMDiS {
     }
 
     nSize = faceInRank.size();
+    SerUtil::serialize(out, nSize);
     for (std::map<DofFace, std::map<int, ElementObjectData> >::iterator it = faceInRank.begin();
 	 it != faceInRank.end(); ++it) {
       SerUtil::serialize(out, it->first);
@@ -147,7 +151,7 @@ namespace AMDiS {
       deserialize(in, data);
       faceElements[face] = data;
     }
-    
+   
 
     SerUtil::deserialize(in, vertexOwner);
     SerUtil::deserialize(in, edgeOwner);
diff --git a/AMDiS/src/parallel/InteriorBoundary.cc b/AMDiS/src/parallel/InteriorBoundary.cc
index abf05375fb23e4da6741a542d2260899a98c7314..61bc95fac073f6eb5d3b4eae15bd71556dcfe4cc 100644
--- a/AMDiS/src/parallel/InteriorBoundary.cc
+++ b/AMDiS/src/parallel/InteriorBoundary.cc
@@ -128,6 +128,8 @@ namespace AMDiS {
 
   void InteriorBoundary::serialize(std::ostream &out)
   {
+    FUNCNAME("InteriorBoundary::serialize()");
+
     int mSize = boundary.size();
     SerUtil::serialize(out, mSize);
     for (RankToBoundMap::iterator it = boundary.begin(); it != boundary.end(); ++it) {
diff --git a/AMDiS/src/parallel/MeshDistributor.cc b/AMDiS/src/parallel/MeshDistributor.cc
index f45cfe9e262d2713791c4ce51e2d5425342bbcbf..eeeda478f1060184ca91e37eab7b885e25d1db33 100644
--- a/AMDiS/src/parallel/MeshDistributor.cc
+++ b/AMDiS/src/parallel/MeshDistributor.cc
@@ -41,6 +41,11 @@ namespace AMDiS {
     return (*dof1 < *dof2);
   }
 
+  bool cmpElement(MacroElement* lhs, MacroElement* rhs)
+  {
+    return lhs->getIndex() < rhs->getIndex();
+  }
+  
 
   MeshDistributor::MeshDistributor(std::string str)
     : probStat(0),
@@ -86,18 +91,27 @@ namespace AMDiS {
     // isRankDofs to all matrices and rhs vector and to remove periodic 
     // boundary conditions (if there are some).
     if (deserialized) {
+      updateMacroElementInfo();
+
       setRankDofs();
+
       removePeriodicBoundaryConditions();
 
       macroElIndexMap.clear();
       macroElIndexTypeMap.clear();
-      std::deque<MacroElement*>& macros = mesh->getMacroElements();
-      for (std::deque<MacroElement*>::iterator it = macros.begin();
-	   it != macros.end(); ++it) {
+
+      std::map<int, bool>& elementInRank = partitioner->getElementInRank();
+      for (std::vector<MacroElement*>::iterator it = allMacroElements.begin();
+	   it != allMacroElements.end(); ++it) {
+	elementInRank[(*it)->getIndex()] = false;
 	macroElIndexMap[(*it)->getIndex()] = (*it)->getElement();
 	macroElIndexTypeMap[(*it)->getIndex()] = (*it)->getElType();
       }
 
+      for (std::deque<MacroElement*>::iterator it = mesh->getMacroElements().begin();
+	   it != mesh->getMacroElements().end(); ++it)
+	elementInRank[(*it)->getIndex()] = true;      
+
       return;
     }
 
@@ -561,10 +575,10 @@ namespace AMDiS {
     nTimestepsAfterLastRepartitioning++;
 
     if (repartitioningAllowed) {
-      if (nTimestepsAfterLastRepartitioning >= 20) {
+      //      if (nTimestepsAfterLastRepartitioning >= 20) {
 	repartitionMesh();
 	nTimestepsAfterLastRepartitioning = 0;
-      }
+	//      }
     }
   }
 
@@ -853,7 +867,7 @@ namespace AMDiS {
     int vecSize = data.size();
     SerUtil::serialize(out, vecSize);
     for (int i = 0; i < vecSize; i++) {
-      int dofIndex = (*data[i]);
+      int dofIndex = *(data[i]);
       SerUtil::serialize(out, dofIndex);
     }
   }
@@ -980,7 +994,9 @@ namespace AMDiS {
     if (repartitioning == 0)
       return;
 
-    if (repartCounter == 0) {
+#if (DEBUG != 0)
+    //    if (repartCounter == 0) {
+ {
       std::stringstream oss;
       oss << "partitioning-" << repartCounter << ".vtu";
       DOFVector<double> tmpa(feSpace, "tmp");
@@ -988,9 +1004,11 @@ namespace AMDiS {
       VtkWriter::writeFile(&tmpa, oss.str());
       
       repartCounter++;
-    }
+ }
 
     MSG("USED-SIZE A: %d\n", mesh->getDofAdmin(0).getUsedDofs());
+#endif
+
 
     elemWeights.clear();
 
@@ -1001,20 +1019,18 @@ namespace AMDiS {
       elInfo = stack.traverseNext(elInfo);
     }
 
+
     partitioner->useLocalGlobalDofMap(&mapLocalGlobalDofs);
     oldPartitionVec = partitionVec;
     partitioner->partition(&elemWeights, ADAPTIVE_REPART, 1000.0);    
 
-    MacroInfo* minfo = mesh->getMacroFileInfo();
-    TEST_EXIT_DBG(minfo)("No macro file info!\n");
-
 
     // === Create map that maps macro element indices to pointers to the ===
     // === macro elements.                                               ===
 
     std::map<int, MacroElement*> elIndexMap;
-    for (std::deque<MacroElement*>::iterator it = minfo->mel.begin();
-	 it != minfo->mel.end(); ++it)
+    for (std::vector<MacroElement*>::iterator it = allMacroElements.begin();
+	 it != allMacroElements.end(); ++it)
       elIndexMap[(*it)->getIndex()] = *it;
 
 
@@ -1028,29 +1044,45 @@ namespace AMDiS {
 	   elIt != it->second.end(); ++elIt) {
 	TEST_EXIT_DBG(elIndexMap.count(*elIt) == 1)
 	  ("Could not find macro element %d\n", *elIt);
+
 	newMacroEl.insert(elIndexMap[*elIt]);
       }
     }
 
-    std::set<MacroElement*> allMacroEl = newMacroEl;
+    bool (*fn_pt)(MacroElement*, MacroElement*) = cmpElement;
+    std::set<MacroElement* , bool(*)(MacroElement*, MacroElement*)> allMacroEl(fn_pt);
+
+    for (std::set<MacroElement*>::iterator it = newMacroEl.begin(); 
+	 it != newMacroEl.end(); ++it)
+      allMacroEl.insert(*it);
+
     for (std::deque<MacroElement*>::iterator it = mesh->firstMacroElement();
 	 it != mesh->endOfMacroElements(); ++it)
       allMacroEl.insert(*it);
 
+
     // === Add new macro elements to mesh. ===
 
     for (std::set<MacroElement*>::iterator it = newMacroEl.begin();
 	 it != newMacroEl.end(); ++it) {
       MacroElement *mel = *it;
+
+      // First, reset all neighbour relations. The correct neighbours will be set later.
       for (int i = 0; i < mesh->getGeo(NEIGH); i++)
 	mel->setNeighbour(i, NULL);
 
+      // Create new DOFs for the macro element.
       for (int i = 0; i < mesh->getGeo(VERTEX); i++)
 	mel->getElement()->setDof(i, mesh->getDof(VERTEX));
 
+      // The macro element now has valid DOF pointers.
+      mel->getElement()->setDofValid(true);
+
+      // Push the macro element to all macro elements in mesh.
       mesh->getMacroElements().push_back(mel);
     }
 
+
     // === Send and receive mesh structure codes. ===
 
     std::map<int, MeshCodeVec> sendCodes;
@@ -1097,6 +1129,7 @@ namespace AMDiS {
       for (std::vector<int>::iterator elIt = it->second.begin();
 	   elIt != it->second.end(); ++elIt) {
 	recvCodes[i].fitMeshToStructure(mesh, refineManager, false, *elIt);
+
 	for (unsigned int k = 0; k < testVec.size(); k++)
 	  recvCodes[i].setMeshStructureValues(mesh, *elIt, testVec[k], recvValues[j++]);    
 	
@@ -1107,16 +1140,27 @@ namespace AMDiS {
 
     // === Set correct neighbour information on macro elements. ===
 
-    for (std::set<MacroElement*>::iterator it = allMacroEl.begin();
+    for (std::set<MacroElement* , bool(*)(MacroElement*, MacroElement*)>::iterator it = allMacroEl.begin();
 	 it != allMacroEl.end(); ++it) {
       int elIndex = (*it)->getIndex();
       for (int i = 0; i < mesh->getGeo(NEIGH); i++) {
+
+	TEST_EXIT_DBG(macroElementNeighbours.count(elIndex) == 1)
+	  ("Should not happen!\n");
+	TEST_EXIT_DBG(static_cast<int>(macroElementNeighbours[elIndex].size()) == 
+		      mesh->getGeo(NEIGH))
+	  ("Should not happen!\n");
+
 	int neighIndex = macroElementNeighbours[elIndex][i];
+
 	if (neighIndex == -1 || 
-	    partitioner->getElementInRank()[neighIndex] == false)
+	    partitioner->getElementInRank()[neighIndex] == false) {
 	  (*it)->setNeighbour(i, NULL);
-	else
+	} else {
+	  TEST_EXIT_DBG(elIndexMap.count(neighIndex) == 1)("Should not happen!\n");
+
 	  (*it)->setNeighbour(i, elIndexMap[neighIndex]);
+	}
       }
     }
 
@@ -1136,18 +1180,20 @@ namespace AMDiS {
       }
     }
 
-    mesh->removeMacroElements(deleteMacroElements, feSpace);
 
-    if (mpiRank == 0)
-      debug::writeElementIndexMesh(mesh, "elementIndexxxx.vtu");
+    // Note that also if there are no macros to be deleted, this function will
+    // update the number of elements, vertices, etc. of the mesh.
+    mesh->removeMacroElements(deleteMacroElements, feSpace);
 
 
     // === Remove double DOFs. ===
+
     MeshManipulation meshManipulation(mesh);
     meshManipulation.deleteDoubleDofs(newMacroEl);
 
     mesh->dofCompress();
 
+#if (DEBUG != 0)
     std::stringstream oss;
     oss << "partitioning-" << repartCounter << ".vtu";    
     DOFVector<double> tmpa(feSpace, "tmp");
@@ -1157,7 +1203,6 @@ namespace AMDiS {
 
     MSG("USED-SIZE B: %d\n", mesh->getDofAdmin(0).getUsedDofs());
 
-#if (DEBUG != 0)
     ParallelDebug::testAllElements(*this);
 #endif
 
@@ -1179,6 +1224,8 @@ namespace AMDiS {
 
     MSG("Debug mode tests finished!\n");
 #endif
+
+    ParallelDebug::printBoundaryInfo(*this);
   }
 
 
@@ -1346,7 +1393,9 @@ namespace AMDiS {
 	if (objData.count(mpiRank) && objData.size() > 1) {
 	  int owner = elObjects.getIterateOwner();
 	  ElementObjectData& rankBoundEl = objData[mpiRank];
-	  
+
+	  TEST_EXIT_DBG(macroElIndexMap[rankBoundEl.elIndex])("Should not happen!\n");
+
 	  AtomicBoundary bound;	    	    
 	  bound.rankObj.el = macroElIndexMap[rankBoundEl.elIndex];
 	  bound.rankObj.elIndex = rankBoundEl.elIndex;
@@ -1646,9 +1695,6 @@ namespace AMDiS {
     // === Traverse on interior boundaries and move all not ranked owned DOFs from ===
     // === rankDofs to boundaryDofs.                                               ===
 
-    RankToDofContainer oldSendDofs = sendDofs;
-    RankToDofContainer oldRecvDofs = recvDofs;
-
     sendDofs.clear();
     recvDofs.clear();
 
@@ -1889,6 +1935,8 @@ namespace AMDiS {
 
     for (std::deque<MacroElement*>::iterator it = mesh->firstMacroElement();
 	 it != mesh->endOfMacroElements(); ++it) {
+      allMacroElements.push_back(*it);
+
       macroElementNeighbours[(*it)->getIndex()].resize(mesh->getGeo(NEIGH));
       for (int i = 0; i < mesh->getGeo(NEIGH); i++)
 	macroElementNeighbours[(*it)->getIndex()][i] = 
@@ -1897,8 +1945,29 @@ namespace AMDiS {
   }
 
 
+  void MeshDistributor::updateMacroElementInfo()
+  {
+    std::deque<MacroElement*>& meshMacros = mesh->getMacroElements();
+
+    for (unsigned int i = 0; i < allMacroElements.size(); i++) {
+      for (std::deque<MacroElement*>::iterator it = meshMacros.begin();
+	   it != meshMacros.end(); ++it) {
+	if ((*it)->getIndex() == allMacroElements[i]->getIndex() &&
+	    *it != allMacroElements[i]) {
+	  delete allMacroElements[i];
+	  allMacroElements[i] = *it;
+	}
+      }
+    }
+  }
+
+
   void MeshDistributor::serialize(std::ostream &out)
   {
+    FUNCNAME("MeshDistributor::serialize()");
+
+    checkMeshChange();
+
     SerUtil::serialize(out, elemWeights);
     SerUtil::serialize(out, partitionVec);
     SerUtil::serialize(out, oldPartitionVec);
@@ -1927,6 +1996,14 @@ namespace AMDiS {
 
     SerUtil::serialize(out, rstart);
     SerUtil::serialize(out, macroElementNeighbours);
+
+    int nSize = allMacroElements.size();
+    SerUtil::serialize(out, nSize);
+    for (int i = 0; i < nSize; i++)
+      allMacroElements[i]->serialize(out);
+
+    SerUtil::serialize(out, nTimestepsAfterLastRepartitioning);
+    SerUtil::serialize(out, repartCounter);
   }
 
 
@@ -1984,6 +2061,21 @@ namespace AMDiS {
     SerUtil::deserialize(in, rstart);
     SerUtil::deserialize(in, macroElementNeighbours);
 
+    int nSize = 0;
+    SerUtil::deserialize(in, nSize);
+    allMacroElements.resize(nSize);
+    for (int i = 0; i < nSize; i++) {
+      // We do not need the neighbour indices, but must still read them.
+      std::vector<int> indices;
+      allMacroElements[i] = new MacroElement(mesh->getDim());
+      allMacroElements[i]->writeNeighboursTo(&indices);
+      allMacroElements[i]->deserialize(in);
+      allMacroElements[i]->getElement()->setMesh(mesh);
+    }
+
+    SerUtil::deserialize(in, nTimestepsAfterLastRepartitioning);
+    SerUtil::deserialize(in, repartCounter);
+
     deserialized = true;
   }
 
diff --git a/AMDiS/src/parallel/MeshDistributor.h b/AMDiS/src/parallel/MeshDistributor.h
index 3020cde5595cdafb2456546801b33ff305af4361..6b9dff24e08ad8e980932e92aabd9d73cfbad4b2 100644
--- a/AMDiS/src/parallel/MeshDistributor.h
+++ b/AMDiS/src/parallel/MeshDistributor.h
@@ -259,6 +259,8 @@ namespace AMDiS {
 
     void createMacroElementInfo();
 
+    void updateMacroElementInfo();
+
     /** \brief
      * Checks for all given interior boundaries if the elements fit together on both
      * sides of the boundaries. If this is not the case, the mesh is adapted. Because
@@ -449,14 +451,14 @@ namespace AMDiS {
     std::map<int, double> elemWeights;
 
     /** \brief
-     * Stores to every coarse element index the number of the partition it 
-     * corresponds to.
+     * Stores to every macro element index the number of the rank that owns this
+     * macro element.
      */
     std::map<int, int> partitionVec;
 
     /** \brief
-     * Stores an old partitioning of elements. To every element index the number
-     * of the parition it corresponds to is stored.
+     * Stores an old partitioning of elements. To every macro element index the
+     * number of the rank it corresponds to is stored.
      */
     std::map<int, int> oldPartitionVec;    
    
@@ -577,6 +579,10 @@ namespace AMDiS {
 
     std::map<int, std::vector<int> > macroElementNeighbours;
 
+    /// Store all macro elements of the overall mesh, i.e., before the macro mesh is
+    /// redistributed for the first time.
+    std::vector<MacroElement*> allMacroElements;
+
     friend class ParallelDebug;
   };
 }
diff --git a/AMDiS/src/parallel/MeshManipulation.cc b/AMDiS/src/parallel/MeshManipulation.cc
index e7516cc74201ee2113ca60405c911d8182128aad..3ac1be257c5102b2ec5afdbe8e04c96b34a522d2 100644
--- a/AMDiS/src/parallel/MeshManipulation.cc
+++ b/AMDiS/src/parallel/MeshManipulation.cc
@@ -1,6 +1,7 @@
 #include "parallel/MeshManipulation.h"
 #include "Mesh.h"
 #include "Traverse.h"
+#include "Debug.h"
 
 namespace AMDiS {
 
@@ -15,11 +16,8 @@ namespace AMDiS {
       elInfo = stack.traverseNext(elInfo);
     }
 
-//     MSG("----------------------------------\n");
     deleteDoubleDofs(leafInMacroEl, newMacroEl, 0);
-//     MSG("----------------------------------\n");
     deleteDoubleDofs(leafInMacroEl, newMacroEl, 1);
-//     MSG("----------------------------------\n");
   }
 
 
@@ -38,10 +36,17 @@ namespace AMDiS {
       Element *el = elInfo->getElement();
 
       for (int i = 0; i < mesh->getGeo(NEIGH); i++) {
+
 	Element *neigh = elInfo->getNeighbour(i);
 
 	if (!neigh)
-	  continue;      
+	  continue;
+
+	TEST_EXIT_DBG(leafInMacroEl.count(el->getIndex()) == 1)
+	  ("Should not happen!\n");
+	TEST_EXIT_DBG(leafInMacroEl.count(neigh->getIndex()) == 1)
+	  ("No macro element for el %d, that is %d-ith neigbour of element %d!\n",
+	   neigh->getIndex(), i, el->getIndex());
 
 	int elMacroIndex = leafInMacroEl[el->getIndex()]->getIndex();
 	int neighMacroIndex = leafInMacroEl[neigh->getIndex()]->getIndex();
@@ -57,8 +62,6 @@ namespace AMDiS {
 	     newMacroEl.count(leafInMacroEl[el->getIndex()]) == 0 &&
 	     newMacroEl.count(leafInMacroEl[neigh->getIndex()]) == 1)) {
 
-// 	  MSG("EL: %d   NEIGH:  %d\n", el->getIndex(), neigh->getIndex());
-
 	  if (el->getEdge(i) != neigh->getEdge(elInfo->getOppVertex(i))) {
 	    std::vector<int> elEdgeIndex(2);
 	    elEdgeIndex[0] = el->getVertexOfEdge(i, 0);
@@ -76,15 +79,13 @@ namespace AMDiS {
 	    neighEdgeDof[0] = neigh->getDof(neighEdgeIndex[0], 0);
 	    neighEdgeDof[1] = neigh->getDof(neighEdgeIndex[1], 0);
 
-// 	    MSG("TEST DOFs %d  %d    <->    %d  %d\n", 
-// 		elEdgeDof[0], elEdgeDof[1], neighEdgeDof[0], neighEdgeDof[1]);
-
 	    if ((elEdgeDof[0] < elEdgeDof[1] && neighEdgeDof[0] > neighEdgeDof[1]) ||
 		(elEdgeDof[0] > elEdgeDof[1] && neighEdgeDof[0] < neighEdgeDof[1])) {
 	      std::swap(neighEdgeIndex[0], neighEdgeIndex[1]);
 	      std::swap(neighEdgeDof[0], neighEdgeDof[1]);
 	    }
 
+
 	    for (int j = 0; j < 2; j++) {
 
 	      if (neighEdgeDof[j] != elEdgeDof[j]) {
@@ -104,14 +105,16 @@ namespace AMDiS {
 		  mapDofsMacro[replaceDof] = elMacroIndex;
 		}
 	      }
-	    }	      
+	    }	     
 	  }
 	}
       }
 
+
       elInfo = stack.traverseNext(elInfo);
     }
 
+
     bool changed = false;
     do {
       changed = false;
@@ -129,24 +132,15 @@ namespace AMDiS {
     elInfo = stack.traverseFirst(mesh, -1, Mesh::CALL_EVERY_EL_PREORDER);
     while (elInfo) {
       for (int i = 0; i < mesh->getGeo(VERTEX); i++)
-	if (mapDelDofs.count(elInfo->getElement()->getDof(i)) == 1) {
-// 	  MSG("SET DOF %d TO DOF %d in EL %d\n",
-// 	      elInfo->getElement()->getDof(i, 0),
-// 	      *(mapDelDofs[elInfo->getElement()->getDof(i)]),
-// 	      elInfo->getElement()->getIndex());
-
+	if (mapDelDofs.count(elInfo->getElement()->getDof(i)) == 1)
 	  elInfo->getElement()->setDof(i, const_cast<DegreeOfFreedom*>(mapDelDofs[elInfo->getElement()->getDof(i)]));
-	}      
 
       elInfo = stack.traverseNext(elInfo);
     }
 
     for (std::map<const DegreeOfFreedom*, const DegreeOfFreedom*>::iterator it = mapDelDofs.begin();
-	 it != mapDelDofs.end(); ++it) {
-      //      MSG("DOF TO DEL: %d\n", (*(it->first)));
-
+	 it != mapDelDofs.end(); ++it)
       mesh->freeDof(const_cast<DegreeOfFreedom*>(it->first), VERTEX);
-    }
   }
 
 }
diff --git a/AMDiS/src/parallel/ParMetisPartitioner.cc b/AMDiS/src/parallel/ParMetisPartitioner.cc
index 30332bbda8ce9ba0cfbfd2d863e5e45e6013919a..183e6760c87bb8c6d8d7d242edc0656d163edaec 100644
--- a/AMDiS/src/parallel/ParMetisPartitioner.cc
+++ b/AMDiS/src/parallel/ParMetisPartitioner.cc
@@ -207,7 +207,8 @@ namespace AMDiS {
     if (parMetisMesh) 
       delete parMetisMesh;
 
-    TEST_EXIT(elementInRank.size() != 0)("Should not happen!\n");
+    TEST_EXIT_DBG(elementInRank.size() != 0)("Should not happen!\n");
+
     parMetisMesh = new ParMetisMesh(mesh, mpiComm, elementInRank, mapLocalGlobal);
 
     int nElements = parMetisMesh->getNumElements();
diff --git a/AMDiS/src/parallel/ParMetisPartitioner.h b/AMDiS/src/parallel/ParMetisPartitioner.h
index 9730560f9761fc98b3a8c10a664bd8f10bda9563..12ad24c85805aee98065aba2565eb34f10091f38 100644
--- a/AMDiS/src/parallel/ParMetisPartitioner.h
+++ b/AMDiS/src/parallel/ParMetisPartitioner.h
@@ -186,10 +186,7 @@ namespace AMDiS {
      */
     void fillCoarsePartitionVec(std::map<int, int> *partitionVec);
 
-    /* \brief
-     * Creates an initial paritioning of the AMDiS mesh by seting the partition status
-     * of all elements to either IN or UNDEFINED.
-     */
+    /// Creates an initial paritioning of the AMDiS mesh.
     void createPartitionData();
 
     void useLocalGlobalDofMap(std::map<DegreeOfFreedom, DegreeOfFreedom> *m)