diff --git a/CMakeLists.txt b/CMakeLists.txt
index a7ecb6035f2e753232dbf758e22330a9e6a46bed..0bb7f48e75ec8392cb92c398eff926df57ba2007 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -4,7 +4,7 @@ project("dumux" C CXX)
 # general stuff
 cmake_minimum_required(VERSION 2.8.6)
 
-if(NOT (dune-common_DIR 
+if(NOT (dune-common_DIR
         OR dune-common_ROOT
         OR "${CMAKE_PREFIX_PATH}" MATCHES ".*dune-common.*"))
     string(REPLACE  ${CMAKE_PROJECT_NAME}
diff --git a/dumux/decoupled/2p/diffusion/fv/fvpressure2p.hh b/dumux/decoupled/2p/diffusion/fv/fvpressure2p.hh
index fd52260b54c45bc1c240b1cf7f9f6844b6d4c37f..3599de3b69aafc637fa1c78a99d5f8a46dc7dcd4 100644
--- a/dumux/decoupled/2p/diffusion/fv/fvpressure2p.hh
+++ b/dumux/decoupled/2p/diffusion/fv/fvpressure2p.hh
@@ -443,8 +443,8 @@ public:
                 (*mobilityNW)[i]  = cellData.mobility(nPhaseIdx);
                 if (compressibility_)
                 {
-                	(*mobilityW)[i] = (*mobilityW)[i]/cellData.density(wPhaseIdx);
-                	(*mobilityNW)[i] = (*mobilityNW)[i]/cellData.density(nPhaseIdx);
+                    (*mobilityW)[i] = (*mobilityW)[i]/cellData.density(wPhaseIdx);
+                    (*mobilityNW)[i] = (*mobilityNW)[i]/cellData.density(nPhaseIdx);
                 }
             }
             if (vtkOutputLevel_ > 1)
diff --git a/dumux/freeflow/stokes/stokeslocalresidual.hh b/dumux/freeflow/stokes/stokeslocalresidual.hh
index e47482ac97997762445a16663c66740a03c74299..2e3043f28cfaa19864f18b138bc1a8b493ed3cd1 100644
--- a/dumux/freeflow/stokes/stokeslocalresidual.hh
+++ b/dumux/freeflow/stokes/stokeslocalresidual.hh
@@ -617,7 +617,7 @@ protected:
 
     /*!
      * \brief Interpolate the pressure at corner points of the grid, thus taking the degree of freedom there.
-     * 		  This is required due to stability reasons.
+     *        This is required due to stability reasons.
      */
     void interpolateCornerPoints_(const BoundaryTypes &bcTypes, const int scvIdx)
     {
diff --git a/dumux/freeflow/stokesnc/stokesnclocalresidual.hh b/dumux/freeflow/stokesnc/stokesnclocalresidual.hh
index a766f08311da74da71f37609e8e381fa7a7e7434..0c1255b81d496d9a4717673a2792b15bbc748210 100644
--- a/dumux/freeflow/stokesnc/stokesnclocalresidual.hh
+++ b/dumux/freeflow/stokesnc/stokesnclocalresidual.hh
@@ -49,7 +49,7 @@ class StokesncLocalResidual : public StokesLocalResidual<TypeTag>
     typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
 
     //dimensions
-	enum {  dim = GridView::dimension };
+    enum {  dim = GridView::dimension };
     //number of equations
     enum {  numEq = GET_PROP_VALUE(TypeTag, NumEq) };
     //number of components
@@ -63,25 +63,25 @@ class StokesncLocalResidual : public StokesLocalResidual<TypeTag>
     //primary variable indices
     enum {  pressureIdx = Indices::pressureIdx };
     //phase employed
-    enum { 	phaseIdx = Indices::phaseIdx };
+    enum {  phaseIdx = Indices::phaseIdx };
     //component indices
     enum {  phaseCompIdx = Indices::phaseCompIdx,
             transportCompIdx = Indices::transportCompIdx };
 
-	typedef typename GET_PROP_TYPE(TypeTag, PrimaryVariables) PrimaryVariables;
+    typedef typename GET_PROP_TYPE(TypeTag, PrimaryVariables) PrimaryVariables;
     typedef typename GET_PROP_TYPE(TypeTag, VolumeVariables) VolumeVariables;
     typedef typename GET_PROP_TYPE(TypeTag, FluxVariables) FluxVariables;
     typedef typename GET_PROP_TYPE(TypeTag, ElementVolumeVariables) ElementVolumeVariables;
     typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
 
-	typedef Dune::FieldVector<Scalar, dim> DimVector;
+    typedef Dune::FieldVector<Scalar, dim> DimVector;
 
-	typedef typename GridView::Intersection Intersection;
+    typedef typename GridView::Intersection Intersection;
     typedef typename GET_PROP_TYPE(TypeTag, BoundaryTypes) BoundaryTypes;
 
-	static const bool calculateNavierStokes = GET_PROP_VALUE(TypeTag, EnableNavierStokes);
+    static const bool calculateNavierStokes = GET_PROP_VALUE(TypeTag, EnableNavierStokes);
 
-	//! property that defines whether mole or mass fractions are used
+    //! property that defines whether mole or mass fractions are used
     static const bool useMoles = GET_PROP_VALUE(TypeTag, UseMoles);
 
 public:
@@ -108,7 +108,7 @@ public:
         // is required to compute the derivative of the storage term
         // using the implicit Euler method.
         const ElementVolumeVariables &elemVolVars = usePrevSol ?
-        			this->prevVolVars_() : this->curVolVars_();
+                    this->prevVolVars_() : this->curVolVars_();
         const VolumeVariables &volVars = elemVolVars[scvIdx];
 
         if (useMoles)
@@ -155,12 +155,12 @@ public:
     void computeAdvectiveFlux(PrimaryVariables &flux,
                               const FluxVariables &fluxVars) const
     {
-      	// call ParentType function
-		ParentType::computeAdvectiveFlux(flux,fluxVars);
+        // call ParentType function
+        ParentType::computeAdvectiveFlux(flux,fluxVars);
 
         // data attached to upstream and the downstream vertices
-		const VolumeVariables &up = this->curVolVars_(fluxVars.upstreamIdx());
-		const VolumeVariables &dn = this->curVolVars_(fluxVars.downstreamIdx());
+        const VolumeVariables &up = this->curVolVars_(fluxVars.upstreamIdx());
+        const VolumeVariables &dn = this->curVolVars_(fluxVars.downstreamIdx());
 
         Scalar tmp = fluxVars.normalVelocity();
 
diff --git a/dumux/freeflow/stokesnc/stokesncmodel.hh b/dumux/freeflow/stokesnc/stokesncmodel.hh
index 225beb7e3ee4b06acba2a0ed4cc7cd3c9b3fdc03..876b19603d1f060bcfeab5310b23083716451f45 100644
--- a/dumux/freeflow/stokesnc/stokesncmodel.hh
+++ b/dumux/freeflow/stokesnc/stokesncmodel.hh
@@ -81,11 +81,11 @@ class StokesncModel : public StokesModel<TypeTag>
     typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
 
     enum {  dim = GridView::dimension,
-			transportCompIdx = Indices::transportCompIdx,
-			phaseIdx = GET_PROP_VALUE(TypeTag, PhaseIdx),
-			useMoles = GET_PROP_VALUE(TypeTag, UseMoles),
-			numComponents = Indices::numComponents
-	};
+            transportCompIdx = Indices::transportCompIdx,
+            phaseIdx = GET_PROP_VALUE(TypeTag, PhaseIdx),
+            useMoles = GET_PROP_VALUE(TypeTag, UseMoles),
+            numComponents = Indices::numComponents
+    };
 
     typedef typename GET_PROP_TYPE(TypeTag, FVElementGeometry) FVElementGeometry;
     typedef typename GET_PROP_TYPE(TypeTag, ElementBoundaryTypes) ElementBoundaryTypes;
@@ -102,7 +102,7 @@ public:
                             MultiWriter &writer)
     {
 
-		typedef Dune::BlockVector<Dune::FieldVector<Scalar, 1> > ScalarField;
+        typedef Dune::BlockVector<Dune::FieldVector<Scalar, 1> > ScalarField;
         typedef Dune::BlockVector<Dune::FieldVector<Scalar, dim> > VelocityField;
 
         const Scalar scale_ = GET_PROP_VALUE(TypeTag, Scaling);
@@ -111,17 +111,17 @@ public:
         unsigned numVertices = this->gridView_().size(dim);
         ScalarField &pN = *writer.allocateManagedBuffer(numVertices);
         ScalarField &delP = *writer.allocateManagedBuffer(numVertices);
-		ScalarField &T = *writer.allocateManagedBuffer(numVertices);
+        ScalarField &T = *writer.allocateManagedBuffer(numVertices);
 
-		ScalarField *moleFraction[numComponents];
+        ScalarField *moleFraction[numComponents];
             for (int i = 0; i < numComponents; ++i)
                 moleFraction[i] = writer.template allocateManagedBuffer<Scalar, 1>(numVertices);
 
-		ScalarField *massFraction[numComponents];
-		for (int i = 0; i < numComponents; ++i)
-			massFraction[i] = writer.template allocateManagedBuffer<Scalar, 1>(numVertices);
+        ScalarField *massFraction[numComponents];
+        for (int i = 0; i < numComponents; ++i)
+            massFraction[i] = writer.template allocateManagedBuffer<Scalar, 1>(numVertices);
 
-		ScalarField &rho = *writer.allocateManagedBuffer(numVertices);
+        ScalarField &rho = *writer.allocateManagedBuffer(numVertices);
         ScalarField &mu = *writer.allocateManagedBuffer(numVertices);
         VelocityField &velocity = *writer.template allocateManagedBuffer<Scalar, dim> (numVertices);
 
@@ -156,23 +156,23 @@ public:
 
                 pN[vIdxGlobal] = volVars.pressure()*scale_;
                 delP[vIdxGlobal] = volVars.pressure()*scale_ - 1e5;
-				for (int compIdx = 0; compIdx < numComponents; ++compIdx)
+                for (int compIdx = 0; compIdx < numComponents; ++compIdx)
                     {
                         (*moleFraction[compIdx])[vIdxGlobal]= volVars.moleFraction(compIdx);
-						(*massFraction[compIdx])[vIdxGlobal]= volVars.massFraction(compIdx);
+                        (*massFraction[compIdx])[vIdxGlobal]= volVars.massFraction(compIdx);
                         Valgrind::CheckDefined((*moleFraction[compIdx])[vIdxGlobal]);
-						Valgrind::CheckDefined((*massFraction[compIdx])[vIdxGlobal]);
-					}
+                        Valgrind::CheckDefined((*massFraction[compIdx])[vIdxGlobal]);
+                    }
 
-				T   [vIdxGlobal] = volVars.temperature();
+                T   [vIdxGlobal] = volVars.temperature();
 
-				rho[vIdxGlobal] = volVars.density()*scale_*scale_*scale_;
+                rho[vIdxGlobal] = volVars.density()*scale_*scale_*scale_;
                 mu[vIdxGlobal] = volVars.dynamicViscosity()*scale_;
                 velocity[vIdxGlobal] = volVars.velocity();
                 velocity[vIdxGlobal] *= 1/scale_;
             }
         }
-		writer.attachVertexData(T, "temperature");
+        writer.attachVertexData(T, "temperature");
         writer.attachVertexData(pN, "P");
         writer.attachVertexData(delP, "delP");
 
diff --git a/dumux/freeflow/stokesnc/stokesncvolumevariables.hh b/dumux/freeflow/stokesnc/stokesncvolumevariables.hh
index 8617d33798f80ebad9d75a00bc6ea4272ff10286..e4238bc6a5086b7e5c1811138539c231b3939ca8 100644
--- a/dumux/freeflow/stokesnc/stokesncvolumevariables.hh
+++ b/dumux/freeflow/stokesnc/stokesncvolumevariables.hh
@@ -57,12 +57,12 @@ class StokesncVolumeVariables : public StokesVolumeVariables<TypeTag>
     enum {  transportCompIdx = Indices::transportCompIdx,
             phaseCompIdx = Indices::phaseCompIdx };
     //number of components
-	enum {	numComponents = Indices::numComponents };
+    enum {  numComponents = Indices::numComponents };
     //employed phase index
-	enum {	phaseIdx = GET_PROP_VALUE(TypeTag, PhaseIdx) };
+    enum {  phaseIdx = GET_PROP_VALUE(TypeTag, PhaseIdx) };
     //primary variable indices
-	enum {	massOrMoleFracIdx = Indices::massOrMoleFracIdx };
-	//equation indices
+    enum {  massOrMoleFracIdx = Indices::massOrMoleFracIdx };
+    //equation indices
     enum {  conti0EqIdx = Indices::conti0EqIdx,
             massBalanceIdx = Indices::massBalanceIdx,
             transportEqIdx = Indices::transportEqIdx };
@@ -82,17 +82,17 @@ public:
                 const bool isOldSol)
     {
 
-		// Model is restricted to 2 components when using mass fractions
-		if (!useMoles && numComponents>2)
-		{
-			DUNE_THROW(Dune::NotImplemented, "This model is restricted to 2 components when using mass fractions!\
-			                                  To use mole fractions set property UseMoles true ...");
-		}
+        // Model is restricted to 2 components when using mass fractions
+        if (!useMoles && numComponents>2)
+        {
+            DUNE_THROW(Dune::NotImplemented, "This model is restricted to 2 components when using mass fractions!\
+                                              To use mole fractions set property UseMoles true ...");
+        }
 
-		// set the mole fractions first
+        // set the mole fractions first
         completeFluidState(priVars, problem, element, fvGeometry, scvIdx, this->fluidState(), isOldSol);
 
-		// update vertex data for the mass and momentum balance
+        // update vertex data for the mass and momentum balance
         ParentType::update(priVars,
                            problem,
                            element,
@@ -108,19 +108,19 @@ public:
 
         for (int compIdx=0; compIdx<numComponents; compIdx++)
         {
-			if (phaseCompIdx!=compIdx)
-			{
-				diffCoeff_[compIdx] = FluidSystem::binaryDiffusionCoefficient(this->fluidState(),
+            if (phaseCompIdx!=compIdx)
+            {
+                diffCoeff_[compIdx] = FluidSystem::binaryDiffusionCoefficient(this->fluidState(),
                                                              paramCache,
                                                              phaseIdx,
                                                              compIdx,
                                                              phaseCompIdx);
-			}
-			else
-				diffCoeff_[compIdx] = 0.0;
+            }
+            else
+                diffCoeff_[compIdx] = 0.0;
 
-			Valgrind::CheckDefined(diffCoeff_[compIdx]);
-		}
+            Valgrind::CheckDefined(diffCoeff_[compIdx]);
+        }
     };
 
     /*!
@@ -163,7 +163,7 @@ public:
 
     /*!
      * \brief Returns the mass fraction of a given component in the
-     * 		  given fluid phase within the control volume.
+     *        given fluid phase within the control volume.
      *
      * \param compIdx The component index
      */
diff --git a/dumux/freeflow/stokesncni/stokesncnimodel.hh b/dumux/freeflow/stokesncni/stokesncnimodel.hh
index 991349d92b04863a0c658751ecf27e1a5101565d..ab319cab9b03a3c2aade2bbc6257c8269986bc44 100644
--- a/dumux/freeflow/stokesncni/stokesncnimodel.hh
+++ b/dumux/freeflow/stokesncni/stokesncnimodel.hh
@@ -94,7 +94,7 @@ class StokesncniModel : public StokesncModel<TypeTag>
     enum { transportCompIdx = Indices::transportCompIdx };
     enum { phaseIdx = GET_PROP_VALUE(TypeTag, PhaseIdx) };
     enum { useMoles = GET_PROP_VALUE(TypeTag, UseMoles) };
-	enum { numComponents = Indices::numComponents };
+    enum { numComponents = Indices::numComponents };
 
    typedef typename GET_PROP_TYPE(TypeTag, FVElementGeometry) FVElementGeometry;
     typedef typename GET_PROP_TYPE(TypeTag, ElementBoundaryTypes) ElementBoundaryTypes;
@@ -118,17 +118,17 @@ public:
         unsigned numVertices = this->gridView_().size(dim);
         ScalarField &pn = *writer.allocateManagedBuffer(numVertices);
         ScalarField &delP = *writer.allocateManagedBuffer(numVertices);
-		ScalarField &T = *writer.allocateManagedBuffer(numVertices);
+        ScalarField &T = *writer.allocateManagedBuffer(numVertices);
         ScalarField &h = *writer.allocateManagedBuffer(numVertices);
-		ScalarField *moleFraction[numComponents];
-		for (int i = 0; i < numComponents; ++i)
+        ScalarField *moleFraction[numComponents];
+        for (int i = 0; i < numComponents; ++i)
         moleFraction[i] = writer.template allocateManagedBuffer<Scalar, 1>(numVertices);
 
-		ScalarField *massFraction[numComponents];
-		for (int i = 0; i < numComponents; ++i)
+        ScalarField *massFraction[numComponents];
+        for (int i = 0; i < numComponents; ++i)
         massFraction[i] = writer.template allocateManagedBuffer<Scalar, 1>(numVertices);
 
-		ScalarField &rho = *writer.allocateManagedBuffer(numVertices);
+        ScalarField &rho = *writer.allocateManagedBuffer(numVertices);
         ScalarField &mu = *writer.allocateManagedBuffer(numVertices);
         VelocityField &velocity = *writer.template allocateManagedBuffer<Scalar, dim> (numVertices);
 
@@ -163,24 +163,24 @@ public:
 
                 pn[vIdxGlobal] = volVars.pressure()*scale_;
                 delP[vIdxGlobal] = volVars.pressure()*scale_ - 1e5;
-				for (int compIdx = 0; compIdx < numComponents; ++compIdx)
-				{
-					(*moleFraction[compIdx])[vIdxGlobal]= volVars.moleFraction(compIdx);
-					(*massFraction[compIdx])[vIdxGlobal]= volVars.massFraction(compIdx);
-					Valgrind::CheckDefined((*moleFraction[compIdx])[vIdxGlobal]);
-					Valgrind::CheckDefined((*massFraction[compIdx])[vIdxGlobal]);
-				}
+                for (int compIdx = 0; compIdx < numComponents; ++compIdx)
+                {
+                    (*moleFraction[compIdx])[vIdxGlobal]= volVars.moleFraction(compIdx);
+                    (*massFraction[compIdx])[vIdxGlobal]= volVars.massFraction(compIdx);
+                    Valgrind::CheckDefined((*moleFraction[compIdx])[vIdxGlobal]);
+                    Valgrind::CheckDefined((*massFraction[compIdx])[vIdxGlobal]);
+                }
 
-				T   [vIdxGlobal] = volVars.temperature();
+                T   [vIdxGlobal] = volVars.temperature();
 
-				rho[vIdxGlobal] = volVars.density()*scale_*scale_*scale_;
+                rho[vIdxGlobal] = volVars.density()*scale_*scale_*scale_;
                 mu[vIdxGlobal] = volVars.dynamicViscosity()*scale_;
                 h[vIdxGlobal] = volVars.enthalpy();
                 velocity[vIdxGlobal] = volVars.velocity();
                 velocity[vIdxGlobal] *= 1/scale_;
             }
         }
-		writer.attachVertexData(T, "temperature");
+        writer.attachVertexData(T, "temperature");
         writer.attachVertexData(pn, "pg");
         writer.attachVertexData(delP, "delP");
 
diff --git a/dumux/implicit/2p/2pgridadaptindicator.hh b/dumux/implicit/2p/2pgridadaptindicator.hh
index 9af37602b3f5000e0e64b3fcd91201aaa7f2bf9e..6065177ddea5faf70857a8b6e594417058b3a500 100644
--- a/dumux/implicit/2p/2pgridadaptindicator.hh
+++ b/dumux/implicit/2p/2pgridadaptindicator.hh
@@ -47,8 +47,8 @@ private:
 
     enum
     {
-    	saturationIdx = Indices::saturationIdx,
-		pressureIdx = Indices::pressureIdx
+        saturationIdx = Indices::saturationIdx,
+        pressureIdx = Indices::pressureIdx
     };
     enum
     {
@@ -96,17 +96,17 @@ public:
             // index of the current leaf-elements
             int globalIdxI = problem_.elementMapper().index(element);
 
-        	Scalar satI = 0.0;
+            Scalar satI = 0.0;
 
-        	if(!isBox)
-        		satI = problem_.model().curSol()[globalIdxI][saturationIdx];
-        	else
-        	{
+            if(!isBox)
+                satI = problem_.model().curSol()[globalIdxI][saturationIdx];
+            else
+            {
                 const LocalFiniteElementCache feCache;
                 const auto geometryI = element.geometry();
-            	Dune::GeometryType geomType = geometryI.type();
+                Dune::GeometryType geomType = geometryI.type();
 
-            	GlobalPosition centerI = geometryI.local(geometryI.center());
+                GlobalPosition centerI = geometryI.local(geometryI.center());
                 const LocalFiniteElement &localFiniteElement = feCache.get(geomType);
                 std::vector<Dune::FieldVector<Scalar, 1> > shapeVal;
                 localFiniteElement.localBasis().evaluateFunction(centerI, shapeVal);
@@ -116,7 +116,7 @@ public:
                      int dofIdxGlobal = problem_.model().dofMapper().subIndex(element, i, dofCodim);
                       satI += shapeVal[i]*problem_.model().curSol()[dofIdxGlobal][saturationIdx];
                   }
-        	}
+            }
 
             globalMin = std::min(satI, globalMin);
             globalMax = std::max(satI, globalMax);
@@ -134,17 +134,17 @@ public:
                     // Visit intersection only once
                     if (element.level() > outside.level() || (element.level() == outside.level() && globalIdxI < globalIdxJ))
                     {
-                    	Scalar satJ = 0.0;
+                        Scalar satJ = 0.0;
 
-                    	if(!isBox)
-                    		satJ = problem_.model().curSol()[globalIdxJ][saturationIdx];
-                    	else
-                    	{
+                        if(!isBox)
+                            satJ = problem_.model().curSol()[globalIdxJ][saturationIdx];
+                        else
+                        {
                             const LocalFiniteElementCache feCache;
                             const auto geometryJ = outside.geometry();
-                        	Dune::GeometryType geomType = geometryJ.type();
+                            Dune::GeometryType geomType = geometryJ.type();
 
-                        	GlobalPosition centerJ = geometryJ.local(geometryJ.center());
+                            GlobalPosition centerJ = geometryJ.local(geometryJ.center());
                             const LocalFiniteElement &localFiniteElement = feCache.get(geomType);
                             std::vector<Dune::FieldVector<Scalar, 1> > shapeVal;
                             localFiniteElement.localBasis().evaluateFunction(centerJ, shapeVal);
@@ -155,7 +155,7 @@ public:
 
                                   satJ += shapeVal[i]*problem_.model().curSol()[dofIdxGlobal][saturationIdx];
                               }
-                    	}
+                        }
 
 
 
diff --git a/dumux/implicit/2pminc/2pmincindices.hh b/dumux/implicit/2pminc/2pmincindices.hh
index 50b4145bb99bce08d2b947b74f1521fce1bcda3c..1ac2e22904c1b2da4264e5b2ee997df8165ac1ea 100644
--- a/dumux/implicit/2pminc/2pmincindices.hh
+++ b/dumux/implicit/2pminc/2pmincindices.hh
@@ -69,10 +69,10 @@ struct TwoPMincCommonIndices
  * \tparam formulation The formulation, either pwsn or pnsw
  * \tparam PVOffset The first index in a primary variable vector.
  */
-template <class TypeTag, 
-          int formulation = TwoPMincFormulation::pwsn, 
+template <class TypeTag,
+          int formulation = TwoPMincFormulation::pwsn,
           int PVOffset = 0 >
-struct TwoPMincIndices 
+struct TwoPMincIndices
 : public TwoPMincCommonIndices<TypeTag>, TwoPMincFormulation
 {
     // Primary variable indices
@@ -82,14 +82,14 @@ struct TwoPMincIndices
     // indices of the primary variables
     static const int pwIdx = PVOffset + 0; //!< index of the wetting phase pressure
     static const int snIdx = PVOffset + 1; //!< index of the nonwetting phase saturation
-    
+
     // indices of the equations
     static const int contiWEqIdx = PVOffset + 0; //!< Index of the continuity equation of the wetting phase
     static const int contiNEqIdx = PVOffset + 1; //!< Index of the continuity equation of the non-wetting phase
-    
+
     static const int pIdxc(int numC) {return pwIdx + 2 *numC;} //!< index of the wetting phase pressure for continuum numC
     static const int sIdxc(int numC) {return snIdx + 2 *numC;} //!< index of the non-wetting phase saturation for continuum numC
-    
+
     static const int contiWEqIdxc(int numC) {return contiWEqIdx + 2 *numC;} //!< Index of the continuity equation of the wetting phase for continuum numC
     static const int contiNEqIdxc(int numC) {return contiNEqIdx + 2 *numC;} //!< Index of the continuity equation of the non-wetting phase for continuum numC
 };
@@ -120,7 +120,7 @@ struct TwoPMincIndices<TypeTag, TwoPMincFormulation::pnsw, PVOffset>
 
     static const int pIdxc(int numC) {return pnIdx + 2 *numC;} //!< index of the nonwetting phase pressure for continuum numC
     static const int sIdxc(int numC) {return swIdx + 2 *numC;} //!< index of the wetting phase saturation for continuum numC
-     
+
     static const int contiWEqIdxc(int numC) {return contiWEqIdx + 2 *numC;} //!< Index of the continuity equation of the wetting phase for continuum numC
     static const int contiNEqIdxc(int numC) {return contiNEqIdx + 2 *numC;} //!< Index of the continuity equation of the non-wetting phase for continuum numC
 };
diff --git a/dumux/implicit/2pminc/2pmincmodel.hh b/dumux/implicit/2pminc/2pmincmodel.hh
index 5e2948685b5e8b3ad7ce68152b02db65e3da921d..924bcb9f142369d4c754603daf95494e71853c01 100644
--- a/dumux/implicit/2pminc/2pmincmodel.hh
+++ b/dumux/implicit/2pminc/2pmincmodel.hh
@@ -253,7 +253,7 @@ public:
             interfaceArea_[nC] = 0.0;
 //            transmissibility_[nC] = 0.0;
             distNestedContinua_[nC] =0.0;
-//        	volumetricFraction_[nC] =0.0;
+//          volumetricFraction_[nC] =0.0;
         }
 
         /*
diff --git a/dumux/implicit/2pminc/2pmincvolumevariables.hh b/dumux/implicit/2pminc/2pmincvolumevariables.hh
index 08be6ba3e5cfa0320a3ef2a98262b3d063f4dbd2..78e2e0cf83c43e82ca9326bf2e536897d1e8e751 100644
--- a/dumux/implicit/2pminc/2pmincvolumevariables.hh
+++ b/dumux/implicit/2pminc/2pmincvolumevariables.hh
@@ -97,7 +97,7 @@ public:
         const MaterialLawParams &materialParamsMatrix =
             problem.spatialParams().materialLawParamsMatrix(element, fvGeometry, scvIdx);
 
-		// relative permeabilities krw/krn for fractures (idx 0) and matrix elements (>= idx 1)
+        // relative permeabilities krw/krn for fractures (idx 0) and matrix elements (>= idx 1)
         for (int cIdx = 0; cIdx < numContinua; ++cIdx)
         {
             Scalar krw;
diff --git a/dumux/implicit/2pnc/2pncfluxvariables.hh b/dumux/implicit/2pnc/2pncfluxvariables.hh
index 165dc8857d7537c259d4416e693764efe256e86e..4a31de496322a6ddb0a7caa3a4f2e6a499e73d80 100644
--- a/dumux/implicit/2pnc/2pncfluxvariables.hh
+++ b/dumux/implicit/2pnc/2pncfluxvariables.hh
@@ -46,7 +46,7 @@ namespace Dumux
 template <class TypeTag>
 class TwoPNCFluxVariables : public GET_PROP_TYPE(TypeTag, BaseFluxVariables)
 {
-	typedef typename GET_PROP_TYPE(TypeTag, BaseFluxVariables) BaseFluxVariables;
+    typedef typename GET_PROP_TYPE(TypeTag, BaseFluxVariables) BaseFluxVariables;
     typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
     typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
 
@@ -105,8 +105,8 @@ public:
             potentialGrad_[phaseIdx] = Scalar(0);
             for (int compIdx = 0; compIdx < numComponents; ++compIdx)
             {
-            	massFractionGrad_[phaseIdx][compIdx] = Scalar(0);
-            	moleFractionGrad_[phaseIdx][compIdx] = Scalar(0);
+                massFractionGrad_[phaseIdx][compIdx] = Scalar(0);
+                moleFractionGrad_[phaseIdx][compIdx] = Scalar(0);
             }
         }
         calculateGradients_(problem, element, elemVolVars);
@@ -269,7 +269,7 @@ protected:
             Scalar tauI =  1.0/(volVarsI.porosity() * volVarsI.porosity()) *
                             pow(volVarsI.porosity() * volVarsI.saturation(phaseIdx), 7.0/3);
 
-            Scalar tauJ =	1.0/(volVarsJ.porosity() * volVarsJ.porosity()) *
+            Scalar tauJ =   1.0/(volVarsJ.porosity() * volVarsJ.porosity()) *
                             pow(volVarsJ.porosity() * volVarsJ.saturation(phaseIdx), 7.0/3);
             // Diffusion coefficient in the porous medium
 
@@ -300,7 +300,7 @@ public:
      * mobility is not yet included here since this would require a
      * decision on the upwinding approach (which is done in the
      * model and/or local residual file).
-     * 
+     *
      *   \param phaseIdx The phase index
      */
     Scalar KmvpNormal(int phaseIdx) const
@@ -309,7 +309,7 @@ public:
     /*!
      * \brief Return the pressure potential multiplied with the
      *        intrinsic permeability as vector (for velocity output)
-     * 
+     *
      *   \param phaseIdx The phase index
      */
     DimVector Kmvp(int phaseIdx) const
@@ -318,7 +318,7 @@ public:
     /*!
      * \brief Return the local index of the upstream control volume
      *        for a given phase.
-     * 
+     *
      *   \param phaseIdx The phase index
      */
     int upstreamIdx(int phaseIdx) const
@@ -327,7 +327,7 @@ public:
     /*!
      * \brief Return the local index of the downstream control volume
      *        for a given phase.
-     * 
+     *
      *   \param phaseIdx The phase index
      */
     int downstreamIdx(int phaseIdx) const
@@ -335,7 +335,7 @@ public:
 
     /*!
      * \brief The binary diffusion coefficient for each fluid phase.
-     * 
+     *
      *   \param phaseIdx The phase index
      *   \param compIdx The component index
      */
@@ -345,7 +345,7 @@ public:
     /*!
      * \brief Return density \f$\mathrm{[kg/m^3]}\f$ of a phase at the integration
      *        point.
-     * 
+     *
      * \param phaseIdx The phase index
      */
     Scalar density(int phaseIdx) const
@@ -354,7 +354,7 @@ public:
     /*!
      * \brief Return molar density \f$\mathrm{[mol/m^3]}\f$ of a phase at the integration
      *        point.
-     * 
+     *
      * \param phaseIdx The phase index
      */
     Scalar molarDensity(int phaseIdx) const
@@ -362,7 +362,7 @@ public:
 
     /*!
      * \brief The concentration gradient of a component in a phase.
-     * 
+     *
      * \param phaseIdx The phase index
      * \param compIdx The component index
      */
@@ -371,7 +371,7 @@ public:
 
     /*!
      * \brief The molar concentration gradient of a component in a phase.
-     * 
+     *
      * \param phaseIdx The phase index
      * \param compIdx The component index
      */
diff --git a/dumux/implicit/2pnc/2pncindices.hh b/dumux/implicit/2pnc/2pncindices.hh
index 671917a14488bd82a5f8d39c79cdd28ba09a5364..16e55dad8d5b2feb21914842b0fc2a9478090a6c 100644
--- a/dumux/implicit/2pnc/2pncindices.hh
+++ b/dumux/implicit/2pnc/2pncindices.hh
@@ -53,10 +53,10 @@ struct TwoPNCFormulation//TODO: This might need to be change similar to 2p2c ind
 template <class TypeTag, int PVOffset = 0>
 class TwoPNCIndices
 {
-	typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
+    typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
 
 public:
-	// Phase indices
+    // Phase indices
     static const int wPhaseIdx = FluidSystem::wPhaseIdx; //!< Index of the wetting phase
     static const int nPhaseIdx = FluidSystem::nPhaseIdx; //!< Index of the non-wetting phase
     // present phases (-> 'pseudo' primary variable)
diff --git a/dumux/implicit/2pnc/2pnclocalresidual.hh b/dumux/implicit/2pnc/2pnclocalresidual.hh
index 731f19090cd86cc79fa2b3099b60d74c6148c4a4..2191d9e4408b044cd70e95c140ebef45826d2cc1 100644
--- a/dumux/implicit/2pnc/2pnclocalresidual.hh
+++ b/dumux/implicit/2pnc/2pnclocalresidual.hh
@@ -241,7 +241,7 @@ public:
          // add advective flux of current component in current
          // phase
             unsigned int eqIdx = conti0EqIdx + compIdx;
-            
+
          if (eqIdx != replaceCompEqIdx)
          {
             // upstream vertex
@@ -327,7 +327,7 @@ public:
                                         this->fvGeometry_(),
                                         scvIdx,
                                             this->curVolVars_());
-        
+
         Valgrind::CheckDefined(source);
     }
 
@@ -337,7 +337,7 @@ protected:
     void evalPhaseStorage_(int phaseIdx)
     {
         // evaluate the storage terms of a single phase
-        for (int i=0; i < this->fvGeometry_().numScv; i++) 
+        for (int i=0; i < this->fvGeometry_().numScv; i++)
         {
             PrimaryVariables &result = this->residual_[i];
             const ElementVolumeVariables &elemVolVars = this->curVolVars_();
diff --git a/dumux/implicit/2pnc/2pncpropertydefaults.hh b/dumux/implicit/2pnc/2pncpropertydefaults.hh
index b2a1c4209194763fbc5aa74af3454a80c3f4e9f3..3d4212c160f0e0e53d88e2a64a270093a1f9308f 100644
--- a/dumux/implicit/2pnc/2pncpropertydefaults.hh
+++ b/dumux/implicit/2pnc/2pncpropertydefaults.hh
@@ -73,7 +73,7 @@ private:
 public:
     static const int value = FluidSystem::numComponents;
 };
-//! The major components belonging to the existing phases are mentioned here e.g., 2 for water and air being the major component in the liquid and gas phases in a 2 phase system 
+//! The major components belonging to the existing phases are mentioned here e.g., 2 for water and air being the major component in the liquid and gas phases in a 2 phase system
 SET_PROP(TwoPNC, NumMajorComponents)
 {
 private:
diff --git a/dumux/implicit/2pnc/2pncvolumevariables.hh b/dumux/implicit/2pnc/2pncvolumevariables.hh
index 6a6bb92d5271405c3e5daf2c1895f7e541f58987..27b7d97c861161af0fdbd0dfb282385baf7343d9 100644
--- a/dumux/implicit/2pnc/2pncvolumevariables.hh
+++ b/dumux/implicit/2pnc/2pncvolumevariables.hh
@@ -264,7 +264,7 @@ public:
 
             Miscible2pNCComposition::solve(fluidState,
                         paramCache,
-                        wPhaseIdx,	//known phaseIdx
+                        wPhaseIdx,  //known phaseIdx
                         /*setViscosity=*/true,
                         /*setInternalEnergy=*/false);
         }
diff --git a/dumux/implicit/2pncmin/2pncminfluxvariables.hh b/dumux/implicit/2pncmin/2pncminfluxvariables.hh
index c48b67bf90c9b541d43f7b41f5040a1cb0b3ffc6..3e0f7940f7aa9cc60e97a6e4437944369146cc95 100644
--- a/dumux/implicit/2pncmin/2pncminfluxvariables.hh
+++ b/dumux/implicit/2pncmin/2pncminfluxvariables.hh
@@ -49,7 +49,7 @@ class TwoPNCMinFluxVariables : public TwoPNCFluxVariables<TypeTag>
 {
     typedef TwoPNCFluxVariables<TypeTag> ParentType;
     typedef TwoPNCMinFluxVariables<TypeTag> ThisType;
-    
+
     typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
     typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
 
@@ -100,7 +100,7 @@ public:
                      const int fIdx,
                      const ElementVolumeVariables &elemVolVars,
                      const bool onBoundary = false)
-    : ParentType(problem, element, fvGeometry, fIdx, elemVolVars, onBoundary) 
+    : ParentType(problem, element, fvGeometry, fIdx, elemVolVars, onBoundary)
     {
         for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
             this->density_[phaseIdx] = Scalar(0);
@@ -108,8 +108,8 @@ public:
             this->potentialGrad_[phaseIdx] = Scalar(0);
             for (int compIdx = 0; compIdx < numComponents; ++compIdx)
             {
-            	this->massFractionGrad_[phaseIdx][compIdx] = Scalar(0);
-            	this->moleFractionGrad_[phaseIdx][compIdx] = Scalar(0);
+                this->massFractionGrad_[phaseIdx][compIdx] = Scalar(0);
+                this->moleFractionGrad_[phaseIdx][compIdx] = Scalar(0);
             }
         }
         this->calculateGradients_(problem, element, elemVolVars);
@@ -117,7 +117,7 @@ public:
         this->calculateporousDiffCoeff_(problem, element, elemVolVars);
     };
 
-protected:    
+protected:
     void calculateVelocities_(const Problem &problem,
                               const Element &element,
                               const ElementVolumeVariables &elemVolVars)
diff --git a/dumux/implicit/2pncmin/2pncminindices.hh b/dumux/implicit/2pncmin/2pncminindices.hh
index dc1b2e9e0cedcf6fcac673409271bb32bc407229..8f4e1036667bb45ae69aa543c38b1642a4138d12 100644
--- a/dumux/implicit/2pncmin/2pncminindices.hh
+++ b/dumux/implicit/2pncmin/2pncminindices.hh
@@ -19,7 +19,7 @@
 
 /*!
  * \file
- * \brief Defines the indices required for the two-phase n-component mineralization 
+ * \brief Defines the indices required for the two-phase n-component mineralization
  *        fully implicit model.
  */
 #ifndef DUMUX_2PNCMIN_INDICES_HH
@@ -37,8 +37,8 @@ namespace Dumux
  * \tparam PVOffset The first index in a primary variable vector.
  */
 template <class TypeTag, int PVOffset = 0>
-	class TwoPNCMinIndices: public TwoPNCIndices<TypeTag, PVOffset>
-{ 
+    class TwoPNCMinIndices: public TwoPNCIndices<TypeTag, PVOffset>
+{
 };
 
 // \}
diff --git a/dumux/implicit/2pncmin/2pncminlocalresidual.hh b/dumux/implicit/2pncmin/2pncminlocalresidual.hh
index c0fc7f9529a46d5a231d9c8f6ce6ba08994bc075..72299478f459cc52f34eb863fdc08e3b991f748b 100644
--- a/dumux/implicit/2pncmin/2pncminlocalresidual.hh
+++ b/dumux/implicit/2pncmin/2pncminlocalresidual.hh
@@ -43,9 +43,9 @@ template<class TypeTag>
 class TwoPNCMinLocalResidual: public TwoPNCLocalResidual<TypeTag>
 {
 protected:
-	typedef TwoPNCLocalResidual<TypeTag> ParentType;
-    typedef TwoPNCMinLocalResidual<TypeTag> ThisType; 
-	typedef typename GET_PROP_TYPE(TypeTag, Problem) Problem;
+    typedef TwoPNCLocalResidual<TypeTag> ParentType;
+    typedef TwoPNCMinLocalResidual<TypeTag> ThisType;
+    typedef typename GET_PROP_TYPE(TypeTag, Problem) Problem;
     typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
     typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
     typedef typename GET_PROP_TYPE(TypeTag, SolutionVector) SolutionVector;
@@ -120,7 +120,7 @@ public:
   {
       //call parenttype function
       ParentType::computeStorage(storage, scvIdx, usePrevSol);
-      
+
       const ElementVolumeVariables &elemVolVars = usePrevSol ? this->prevVolVars_()
       : this->curVolVars_();
     const VolumeVariables &volVars = elemVolVars[scvIdx];
diff --git a/dumux/implicit/2pncmin/2pncminpropertydefaults.hh b/dumux/implicit/2pncmin/2pncminpropertydefaults.hh
index ad4412120d2e2c7b3afb6364d03546720e98020a..fa29ac3f81627793d536daddfff28a5da671b5c3 100644
--- a/dumux/implicit/2pncmin/2pncminpropertydefaults.hh
+++ b/dumux/implicit/2pncmin/2pncminpropertydefaults.hh
@@ -49,8 +49,8 @@ namespace Properties {
 
 /*!
  * \brief Set the property for the number of secondary components.
- * Secondary components are components calculated from 
- * primary components by equilibrium relations and 
+ * Secondary components are components calculated from
+ * primary components by equilibrium relations and
  * do not have mass balance equation on their own.
  * These components are important in the context of bio-mineralization applications.
  * We just forward the number from the fluid system
@@ -81,7 +81,7 @@ public:
 };
 
 /*!
- * \brief Set the property for the number of equations. 
+ * \brief Set the property for the number of equations.
  * For each component and each precipitated mineral/solid phase one equation has to
  * be solved.
  */
diff --git a/dumux/implicit/2pncmin/2pncminvolumevariables.hh b/dumux/implicit/2pncmin/2pncminvolumevariables.hh
index 1b7391af03a4a0b7001105446ce53abb855a0a6b..f8e3b7a060eba1125c7236efbcc705abc22772f9 100644
--- a/dumux/implicit/2pncmin/2pncminvolumevariables.hh
+++ b/dumux/implicit/2pncmin/2pncminvolumevariables.hh
@@ -63,7 +63,7 @@ class TwoPNCMinVolumeVariables : public TwoPNCVolumeVariables<TypeTag>
     typedef typename GET_PROP_TYPE(TypeTag, MaterialLawParams) MaterialLawParams;
     typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
 
-    enum 
+    enum
     {
         dim = GridView::dimension,
         dimWorld=GridView::dimensionworld,
@@ -107,7 +107,7 @@ class TwoPNCMinVolumeVariables : public TwoPNCVolumeVariables<TypeTag>
     enum { isBox = GET_PROP_VALUE(TypeTag, ImplicitIsBox) };
     enum { dofCodim = isBox ? dim : 0 };
 public:
-  
+
       //! The type of the object returned by the fluidState() method
       typedef CompositionalFluidState<Scalar, FluidSystem> FluidState;
     /*!
@@ -144,10 +144,10 @@ public:
        precipitateVolumeFraction_[sPhaseIdx] = priVars[numComponents + sPhaseIdx];
        sumPrecipitates_+= precipitateVolumeFraction_[sPhaseIdx];
     }
-    
+
 //         for(int sPhaseIdx = 0; sPhaseIdx < numSPhases; ++sPhaseIdx)
 //     {
-//         Chemistry chemistry; // the non static functions can not be called without abject 
+//         Chemistry chemistry; // the non static functions can not be called without abject
 //         saturationIdx_[sPhaseIdx] = chemistry.omega(sPhaseIdx);
 //     }
 // TODO/FIXME: The salt crust porosity is not clearly defined. However form literature review it is
@@ -162,7 +162,7 @@ public:
 
    salinity_= 0.0;
    moleFractionSalinity_ = 0.0;
-   for (int compIdx = numMajorComponents; compIdx< numComponents; compIdx++)	//sum of the mass fraction of the components
+   for (int compIdx = numMajorComponents; compIdx< numComponents; compIdx++)    //sum of the mass fraction of the components
    {
        if(this->fluidState_.moleFraction(wPhaseIdx, compIdx)> 0)
        {
@@ -283,12 +283,12 @@ public:
             // can be used by the Miscible2pNcComposition constraint solver
             for (int compIdx=numMajorComponents; compIdx<numComponents; ++compIdx)
             {
-            	fluidState.setMoleFraction(wPhaseIdx, compIdx, priVars[compIdx]);
+                fluidState.setMoleFraction(wPhaseIdx, compIdx, priVars[compIdx]);
             }
 
             Miscible2pNCComposition::solve(fluidState,
                                             paramCache,
-                                            wPhaseIdx,	//known phaseIdx
+                                            wPhaseIdx,  //known phaseIdx
                                             /*setViscosity=*/true,
                                             /*setInternalEnergy=*/false);
         }
@@ -312,12 +312,12 @@ public:
             for (int compIdx=numMajorComponents; compIdx<numComponents; ++compIdx)
                 moleFrac[compIdx] = (priVars[compIdx]*fugCoeffL[compIdx]*fluidState.pressure(wPhaseIdx))
                 /(fugCoeffG[compIdx]*fluidState.pressure(nPhaseIdx));
-            
+
             moleFrac[wCompIdx] =  priVars[switchIdx];
             Scalar sumMoleFracNotGas = 0;
             for (int compIdx=numMajorComponents; compIdx<numComponents; ++compIdx)
             {
-            		sumMoleFracNotGas+=moleFrac[compIdx];
+                    sumMoleFracNotGas+=moleFrac[compIdx];
             }
             sumMoleFracNotGas += moleFrac[wCompIdx];
             moleFrac[nCompIdx] = 1 - sumMoleFracNotGas;
@@ -329,7 +329,7 @@ public:
             // Set fluid state mole fractions
             for (int compIdx=0; compIdx<numComponents; ++compIdx)
             {
-            	fluidState.setMoleFraction(nPhaseIdx, compIdx, moleFrac[compIdx]);
+                fluidState.setMoleFraction(nPhaseIdx, compIdx, moleFrac[compIdx]);
             }
 
             // calculate the composition of the remaining phases (as
@@ -344,7 +344,7 @@ public:
 
             }
         else if (phasePresence == wPhaseOnly){
-        
+
         // only the liquid phase is present, i.e. liquid phase
         // composition is stored explicitly.
         // extract _mass_ fractions in the gas phase
@@ -352,13 +352,13 @@ public:
 
             for (int compIdx=numMajorComponents; compIdx<numComponents; ++compIdx)
             {
-            	moleFrac[compIdx] = priVars[compIdx];
+                moleFrac[compIdx] = priVars[compIdx];
             }
             moleFrac[nCompIdx] = priVars[switchIdx];
             Scalar sumMoleFracNotWater = 0;
             for (int compIdx=numMajorComponents; compIdx<numComponents; ++compIdx)
             {
-            		sumMoleFracNotWater+=moleFrac[compIdx];
+                    sumMoleFracNotWater+=moleFrac[compIdx];
             }
             sumMoleFracNotWater += moleFrac[nCompIdx];
             moleFrac[wCompIdx] = 1 -sumMoleFracNotWater;
@@ -387,54 +387,54 @@ public:
         }
     }
     /*!
-     * \brief Returns the volume fraction of the precipitate (solid phase) 
+     * \brief Returns the volume fraction of the precipitate (solid phase)
      * for the given phaseIdx
-     * 
+     *
      * \param phaseIdx the index of the solid phase
      */
     Scalar precipitateVolumeFraction(int phaseIdx) const
     { return precipitateVolumeFraction_[phaseIdx - numPhases]; }
-   
+
     /*!
-     * \brief Returns the inital porosity of the 
+     * \brief Returns the inital porosity of the
      * pure, precipitate-free porous medium
      */
     Scalar initialPorosity() const
     { return initialPorosity_;}
-    
+
     /*!
-     * \brief Returns the inital permeability of the 
+     * \brief Returns the inital permeability of the
      * pure, precipitate-free porous medium
-     */    
+     */
     Scalar initialPermeability() const
     { return initialPermeability_;}
-    
+
     /*!
-     * \brief Returns the factor for the reduction of the initial permeability 
+     * \brief Returns the factor for the reduction of the initial permeability
      * due precipitates in the porous medium
-     */ 	
+     */
     Scalar permeabilityFactor() const
     { return permeabilityFactor_; }
-    
+
     /*!
      * \brief Returns the mole fraction of a component in the phase
-     * 
+     *
      * \param phaseIdx the index of the fluid phase
      * \param compIdx the index of the component
-     */ 
+     */
     Scalar moleFraction(int phaseIdx, int compIdx) const
     {
        return this->fluidState_.moleFraction(phaseIdx, compIdx);
     }
-    
+
     /*!
      * \brief Returns the mole fraction of the salinity in the liquid phase
      */
     Scalar moleFracSalinity() const
     {
-    	return moleFractionSalinity_;
+        return moleFractionSalinity_;
     }
-    
+
     /*!
      * \brief Returns the salinity (mass fraction) in the liquid phase
      */
@@ -442,12 +442,12 @@ public:
     {
         return salinity_;
     }
-    
+
     /*!
-     * \brief Returns the density of the phase for all fluid and solid phases 
-     * 
+     * \brief Returns the density of the phase for all fluid and solid phases
+     *
      * \param phaseIdx the index of the fluid phase
-     */ 
+     */
     Scalar density(int phaseIdx) const
     {
         if (phaseIdx < numPhases)
@@ -472,13 +472,13 @@ public:
         else
             DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx);
     }
-    
+
     /*!
      * \brief Returns the molality of a component in the phase
-     * 
+     *
      * \param phaseIdx the index of the fluid phase
      * \param compIdx the index of the component
-     */ 
+     */
      Scalar molality(int phaseIdx, int compIdx) const // [moles/Kg]
     { return this->fluidState_.moleFraction(phaseIdx, compIdx)/FluidSystem::molarMass(compIdx);}
 
@@ -492,7 +492,7 @@ protected:
     {
         return problem.temperatureAtPos(fvGeometry.subContVol[scvIdx].global);
     }
-    
+
     template<class ParameterCache>
     static Scalar enthalpy_(const FluidState& fluidState,
                             const ParameterCache& paramCache,
diff --git a/dumux/implicit/3p3c/3p3cfluxvariables.hh b/dumux/implicit/3p3c/3p3cfluxvariables.hh
index 80f3f46b272e696d27e92f97d0b9e99a31c0b52b..f3c3d62d8d236af9b8f72fe1616b9efbbcc7e21e 100644
--- a/dumux/implicit/3p3c/3p3cfluxvariables.hh
+++ b/dumux/implicit/3p3c/3p3cfluxvariables.hh
@@ -66,7 +66,7 @@ class ThreePThreeCFluxVariables : public GET_PROP_TYPE(TypeTag, BaseFluxVariable
 
     typedef typename GET_PROP_TYPE(TypeTag, FVElementGeometry) FVElementGeometry;
 
-    typedef Dune::FieldVector<Scalar, dim> 	DimVector;
+    typedef Dune::FieldVector<Scalar, dim>  DimVector;
     typedef Dune::FieldVector<Scalar, dimWorld> GlobalPosition;
 
     typedef typename GET_PROP_TYPE(TypeTag, Indices) Indices;
diff --git a/dumux/implicit/adaptive/adaptionhelper.hh b/dumux/implicit/adaptive/adaptionhelper.hh
index f15f75eea7df5b38718e0cd21b119bbafc2b487d..9fb9be811cb741a537709c0093e7760642ef8178 100644
--- a/dumux/implicit/adaptive/adaptionhelper.hh
+++ b/dumux/implicit/adaptive/adaptionhelper.hh
@@ -54,7 +54,7 @@ private:
 
     struct AdaptedValues
     {
-    	PrimaryVariables u;
+        PrimaryVariables u;
         int count;
         AdaptedValues()
         {
@@ -98,7 +98,7 @@ public:
      *  @param gridView a DUNE gridview object corresponding to diffusion and transport equation
      */
     AdaptionHelper(const GridView& gridView) :
-    	gridView_(gridView), grid_(gridView.grid()), adaptionMap_(grid_, dofCodim)
+        gridView_(gridView), grid_(gridView.grid()), adaptionMap_(grid_, dofCodim)
     {}
 
 
@@ -125,46 +125,46 @@ public:
 
             if(!isBox)
             {
-				for (const auto& element : Dune::elements(levelView))
-				{
-					//get your map entry
-					AdaptedValues &adaptedValues = adaptionMap_[element];
-
-					// put your value in the map
-					if (element.isLeaf())
-					{
-						// get index
-						int indexI = this->elementIndex(problem, element);
-
-						storeAdaptionValues(adaptedValues, problem.model().curSol()[indexI]);
-
-						adaptedValues.count = 1;
-					}
-					//Average in father
-					if (element.level() > 0)
-					{
-						AdaptedValues& adaptedValuesFather = adaptionMap_[element.father()];
-						adaptedValuesFather.count += 1;
-						storeAdaptionValues(adaptedValues, adaptedValuesFather);
-					}
-
-				}
-        	}
+                for (const auto& element : Dune::elements(levelView))
+                {
+                    //get your map entry
+                    AdaptedValues &adaptedValues = adaptionMap_[element];
+
+                    // put your value in the map
+                    if (element.isLeaf())
+                    {
+                        // get index
+                        int indexI = this->elementIndex(problem, element);
+
+                        storeAdaptionValues(adaptedValues, problem.model().curSol()[indexI]);
+
+                        adaptedValues.count = 1;
+                    }
+                    //Average in father
+                    if (element.level() > 0)
+                    {
+                        AdaptedValues& adaptedValuesFather = adaptionMap_[element.father()];
+                        adaptedValuesFather.count += 1;
+                        storeAdaptionValues(adaptedValues, adaptedValuesFather);
+                    }
+
+                }
+            }
             else
             {
-				for (const auto& entity : Dune::entities(levelView, Dune::Codim<dofCodim>()))
-				{
-					//get your map entry
-					AdaptedValues &adaptedValues = adaptionMap_[entity];
+                for (const auto& entity : Dune::entities(levelView, Dune::Codim<dofCodim>()))
+                {
+                    //get your map entry
+                    AdaptedValues &adaptedValues = adaptionMap_[entity];
 
-					// put your value in the map
-					int indexI = this->dofIndex(problem, entity);
+                    // put your value in the map
+                    int indexI = this->dofIndex(problem, entity);
 
-					storeAdaptionValues(adaptedValues, problem.model().curSol()[indexI]);
+                    storeAdaptionValues(adaptedValues, problem.model().curSol()[indexI]);
 
-					adaptedValues.count = 1;
+                    adaptedValues.count = 1;
 
-				}
+                }
             }
         }
     }
@@ -201,27 +201,27 @@ public:
                     //entry is in map, write in leaf
                     if (element.isLeaf())
                     {
-                    	if(!isBox)
-                    	{
-							AdaptedValues &adaptedValues = adaptionMap_[element];
-							int newIdxI = this->elementIndex(problem, element);
-
-							setAdaptionValues(adaptedValues, problem.model().curSol()[newIdxI]);
-						}
-                    	else
-                    	{
+                        if(!isBox)
+                        {
+                            AdaptedValues &adaptedValues = adaptionMap_[element];
+                            int newIdxI = this->elementIndex(problem, element);
+
+                            setAdaptionValues(adaptedValues, problem.model().curSol()[newIdxI]);
+                        }
+                        else
+                        {
                             unsigned int numSubEntities = element.subEntities(dofCodim);
 
-                        	for(unsigned int i = 0; i < numSubEntities; i++)
-                        	{
+                            for(unsigned int i = 0; i < numSubEntities; i++)
+                            {
                                 auto subEntity = element.template subEntity <dofCodim>(i);
                                 AdaptedValues &adaptedValues = adaptionMap_[subEntity];
-    							int newIdxI = this->dofIndex(problem, subEntity);
+                                int newIdxI = this->dofIndex(problem, subEntity);
 
-    							setAdaptionValues(adaptedValues, problem.model().curSol()[newIdxI]);
+                                setAdaptionValues(adaptedValues, problem.model().curSol()[newIdxI]);
 
-                        	}
-                    	}
+                            }
+                        }
                     }
                 }
                 else
@@ -233,59 +233,59 @@ public:
 
                         if(!isBox)
                         {
-							// create new entry: reconstruct from adaptionMap_[*father] to a new
-							// adaptionMap_[*son]
-							reconstructAdaptionValues(adaptionMap_, epFather, element, problem);
-
-							// access new son
-							AdaptedValues& adaptedValues = adaptionMap_[element];
-							adaptedValues.count = 1;
-
-							// if we are on leaf, store reconstructed values of son in CellData object
-							if (element.isLeaf())
-							{
-								// acess new CellData object
-								int newIdxI = this->elementIndex(problem, element);
-
-								setAdaptionValues(adaptedValues, problem.model().curSol()[newIdxI]);
-							}
+                            // create new entry: reconstruct from adaptionMap_[*father] to a new
+                            // adaptionMap_[*son]
+                            reconstructAdaptionValues(adaptionMap_, epFather, element, problem);
+
+                            // access new son
+                            AdaptedValues& adaptedValues = adaptionMap_[element];
+                            adaptedValues.count = 1;
+
+                            // if we are on leaf, store reconstructed values of son in CellData object
+                            if (element.isLeaf())
+                            {
+                                // acess new CellData object
+                                int newIdxI = this->elementIndex(problem, element);
+
+                                setAdaptionValues(adaptedValues, problem.model().curSol()[newIdxI]);
+                            }
                         }
                         else
                         {
                             unsigned int numSubEntities= element.subEntities(dofCodim);
-                    		const auto geometryI = element.geometry();
-
-                        	for(unsigned int i = 0; i < numSubEntities; i++)
-                        	{
-                        		auto subEntity = element.template subEntity <dofCodim>(i);
-    							AdaptedValues &adaptedValues = adaptionMap_[subEntity];
-
-    							if(adaptedValues.count == 0){
-									LocalPosition dofCenterPos = geometryI.local(subEntity.geometry().center());
-									const LocalFiniteElementCache feCache;
-									Dune::GeometryType geomType = epFather.geometry().type();
-
-									const LocalFiniteElement &localFiniteElement = feCache.get(geomType);
-									std::vector<Dune::FieldVector<Scalar, 1> > shapeVal;
-									localFiniteElement.localBasis().evaluateFunction(dofCenterPos, shapeVal);
-									PrimaryVariables u(0);
-									for (int j = 0; j < shapeVal.size(); ++j)
-									{
-										AdaptedValues & adaptedValuesFather = adaptionMap_[epFather.template subEntity <dofCodim>(j)];
-										u.axpy(shapeVal[j], adaptedValuesFather.u);
-									}
-
-									adaptedValues.u = u;
-									adaptedValues.count = 1;
-    							}
-
-    							if (element.isLeaf())
-    							{
-        							int newIdxI = this->dofIndex(problem, subEntity);
-    								setAdaptionValues(adaptedValues, problem.model().curSol()[newIdxI]);
-    							}
-
-                        	}
+                            const auto geometryI = element.geometry();
+
+                            for(unsigned int i = 0; i < numSubEntities; i++)
+                            {
+                                auto subEntity = element.template subEntity <dofCodim>(i);
+                                AdaptedValues &adaptedValues = adaptionMap_[subEntity];
+
+                                if(adaptedValues.count == 0){
+                                    LocalPosition dofCenterPos = geometryI.local(subEntity.geometry().center());
+                                    const LocalFiniteElementCache feCache;
+                                    Dune::GeometryType geomType = epFather.geometry().type();
+
+                                    const LocalFiniteElement &localFiniteElement = feCache.get(geomType);
+                                    std::vector<Dune::FieldVector<Scalar, 1> > shapeVal;
+                                    localFiniteElement.localBasis().evaluateFunction(dofCenterPos, shapeVal);
+                                    PrimaryVariables u(0);
+                                    for (int j = 0; j < shapeVal.size(); ++j)
+                                    {
+                                        AdaptedValues & adaptedValuesFather = adaptionMap_[epFather.template subEntity <dofCodim>(j)];
+                                        u.axpy(shapeVal[j], adaptedValuesFather.u);
+                                    }
+
+                                    adaptedValues.u = u;
+                                    adaptedValues.count = 1;
+                                }
+
+                                if (element.isLeaf())
+                                {
+                                    int newIdxI = this->dofIndex(problem, subEntity);
+                                    setAdaptionValues(adaptedValues, problem.model().curSol()[newIdxI]);
+                                }
+
+                            }
                         }
                     }
 
@@ -334,15 +334,15 @@ public:
     static void storeAdaptionValues(AdaptedValues& adaptedValues,
                                     AdaptedValues& adaptedValuesFather)
     {
-    	if(!isBox)
-    	{
-			adaptedValuesFather.u += adaptedValues.u;
-			adaptedValuesFather.u /= adaptedValues.count;
-    	}
-    	else
-    	{
-    		adaptedValuesFather.u = adaptedValues.u;
-    	}
+        if(!isBox)
+        {
+            adaptedValuesFather.u += adaptedValues.u;
+            adaptedValuesFather.u /= adaptedValues.count;
+        }
+        else
+        {
+            adaptedValuesFather.u = adaptedValues.u;
+        }
     }
     //! Set adapted values in CellData
     /**
@@ -355,10 +355,10 @@ public:
      */
     static void setAdaptionValues(AdaptedValues& adaptedValues, PrimaryVariables& u)
     {
-    	PrimaryVariables uNew = adaptedValues.u;
-    	uNew /= adaptedValues.count;
+        PrimaryVariables uNew = adaptedValues.u;
+        uNew /= adaptedValues.count;
 
-    	u = uNew;
+        u = uNew;
     }
 
     //! Reconstructs sons entries from data of father cell
@@ -373,13 +373,13 @@ public:
      * \param problem The problem
      */
     static void reconstructAdaptionValues(Dune::PersistentContainer<Grid, AdaptedValues>& adaptionMap,
-    		const Element& father, const Element& son, const Problem& problem)
+            const Element& father, const Element& son, const Problem& problem)
     {
-		AdaptedValues& adaptedValues = adaptionMap[son];
-		AdaptedValues& adaptedValuesFather = adaptionMap[father];
+        AdaptedValues& adaptedValues = adaptionMap[son];
+        AdaptedValues& adaptedValuesFather = adaptionMap[father];
 
-		adaptedValues.u = adaptedValuesFather.u;
-		adaptedValues.u /= adaptedValuesFather.count;
+        adaptedValues.u = adaptedValuesFather.u;
+        adaptedValues.u /= adaptedValuesFather.count;
     }
 
     int dofIndex(const Problem& problem, const DofEntity& entity) const
diff --git a/dumux/implicit/adaptive/gridadaptinitializationindicator.hh b/dumux/implicit/adaptive/gridadaptinitializationindicator.hh
index afad4a6e886795e28913c56c730d5138958ca0eb..04e8d08625ad30432e446ae707e3ba63a15167d6 100644
--- a/dumux/implicit/adaptive/gridadaptinitializationindicator.hh
+++ b/dumux/implicit/adaptive/gridadaptinitializationindicator.hh
@@ -224,7 +224,7 @@ public:
             }
 
             if (!(refineAtSource_ || refineAtFluxBC_ || refineAtDirichletBC_))
-            	continue;
+                continue;
 
             // get the fvGeometry and elementVolVars needed for the bc and source interfaces
             FVElementGeometry fvGeometry;
diff --git a/dumux/implicit/box/intersectiontovertexbc.hh b/dumux/implicit/box/intersectiontovertexbc.hh
index 27ba5658534cf629599f71379a75c3d63cbe8a61..9c5f7fdd2420ae34367ab4adc50b9bf8c717285e 100644
--- a/dumux/implicit/box/intersectiontovertexbc.hh
+++ b/dumux/implicit/box/intersectiontovertexbc.hh
@@ -91,7 +91,7 @@ public:
         values.setAllNeumann();
         int vIdxGlobal = problem_.vertexMapper().index(vertex);
 
-	const BoundaryTypes& bcTypes = vertexBC[vIdxGlobal];
+        const BoundaryTypes& bcTypes = vertexBC[vIdxGlobal];
 
         for (int eqIdx = 0; eqIdx < numEq; eqIdx++)
             if (bcTypes.isDirichlet(eqIdx))
diff --git a/dumux/implicit/co2/co2model.hh b/dumux/implicit/co2/co2model.hh
index 371c7009e8733f13e5bb5c45bb3863fe39f06088..5c88e86e2c93a71a6c716dd9a1dfaad5786c3791 100644
--- a/dumux/implicit/co2/co2model.hh
+++ b/dumux/implicit/co2/co2model.hh
@@ -42,7 +42,7 @@ namespace Dumux
  *   case in the 2p2c model. Instead mole fractions are calculated in the FluidSystem with a given
  *   temperature, pressurem and salinity.
  *   The model is able to use either mole or mass fractions. The property useMoles can be set to either true or false in the
- * 	 problem file. Make sure that the according units are used in the problem setup. useMoles is set to false by default.
+ *   problem file. Make sure that the according units are used in the problem setup. useMoles is set to false by default.
  *
  */
 
diff --git a/dumux/implicit/co2/co2volumevariables.hh b/dumux/implicit/co2/co2volumevariables.hh
index 55b2c1fd65df7b9ac8c06618718d77c5aff5c8d1..c50f4489ece29124cdbc1e9f509047fa76f7ec3d 100644
--- a/dumux/implicit/co2/co2volumevariables.hh
+++ b/dumux/implicit/co2/co2volumevariables.hh
@@ -232,23 +232,23 @@ public:
                 Scalar xwH2O = 1 - xwCO2;
                 ParentType::fluidState_.setMoleFraction(wPhaseIdx, nCompIdx, xwCO2);
                 ParentType::fluidState_.setMoleFraction(wPhaseIdx, wCompIdx, xwH2O);
-			}
+            }
 
-			//Get the phase densities from the FluidSystem and set them in the fluidState
+            //Get the phase densities from the FluidSystem and set them in the fluidState
 
-			Scalar rhoW = FluidSystem::density(ParentType::fluidState_, paramCache, wPhaseIdx);
-			Scalar rhoN = FluidSystem::density(ParentType::fluidState_, paramCache, nPhaseIdx);
+            Scalar rhoW = FluidSystem::density(ParentType::fluidState_, paramCache, wPhaseIdx);
+            Scalar rhoN = FluidSystem::density(ParentType::fluidState_, paramCache, nPhaseIdx);
 
-			ParentType::fluidState_.setDensity(wPhaseIdx, rhoW);
-			ParentType::fluidState_.setDensity(nPhaseIdx, rhoN);
+            ParentType::fluidState_.setDensity(wPhaseIdx, rhoW);
+            ParentType::fluidState_.setDensity(nPhaseIdx, rhoN);
 
-			//Get the phase viscosities from the FluidSystem and set them in the fluidState
+            //Get the phase viscosities from the FluidSystem and set them in the fluidState
 
-			Scalar muW = FluidSystem::viscosity(ParentType::fluidState_, paramCache, wPhaseIdx);
-			Scalar muN = FluidSystem::viscosity(ParentType::fluidState_, paramCache, nPhaseIdx);
+            Scalar muW = FluidSystem::viscosity(ParentType::fluidState_, paramCache, wPhaseIdx);
+            Scalar muN = FluidSystem::viscosity(ParentType::fluidState_, paramCache, nPhaseIdx);
 
-			ParentType::fluidState_.setViscosity(wPhaseIdx, muW);
-			ParentType::fluidState_.setViscosity(nPhaseIdx, muN);
+            ParentType::fluidState_.setViscosity(wPhaseIdx, muW);
+            ParentType::fluidState_.setViscosity(nPhaseIdx, muN);
           }
           else if (phasePresence == wPhaseOnly) {
                // only the wetting phase is present, i.e. wetting phase
diff --git a/dumux/implicit/common/implicitdarcyfluxvariables.hh b/dumux/implicit/common/implicitdarcyfluxvariables.hh
index fd80f14f30b605baf8a0d4260bd8c17766365c73..30b53e56ad2283cce8172c922322c5c42c4aae65 100644
--- a/dumux/implicit/common/implicitdarcyfluxvariables.hh
+++ b/dumux/implicit/common/implicitdarcyfluxvariables.hh
@@ -313,8 +313,8 @@ protected:
         }// loop all phases
     }
 
-    const FVElementGeometry &fvGeometry_;   	//!< Information about the geometry of discretization
-    const unsigned int faceIdx_;            	//!< The index of the sub control volume face
+    const FVElementGeometry &fvGeometry_;       //!< Information about the geometry of discretization
+    const unsigned int faceIdx_;                //!< The index of the sub control volume face
     const bool      onBoundary_;                //!< Specifying whether we are currently on the boundary of the simulation domain
     unsigned int    upstreamIdx_[numPhases] , downstreamIdx_[numPhases]; //!< local index of the upstream / downstream vertex
     Scalar          volumeFlux_[numPhases] ;    //!< Velocity multiplied with normal (magnitude=area)
diff --git a/dumux/implicit/common/implicitforchheimerfluxvariables.hh b/dumux/implicit/common/implicitforchheimerfluxvariables.hh
index 18e3bb9470b4c05119e2bdea8102efaf7976ca08..54ad6efa9d5eaa80841b4c37d6c5d1feabd95885 100644
--- a/dumux/implicit/common/implicitforchheimerfluxvariables.hh
+++ b/dumux/implicit/common/implicitforchheimerfluxvariables.hh
@@ -260,13 +260,13 @@ protected:
      * \param phaseIdx The index of the currently considered phase
      */
      void forchheimerResidual_(GlobalPosition & residual,
-    		 	 	 	 	 const Scalar forchCoeff,
-    		 	 	 	 	 const DimWorldMatrix & sqrtK,
-    		 	 	 	 	 const DimWorldMatrix & K,
-    		 	 	 	 	 const GlobalPosition & velocity,
-    		 	 	 	 	 const ElementVolumeVariables & elemVolVars,
-    		 	 	 	 	 const GlobalPosition & potentialGrad,
-    		 	 	 	 	 const unsigned int phaseIdx) const
+                             const Scalar forchCoeff,
+                             const DimWorldMatrix & sqrtK,
+                             const DimWorldMatrix & K,
+                             const GlobalPosition & velocity,
+                             const ElementVolumeVariables & elemVolVars,
+                             const GlobalPosition & potentialGrad,
+                             const unsigned int phaseIdx) const
      {
          const VolumeVariables upVolVars    = elemVolVars[this->upstreamIdx(phaseIdx)];
          const VolumeVariables downVolVars  = elemVolVars[this->downstreamIdx(phaseIdx)];
@@ -330,11 +330,11 @@ protected:
       * \param phaseIdx The index of the currently considered phase
       */
      void forchheimerDerivative_(DimWorldMatrix & derivative,
-    		 	 	 	 	 	 const Scalar forchCoeff,
-    		 	 	 	 	 	 const DimWorldMatrix & sqrtK,
-    		 	 	 	 	 	 const GlobalPosition & velocity,
-    		 	 	 	 	 	 const ElementVolumeVariables & elemVolVars,
-    		 	 	 	 	 	 const unsigned int phaseIdx) const
+                                 const Scalar forchCoeff,
+                                 const DimWorldMatrix & sqrtK,
+                                 const GlobalPosition & velocity,
+                                 const ElementVolumeVariables & elemVolVars,
+                                 const unsigned int phaseIdx) const
      {
          const VolumeVariables upVolVars    = elemVolVars[this->upstreamIdx(phaseIdx)];
          const VolumeVariables downVolVars  = elemVolVars[this->downstreamIdx(phaseIdx)];
diff --git a/dumux/implicit/common/implicitmodel.hh b/dumux/implicit/common/implicitmodel.hh
index a0f09fa37ffeb5259aba99ab56cabbcee63130c7..24e75b700ee0a590de757f428f169ad98c026a01 100644
--- a/dumux/implicit/common/implicitmodel.hh
+++ b/dumux/implicit/common/implicitmodel.hh
@@ -313,7 +313,7 @@ public:
 
     void adaptVariableSize()
     {
-    	uCur_.resize(numDofs());
+        uCur_.resize(numDofs());
     }
 
     /*!
@@ -456,9 +456,9 @@ public:
      */
     void updateBegin()
     {
-    	if(GET_PROP_VALUE(TypeTag, AdaptiveGrid) && problem_().gridAdapt().wasAdapted())
-    	{
-    		uPrev_ = uCur_;
+        if(GET_PROP_VALUE(TypeTag, AdaptiveGrid) && problem_().gridAdapt().wasAdapted())
+        {
+            uPrev_ = uCur_;
 
             updateBoundaryIndices_();
 
@@ -480,7 +480,7 @@ public:
                           false);
             }
 
-    	}
+        }
 
     }
 
diff --git a/dumux/implicit/cornerpoint/cpdarcyfluxvariables.hh b/dumux/implicit/cornerpoint/cpdarcyfluxvariables.hh
index 682ccef2a236e3e1846c9d779368bf9ace5e46b2..2532e5928fb384db2ff3f0a7352db786a78282ac 100644
--- a/dumux/implicit/cornerpoint/cpdarcyfluxvariables.hh
+++ b/dumux/implicit/cornerpoint/cpdarcyfluxvariables.hh
@@ -39,7 +39,7 @@ namespace Dumux
 
 namespace Properties
 {
-// forward declaration of properties 
+// forward declaration of properties
 NEW_PROP_TAG(ImplicitMobilityUpwindWeight);
 NEW_PROP_TAG(SpatialParams);
 NEW_PROP_TAG(NumPhases);
@@ -48,7 +48,7 @@ NEW_PROP_TAG(ProblemEnableGravity);
 
 /*!
  * \ingroup ImplicitFluxVariables
- * \brief Evaluates the normal component of the Darcy velocity 
+ * \brief Evaluates the normal component of the Darcy velocity
  * on a (sub)control volume face.
  */
 template <class TypeTag>
@@ -58,7 +58,7 @@ class CpDarcyFluxVariables
     typedef typename GET_PROP_TYPE(TypeTag, SpatialParams) SpatialParams;
     typedef typename GET_PROP_TYPE(TypeTag, ElementVolumeVariables) ElementVolumeVariables;
     typedef typename GET_PROP_TYPE(TypeTag, VolumeVariables) VolumeVariables;
-    
+
     typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
     typedef typename GridView::template Codim<0>::Entity Element;
 
@@ -131,7 +131,7 @@ public:
      */
     const unsigned int downstreamIdx(const unsigned phaseIdx) const
     { return downstreamIdx_[phaseIdx]; }
-    
+
     /*!
      * \brief Return the local index of the upstream control volume
      *        for a given phase.
@@ -264,8 +264,8 @@ protected:
         } // over loop all phases
     }
 
-    const FVElementGeometry &fvGeometry_;   	//!< Information about the geometry of discretization
-    const unsigned int faceIdx_;            	//!< The index of the sub control volume face
+    const FVElementGeometry &fvGeometry_;       //!< Information about the geometry of discretization
+    const unsigned int faceIdx_;                //!< The index of the sub control volume face
     const bool      onBoundary_;                //!< Specifying whether we are currently on the boundary of the simulation domain
     unsigned int    upstreamIdx_[numPhases] , downstreamIdx_[numPhases]; //!< local index of the upstream / downstream vertex
     Scalar          volumeFlux_[numPhases] ;    //!< Velocity multiplied with normal (magnitude=area)
diff --git a/dumux/implicit/cornerpoint/cpelementvolumevariables.hh b/dumux/implicit/cornerpoint/cpelementvolumevariables.hh
index 2751b3fd53b9e36934e1ee13871f68f64910ac2a..00913e6ad6c9667ba8fa4b877005c3c44467fb81 100644
--- a/dumux/implicit/cornerpoint/cpelementvolumevariables.hh
+++ b/dumux/implicit/cornerpoint/cpelementvolumevariables.hh
@@ -94,11 +94,11 @@ public:
         // only treat boundary if current solution is evaluated
         if (!oldSol)
         {
-            // check if element intersects with the boundary 
+            // check if element intersects with the boundary
             ElementBoundaryTypes elemBCTypes;
             elemBCTypes.update(problem, element);
-            if (elemBCTypes.hasDirichlet() 
-                || elemBCTypes.hasNeumann() 
+            if (elemBCTypes.hasDirichlet()
+                || elemBCTypes.hasNeumann()
                 || elemBCTypes.hasOutflow())
             {
                 const int numFaces = 6;
@@ -127,7 +127,7 @@ public:
                                                          /*scvIdx=*/0,
                                                          oldSol);
                     }
-                    else 
+                    else
                     {
                         (*this)[indexInVariables] = (*this)[0];
                     }
diff --git a/dumux/implicit/cornerpoint/cpfvelementgeometry.hh b/dumux/implicit/cornerpoint/cpfvelementgeometry.hh
index 101a539659a73a9567f2e8a67149e54207fa5853..d60d0604a883d5b82b79cd371547e9f482732a3a 100644
--- a/dumux/implicit/cornerpoint/cpfvelementgeometry.hh
+++ b/dumux/implicit/cornerpoint/cpfvelementgeometry.hh
@@ -91,7 +91,7 @@ public:
     SubControlVolumeFace subContVolFace[maxNE]; //!< data of the sub control volume faces
     BoundaryFace boundaryFace[maxBF]; //!< data of the boundary faces
     int numScv; //!< number of subcontrol volumes
-    int numScvf; //!< number of inner-domain subcontrolvolume faces 
+    int numScvf; //!< number of inner-domain subcontrolvolume faces
     int numNeighbors; //!< number of neighboring elements including the element itself
     std::vector<Element> neighbors; //!< stores the neighboring elements
 
@@ -111,7 +111,7 @@ public:
         subContVol[0].inner = true;
         subContVol[0].volume = elementVolume;
 
-        // initialize neighbors list with self: 
+        // initialize neighbors list with self:
         numNeighbors = 1;
         neighbors.clear();
         neighbors.reserve(maxNE);
diff --git a/dumux/implicit/mpnc/energy/mpncfluxvariablesenergy.hh b/dumux/implicit/mpnc/energy/mpncfluxvariablesenergy.hh
index 75dfbc3c17a66929f4d78970cb5a3b70eab71ad0..066cb55620e90b1235d42d001a506fba73faf8ea 100644
--- a/dumux/implicit/mpnc/energy/mpncfluxvariablesenergy.hh
+++ b/dumux/implicit/mpnc/energy/mpncfluxvariablesenergy.hh
@@ -150,10 +150,10 @@ public:
 
         lambdaEff_ = 0;
         calculateEffThermalConductivity_(problem,
-        								 element,
-        								 fvGeometry,
-        								 face,
-        								 elemVolVars);
+                                         element,
+                                         fvGeometry,
+                                         face,
+                                         elemVolVars);
     }
 
     /*!
@@ -172,13 +172,13 @@ public:
 
 protected:
     /*!
-	 * \brief Calculate the effective thermal conductivity of
-	 * 		  the porous medium plus residing phases \f$[W/mK]\f$.
-	 *		  This basically means to access the model for averaging
-	 *		  the individual conductivities, set by the property ThermalConductivityModel.
-	 *		  Except the adapted arguments, this is the same function
-	 *		  as used in the implicit TwoPTwoCNIFluxVariables.
-	 */
+     * \brief Calculate the effective thermal conductivity of
+     *        the porous medium plus residing phases \f$[W/mK]\f$.
+     *        This basically means to access the model for averaging
+     *        the individual conductivities, set by the property ThermalConductivityModel.
+     *        Except the adapted arguments, this is the same function
+     *        as used in the implicit TwoPTwoCNIFluxVariables.
+     */
     void calculateEffThermalConductivity_(const Problem &problem,
                                           const Element &element,
                                           const FVElementGeometry & fvGeometry,
diff --git a/dumux/implicit/mpnc/energy/mpncfluxvariablesenergykinetic.hh b/dumux/implicit/mpnc/energy/mpncfluxvariablesenergykinetic.hh
index 7dc2e668ac606bc02fcf4d2d65cbc7cec8c6b7e2..67b6e1daa62dbbbab7a72e50f4c16cdfa3556aa4 100644
--- a/dumux/implicit/mpnc/energy/mpncfluxvariablesenergykinetic.hh
+++ b/dumux/implicit/mpnc/energy/mpncfluxvariablesenergykinetic.hh
@@ -199,10 +199,10 @@ public:
 
         lambdaEff_ = 0;
         calculateEffThermalConductivity_(problem,
-        								 element,
-        								 fvGeometry,
-        								 face,
-        								 elemVolVars);
+                                         element,
+                                         fvGeometry,
+                                         face,
+                                         elemVolVars);
 
 
     }
@@ -226,13 +226,13 @@ public:
 
 protected:
     /*!
-	 * \brief Calculate the effective thermal conductivity of
-	 * 		  the porous medium plus residing phases \f$[W/mK]\f$.
-	 *		  This basically means to access the model for averaging
-	 *		  the individual conductivities, set by the property ThermalConductivityModel.
-	 *		  Except the adapted arguments, this is the same function
-	 *		  as used in the implicit TwoPTwoCNIFluxVariables.
-	 */
+     * \brief Calculate the effective thermal conductivity of
+     *        the porous medium plus residing phases \f$[W/mK]\f$.
+     *        This basically means to access the model for averaging
+     *        the individual conductivities, set by the property ThermalConductivityModel.
+     *        Except the adapted arguments, this is the same function
+     *        as used in the implicit TwoPTwoCNIFluxVariables.
+     */
     void calculateEffThermalConductivity_(const Problem &problem,
                                           const Element &element,
                                           const FVElementGeometry & fvGeometry,
diff --git a/dumux/implicit/mpnc/energy/mpncindicesenergykinetic.hh b/dumux/implicit/mpnc/energy/mpncindicesenergykinetic.hh
index f69892fdd2928a9e29f7378c2dfd5905ea2fd6ae..bf7abb03a495cbc7f8d6a829bea70db1845531f4 100644
--- a/dumux/implicit/mpnc/energy/mpncindicesenergykinetic.hh
+++ b/dumux/implicit/mpnc/energy/mpncindicesenergykinetic.hh
@@ -30,7 +30,7 @@ namespace Dumux
 
 /*!
  * \brief The indices required for the energy equation. Specialization for the case of
- * 		  *3* energy balance equations.
+ *        *3* energy balance equations.
  */
 template <int PVOffset>
 struct MPNCEnergyIndices<PVOffset, /*enableEnergy=*/true, /*numEnergyEquations=*/3>
@@ -63,7 +63,7 @@ public:
 
 /*!
  * \brief The indices required for the energy equation. Specialization for the case of
- * 		  *2* energy balance equations.
+ *        *2* energy balance equations.
  */
 template <int PVOffset>
 struct MPNCEnergyIndices<PVOffset, /*enableEnergy=*/true, /*numEnergyEquations=*/2>
diff --git a/dumux/implicit/mpnc/energy/mpnclocalresidualenergykinetic.hh b/dumux/implicit/mpnc/energy/mpnclocalresidualenergykinetic.hh
index bd5bac4e84beb881c36ee6f5c5b8ddc64b8737d1..bfe97a2a93911326618d85eed9b0769a246105d5 100644
--- a/dumux/implicit/mpnc/energy/mpnclocalresidualenergykinetic.hh
+++ b/dumux/implicit/mpnc/energy/mpnclocalresidualenergykinetic.hh
@@ -127,10 +127,10 @@ public:
         }
         else
             DUNE_THROW(Dune::NotImplemented,
-            		"wrong index");
+                    "wrong index");
 
         if (!std::isfinite(storage[energyEq0Idx+phaseIdx]))
-        	DUNE_THROW(NumericalProblem, "Calculated non-finite storage");
+            DUNE_THROW(NumericalProblem, "Calculated non-finite storage");
 
     }
     /*!
@@ -166,7 +166,7 @@ public:
                                   elemVolVars,
                                   energyEqIdx);
             if (!std::isfinite(flux[energyEq0Idx + energyEqIdx]))
-            	DUNE_THROW(NumericalProblem, "Calculated non-finite flux in phase " << energyEqIdx);
+                DUNE_THROW(NumericalProblem, "Calculated non-finite flux in phase " << energyEqIdx);
         }
     }
     /*!
@@ -215,9 +215,9 @@ public:
 
         const Scalar massUpwindWeight_ = GET_PARAM_FROM_GROUP(TypeTag, Scalar, Implicit, MassUpwindWeight);
         flux[energyEq0Idx + phaseIdx] += massFlux *
-        									(massUpwindWeight_ * transportedThingUp
-        											+
-        									(1.-massUpwindWeight_) * transportedThingDn ) ;
+                                            (massUpwindWeight_ * transportedThingUp
+                                                    +
+                                            (1.-massUpwindWeight_) * transportedThingDn ) ;
     }
     /*!
         * \brief The heat conduction in the phase
@@ -350,7 +350,7 @@ public:
 
 
             if (!std::isfinite(source[energyEq0Idx + phaseIdx]))
-            	DUNE_THROW(NumericalProblem, "Calculated non-finite source, " << "Tw="<< Tw << " Tn="<< Tn<< " Ts="<< Ts);
+                DUNE_THROW(NumericalProblem, "Calculated non-finite source, " << "Tw="<< Tw << " Tn="<< Tn<< " Ts="<< Ts);
         }// end phases
 
 #define MASS_ENERGY_TRANSPORT 1
@@ -468,7 +468,7 @@ public:
             * volVars.solidHeatCapacity();
 
         if (!std::isfinite(storage[energyEqSolidIdx]))
-        	DUNE_THROW(NumericalProblem, "Calculated non-finite storage");
+            DUNE_THROW(NumericalProblem, "Calculated non-finite storage");
     }
 
     /*! \brief Calculate the storage for all mass balance equations
@@ -494,10 +494,10 @@ public:
         }
         else
             DUNE_THROW(Dune::NotImplemented,
-            		"wrong index");
+                    "wrong index");
 
         if (!std::isfinite(storage[energyEq0Idx]))
-        	DUNE_THROW(NumericalProblem, "Calculated non-finite storage");
+            DUNE_THROW(NumericalProblem, "Calculated non-finite storage");
     }
 
     /*!\brief Evaluates the total flux of all conservation quantities
@@ -533,17 +533,17 @@ public:
                                   energyEqIdx);
 
             if (!std::isfinite(flux[energyEq0Idx + energyEqIdx]))
-            	DUNE_THROW(NumericalProblem, "Calculated non-finite flux in phase " << energyEqIdx);
+                DUNE_THROW(NumericalProblem, "Calculated non-finite flux in phase " << energyEqIdx);
         }
     }
 
     /*! \brief the advective Flux of the enthalpy
-	  *        \param flux The flux over the SCV (sub-control-volume) face for each component
-	  *        \param fluxVars The flux Variables
-	  *        \param elemVolVars The volume variables of the current element
-	  *        \param phaseIdx The local index of the phases
-	  *        \param molarComponentValuesMassTransport
-	  */
+      *        \param flux The flux over the SCV (sub-control-volume) face for each component
+      *        \param fluxVars The flux Variables
+      *        \param elemVolVars The volume variables of the current element
+      *        \param phaseIdx The local index of the phases
+      *        \param molarComponentValuesMassTransport
+      */
     static void computePhaseEnthalpyFlux(PrimaryVariables & flux,
                                          const FluxVariables & fluxVars,
                                          const ElementVolumeVariables & elemVolVars,
@@ -579,9 +579,9 @@ public:
 
         const Scalar massUpwindWeight_ = GET_PARAM_FROM_GROUP(TypeTag, Scalar, Implicit, MassUpwindWeight);
         flux[energyEq0Idx] += massFlux *
-        									(massUpwindWeight_ * transportedThingUp
-        											+
-        									(1.-massUpwindWeight_) * transportedThingDn ) ;
+                                            (massUpwindWeight_ * transportedThingUp
+                                                    +
+                                            (1.-massUpwindWeight_) * transportedThingDn ) ;
     }
 
     /*! \brief The heat conduction in the phase
@@ -629,7 +629,7 @@ public:
                                "wrong index");
     }
 
-	/*! \brief Calculate the source term of the equation
+    /*! \brief Calculate the source term of the equation
    *
    * \param source The source/sink in the sub-control volume for each component
    * \param volVars The volume variables
@@ -664,119 +664,119 @@ public:
    */
     static Scalar qsf(const VolumeVariables & volVars)
     {
-    	const FluidState & fs = volVars.fluidState() ;
-		const Scalar characteristicLength   = volVars.characteristicLength()  ;
+        const FluidState & fs = volVars.fluidState() ;
+        const Scalar characteristicLength   = volVars.characteristicLength()  ;
 
-		// Shi & Wang, Transport in porous media (2011)
-		const Scalar as = 6.0 * (1.0-volVars.porosity()) / characteristicLength ;
+        // Shi & Wang, Transport in porous media (2011)
+        const Scalar as = 6.0 * (1.0-volVars.porosity()) / characteristicLength ;
 
-		const Scalar TFluid 	= volVars.temperature(temperatureFluidIdx);
-		const Scalar TSolid 	= volVars.temperature(temperatureSolidIdx);
+        const Scalar TFluid     = volVars.temperature(temperatureFluidIdx);
+        const Scalar TSolid     = volVars.temperature(temperatureSolidIdx);
 
-		const Scalar satW 		= fs.saturation(wPhaseIdx) ;
-		const Scalar satN 		= fs.saturation(nPhaseIdx) ;
+        const Scalar satW       = fs.saturation(wPhaseIdx) ;
+        const Scalar satN       = fs.saturation(nPhaseIdx) ;
 
-		const Scalar eps = 1e-6 ;
-		Scalar solidToFluidEnergyExchange ;
+        const Scalar eps = 1e-6 ;
+        Scalar solidToFluidEnergyExchange ;
 
-		Scalar fluidConductivity ;
-		if (satW < 1.0 - eps)
-			fluidConductivity = volVars.thermalConductivity(nPhaseIdx) ;
-		else if (satW >= 1.0 - eps)
-			fluidConductivity = volVars.thermalConductivity(wPhaseIdx) ;
-		else
-			DUNE_THROW(Dune::NotImplemented,
-					   "wrong range");
+        Scalar fluidConductivity ;
+        if (satW < 1.0 - eps)
+            fluidConductivity = volVars.thermalConductivity(nPhaseIdx) ;
+        else if (satW >= 1.0 - eps)
+            fluidConductivity = volVars.thermalConductivity(wPhaseIdx) ;
+        else
+            DUNE_THROW(Dune::NotImplemented,
+                       "wrong range");
 
-		const Scalar factorEnergyTransfer   = volVars.factorEnergyTransfer()  ;
+        const Scalar factorEnergyTransfer   = volVars.factorEnergyTransfer()  ;
 
-		solidToFluidEnergyExchange = factorEnergyTransfer * (TSolid - TFluid) / characteristicLength * as * fluidConductivity ;
+        solidToFluidEnergyExchange = factorEnergyTransfer * (TSolid - TFluid) / characteristicLength * as * fluidConductivity ;
 
-		const Scalar epsRegul = 1e-3 ;
+        const Scalar epsRegul = 1e-3 ;
 
-		if (satW < (0 + eps) )
-			solidToFluidEnergyExchange *=  volVars.nusseltNumber(nPhaseIdx) ;
+        if (satW < (0 + eps) )
+            solidToFluidEnergyExchange *=  volVars.nusseltNumber(nPhaseIdx) ;
 
-		else if ( (satW >= 0 + eps) and (satW < 1.0-eps) ){
-			solidToFluidEnergyExchange *=  (volVars.nusseltNumber(nPhaseIdx) * satN );
-			Scalar qBoil ;
+        else if ( (satW >= 0 + eps) and (satW < 1.0-eps) ){
+            solidToFluidEnergyExchange *=  (volVars.nusseltNumber(nPhaseIdx) * satN );
+            Scalar qBoil ;
 
-			if (satW<=epsRegul){// regularize
-				Spline sp(0.0, 						epsRegul, 							// x1, x2
-						  QBoilFunc(volVars, 0.0), 	QBoilFunc(volVars, epsRegul), 		// y1, y2
-						  0.0, 						dQBoil_dSw(volVars, epsRegul) ); 	// m1, m2
+            if (satW<=epsRegul){// regularize
+                Spline sp(0.0,                      epsRegul,                           // x1, x2
+                          QBoilFunc(volVars, 0.0),  QBoilFunc(volVars, epsRegul),       // y1, y2
+                          0.0,                      dQBoil_dSw(volVars, epsRegul) );    // m1, m2
 
-				qBoil = sp.eval(satW) ;
-			}
+                qBoil = sp.eval(satW) ;
+            }
 
-			else if (satW>= (1.0-epsRegul) ){// regularize
-				Spline sp(1.0-epsRegul, 									1.0, 	// x1, x2
-						  QBoilFunc(volVars, 1.0-epsRegul), 				0.0, 	// y1, y2
-						  dQBoil_dSw(volVars, 1.0-epsRegul),					0.0 ); 		// m1, m2
+            else if (satW>= (1.0-epsRegul) ){// regularize
+                Spline sp(1.0-epsRegul,                                     1.0,    // x1, x2
+                          QBoilFunc(volVars, 1.0-epsRegul),                 0.0,    // y1, y2
+                          dQBoil_dSw(volVars, 1.0-epsRegul),                    0.0 );      // m1, m2
 
-				qBoil = sp.eval(satW) ;
-			}
-			else
-				qBoil = QBoilFunc(volVars, satW)  ;
+                qBoil = sp.eval(satW) ;
+            }
+            else
+                qBoil = QBoilFunc(volVars, satW)  ;
 
-			solidToFluidEnergyExchange += qBoil;
-		}
-		else if (satW >= 1.0-eps)
-			solidToFluidEnergyExchange *=  volVars.nusseltNumber(wPhaseIdx) ;
-		else
-			DUNE_THROW(Dune::NotImplemented,
-					   "wrong range");
+            solidToFluidEnergyExchange += qBoil;
+        }
+        else if (satW >= 1.0-eps)
+            solidToFluidEnergyExchange *=  volVars.nusseltNumber(wPhaseIdx) ;
+        else
+            DUNE_THROW(Dune::NotImplemented,
+                       "wrong range");
 
-		if (!std::isfinite(solidToFluidEnergyExchange))
-		            	DUNE_THROW(NumericalProblem, "Calculated non-finite source, " << "TFluid="<< TFluid << " TSolid="<< TSolid  );
+        if (!std::isfinite(solidToFluidEnergyExchange))
+                        DUNE_THROW(NumericalProblem, "Calculated non-finite source, " << "TFluid="<< TFluid << " TSolid="<< TSolid  );
 
-		return solidToFluidEnergyExchange ;
+        return solidToFluidEnergyExchange ;
     }
 
-	/*! \brief Calculate the energy transfer during boiling, i.e. latent heat
+    /*! \brief Calculate the energy transfer during boiling, i.e. latent heat
    *
    * \param volVars The volume variables
    * \param satW The wetting phase saturation. Not taken from volVars, because we regularize.
    */
     static Scalar QBoilFunc(const VolumeVariables & volVars,
-    						const  Scalar satW)
+                            const  Scalar satW)
     {
-    	// using saturation as input (instead of from volVars)
-    	// in order to make regularization (evaluation at different points) easyer
+        // using saturation as input (instead of from volVars)
+        // in order to make regularization (evaluation at different points) easyer
         const FluidState & fs = volVars.fluidState() ;
-	const Scalar g( 9.81 ) ;
-	const Scalar gamma(0.0589) ;
-        const Scalar TSolid 	= volVars.temperature(temperatureSolidIdx);
+    const Scalar g( 9.81 ) ;
+    const Scalar gamma(0.0589) ;
+        const Scalar TSolid     = volVars.temperature(temperatureSolidIdx);
         const Scalar characteristicLength   = volVars.characteristicLength()  ;
 
-	const Scalar as = 6.0 * (1.0-volVars.porosity()) / characteristicLength ;
-	const Scalar mul = fs.viscosity(wPhaseIdx) ;
-    	const Scalar deltahv = fs.enthalpy(nPhaseIdx) - fs.enthalpy(wPhaseIdx);
-    	const Scalar deltaRho = fs.density(wPhaseIdx) - fs.density(nPhaseIdx) ;
-	const Scalar firstBracket = std::pow(g * deltaRho / gamma, 0.5);
-	const Scalar cp = FluidSystem::heatCapacity(fs, wPhaseIdx) ;
+    const Scalar as = 6.0 * (1.0-volVars.porosity()) / characteristicLength ;
+    const Scalar mul = fs.viscosity(wPhaseIdx) ;
+        const Scalar deltahv = fs.enthalpy(nPhaseIdx) - fs.enthalpy(wPhaseIdx);
+        const Scalar deltaRho = fs.density(wPhaseIdx) - fs.density(nPhaseIdx) ;
+    const Scalar firstBracket = std::pow(g * deltaRho / gamma, 0.5);
+    const Scalar cp = FluidSystem::heatCapacity(fs, wPhaseIdx) ;
     // This use of Tsat is only justified if the fluid is always boiling (tsat equals boiling conditions)
     // If a different state is to be simulated, please use the actual fluid temperature instead.
-	const Scalar Tsat = FluidSystem::vaporTemperature(fs, nPhaseIdx ) ;
-	const Scalar deltaT = TSolid - Tsat ;
-	const Scalar secondBracket = std::pow( (cp *deltaT / (0.006 * deltahv)  ) , 3.0 ) ;
-	const Scalar Prl = volVars.prandtlNumber(wPhaseIdx) ;
-	const Scalar thirdBracket = std::pow( 1/Prl , (1.7/0.33) );
-	const Scalar QBoil = satW * as * mul * deltahv * firstBracket * secondBracket * thirdBracket ;
-    	return QBoil;
+    const Scalar Tsat = FluidSystem::vaporTemperature(fs, nPhaseIdx ) ;
+    const Scalar deltaT = TSolid - Tsat ;
+    const Scalar secondBracket = std::pow( (cp *deltaT / (0.006 * deltahv)  ) , 3.0 ) ;
+    const Scalar Prl = volVars.prandtlNumber(wPhaseIdx) ;
+    const Scalar thirdBracket = std::pow( 1/Prl , (1.7/0.33) );
+    const Scalar QBoil = satW * as * mul * deltahv * firstBracket * secondBracket * thirdBracket ;
+        return QBoil;
     }
 
-	/*! \brief Calculate the derivative of the energy transfer function during boiling. Needed for regularization.
+    /*! \brief Calculate the derivative of the energy transfer function during boiling. Needed for regularization.
    *
    * \param volVars The volume variables
    * \param satW The wetting phase saturation. Not taken from volVars, because we regularize.
    */
     static Scalar dQBoil_dSw(const VolumeVariables & volVars,
-    							const Scalar satW)
+                                const Scalar satW)
     {
-    	// on the fly derive w.r.t. Sw.
-    	// Only linearly depending on it (directly)
-    	return (QBoilFunc(volVars, satW) / satW ) ;
+        // on the fly derive w.r.t. Sw.
+        // Only linearly depending on it (directly)
+        return (QBoilFunc(volVars, satW) / satW ) ;
     }
 };
 
diff --git a/dumux/implicit/mpnc/energy/mpncvolumevariablesenergykinetic.hh b/dumux/implicit/mpnc/energy/mpncvolumevariablesenergykinetic.hh
index f9bb56539b20a330e042c67f5a697805b4ee5aa2..9e62d336b5f7b8de6a91c152bdf56447e233b415 100644
--- a/dumux/implicit/mpnc/energy/mpncvolumevariablesenergykinetic.hh
+++ b/dumux/implicit/mpnc/energy/mpncvolumevariablesenergykinetic.hh
@@ -114,7 +114,7 @@ public:
 
         temperature_[sPhaseIdx] = priVars[temperature0Idx + sPhaseIdx];
 
-	  Valgrind::CheckDefined(temperature_);
+      Valgrind::CheckDefined(temperature_);
     }
 
     /*!
@@ -141,7 +141,7 @@ public:
 
         for(int phaseIdx =0; phaseIdx<numPhases; ++phaseIdx){
             fluidThermalConductivity_[phaseIdx] =
-	      FluidSystem::thermalConductivity(fluidState, paramCache, phaseIdx);
+          FluidSystem::thermalConductivity(fluidState, paramCache, phaseIdx);
         }
         Valgrind::CheckDefined(fluidThermalConductivity_);
 
@@ -319,7 +319,7 @@ public:
 
         fluidState.setTemperature(priVars[temperature0Idx]);
 
-	  Valgrind::CheckDefined(temperature_);
+      Valgrind::CheckDefined(temperature_);
     }
 
     /*!
@@ -346,7 +346,7 @@ public:
 
         for(int phaseIdx =0; phaseIdx<numPhases; ++phaseIdx){
             fluidThermalConductivity_[phaseIdx] =
-            		FluidSystem::thermalConductivity(fluidState, paramCache, phaseIdx);
+                    FluidSystem::thermalConductivity(fluidState, paramCache, phaseIdx);
         }
         Valgrind::CheckDefined(fluidThermalConductivity_);
 
diff --git a/dumux/implicit/mpnc/energy/mpncvtkwriterenergykinetic.hh b/dumux/implicit/mpnc/energy/mpncvtkwriterenergykinetic.hh
index 37a5ab4a1012b0e9a1a7e4b7ddadaece5492caab..ab846374431e9b0850f5b3a359d2a5e861f26482 100644
--- a/dumux/implicit/mpnc/energy/mpncvtkwriterenergykinetic.hh
+++ b/dumux/implicit/mpnc/energy/mpncvtkwriterenergykinetic.hh
@@ -377,10 +377,10 @@ public:
             const unsigned int vIdxGlobal = this->problem_.vertexMapper().subIndex(element, localVertexIdx, dim);
             const VolumeVariables &volVars = elemVolVars[localVertexIdx];
 
-        	qBoil_[vIdxGlobal] = LocalResidual::QBoilFunc(volVars, volVars.fluidState().saturation(wPhaseIdx));
-        	qsf_[vIdxGlobal] = LocalResidual::qsf(volVars);
+            qBoil_[vIdxGlobal] = LocalResidual::QBoilFunc(volVars, volVars.fluidState().saturation(wPhaseIdx));
+            qsf_[vIdxGlobal] = LocalResidual::qsf(volVars);
 
-        	for (int phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
+            for (int phaseIdx = 0; phaseIdx < numPhases; ++ phaseIdx) {
                 enthalpy_[phaseIdx][vIdxGlobal]          = volVars.fluidState().enthalpy(phaseIdx);
                 internalEnergy_[phaseIdx][vIdxGlobal]    = volVars.fluidState().internalEnergy(phaseIdx);
                 reynoldsNumber_[phaseIdx][vIdxGlobal]    = volVars.reynoldsNumber(phaseIdx);
@@ -470,10 +470,10 @@ private:
                                    EnergyEqVector & buffer,
                                    bool vertexCentered = true)
     {
-    	static const char *name[] = {
-    	            "fluid",
-    	            "solid"
-    	        };
+        static const char *name[] = {
+                    "fluid",
+                    "solid"
+                };
 
         for (int energyEqIdx = 0; energyEqIdx < numEnergyEqs; ++energyEqIdx) {
             std::ostringstream oss;
diff --git a/dumux/implicit/mpnc/mass/mpnclocalresidualmass.hh b/dumux/implicit/mpnc/mass/mpnclocalresidualmass.hh
index c2bdf90d5985214e01cd420bda65abf6727d5ff6..1f5743739f6a53101dbb04e171826f7182c482fd 100644
--- a/dumux/implicit/mpnc/mass/mpnclocalresidualmass.hh
+++ b/dumux/implicit/mpnc/mass/mpnclocalresidualmass.hh
@@ -183,9 +183,9 @@ if (!std::isfinite(volumeFlux))
                         ((     massUpwindWeight)*up.fluidState().molarity(phaseIdx, compIdx)
                                 +
                         (  1. - massUpwindWeight)*dn.fluidState().molarity(phaseIdx, compIdx) );
-						if (!std::isfinite(flux[compIdx]))
-							DUNE_THROW(NumericalProblem, "Calculated non-finite normal flux in phase " <<  phaseIdx << " comp " << compIdx << "T: "<<  up.fluidState().temperature(phaseIdx) << "S "<<up.fluidState().saturation(phaseIdx)  ) ;
-			}
+                        if (!std::isfinite(flux[compIdx]))
+                            DUNE_THROW(NumericalProblem, "Calculated non-finite normal flux in phase " <<  phaseIdx << " comp " << compIdx << "T: "<<  up.fluidState().temperature(phaseIdx) << "S "<<up.fluidState().saturation(phaseIdx)  ) ;
+            }
         }
     }
 
@@ -371,10 +371,10 @@ public:
     static void computeSource(PrimaryVariables &source,
                               const VolumeVariables &volVars)
     {
-//    	static_assert(not enableKineticEnergy, // enableKinetic is disabled, in this specialization
-//    			"In the case of kinetic energy transfer the advective energy transport between the phases has to be considered. "
-//    			"It is hard (technically) to say how much mass got transfered in the case of chemical equilibrium. "
-//    			"Therefore, kineticEnergy and no kinetic mass does not fit (yet).");
+//      static_assert(not enableKineticEnergy, // enableKinetic is disabled, in this specialization
+//              "In the case of kinetic energy transfer the advective energy transport between the phases has to be considered. "
+//              "It is hard (technically) to say how much mass got transfered in the case of chemical equilibrium. "
+//              "Therefore, kineticEnergy and no kinetic mass does not fit (yet).");
 
 
         // mass transfer is not considered in this mass module
diff --git a/dumux/implicit/mpnc/mass/mpnclocalresidualmasskinetic.hh b/dumux/implicit/mpnc/mass/mpnclocalresidualmasskinetic.hh
index 61b4dbd1e11f25901444bd6291509e60882b6f2a..09400252acfad48805f553fc6b383b7d4acdbb17 100644
--- a/dumux/implicit/mpnc/mass/mpnclocalresidualmasskinetic.hh
+++ b/dumux/implicit/mpnc/mass/mpnclocalresidualmasskinetic.hh
@@ -106,11 +106,11 @@ public:
 
         // copy to the primary variables
         for (int compIdx = 0; compIdx < numComponents; ++compIdx){
-        	storage[conti0EqIdx + phaseIdx*numComponents + compIdx]
-        	        += phaseComponentValues[compIdx];
+            storage[conti0EqIdx + phaseIdx*numComponents + compIdx]
+                    += phaseComponentValues[compIdx];
 
             if (!std::isfinite(storage[conti0EqIdx + phaseIdx*numComponents + compIdx]))
-            	DUNE_THROW(NumericalProblem, "Calculated non-finite storage");
+                DUNE_THROW(NumericalProblem, "Calculated non-finite storage");
         }
 
        Valgrind::CheckDefined(storage);
@@ -145,7 +145,7 @@ public:
                 Valgrind::CheckDefined(flux);
 
                 if (!std::isfinite(flux[conti0EqIdx + phaseIdx*numComponents + compIdx]))
-                	DUNE_THROW(NumericalProblem, "Calculated non-finite flux");
+                    DUNE_THROW(NumericalProblem, "Calculated non-finite flux");
             }
 
             // Right now I think that adding the two contributions individually into the flux is best for debugging and understanding.
@@ -276,7 +276,7 @@ public:
                         source[eqIdx] += componentIntoPhaseMassTransfer[phaseIdx][compIdx] ;
 
                         if (!std::isfinite(source[eqIdx]))
-                        	DUNE_THROW(NumericalProblem, "Calculated non-finite source");
+                            DUNE_THROW(NumericalProblem, "Calculated non-finite source");
             }
         }
         Valgrind::CheckDefined(source);
diff --git a/dumux/implicit/mpnc/mpnclocalresidual.hh b/dumux/implicit/mpnc/mpnclocalresidual.hh
index ac9afb34069052cc15765464fd0236c70437c392..43cfb6dac84f781126bf6e25086d82e5083c8279 100644
--- a/dumux/implicit/mpnc/mpnclocalresidual.hh
+++ b/dumux/implicit/mpnc/mpnclocalresidual.hh
@@ -230,10 +230,10 @@ public:
     /*!
      * \brief Evaluate the local residual.
      *
-     * 		\param element The finite element
-     * 		\param fvGeometry The finite-volume geometry in the fully implicit scheme
-     * 		\param prevElemVolVars The element volume variables of the previous timestep
-     * 		\param curElemVolVars The element volume variables of the current timestep
+     *      \param element The finite element
+     *      \param fvGeometry The finite-volume geometry in the fully implicit scheme
+     *      \param prevElemVolVars The element volume variables of the previous timestep
+     *      \param curElemVolVars The element volume variables of the current timestep
      *      \param bcType The types of the boundary conditions for all vertices of the element
      */
     void eval(const Element &element,
@@ -273,8 +273,8 @@ public:
      * Sets the Dirichlet conditions in a strong sense, in contrast to
      * the general handling in CCLocalResidual.
      *
-     * 		\param isIt
-     * 		\param bcTypes
+     *      \param isIt
+     *      \param bcTypes
      */
     template <class IntersectionIterator>
     void evalDirichletSegment_(const IntersectionIterator &isIt,
diff --git a/dumux/implicit/mpnc/mpncmodelkinetic.hh b/dumux/implicit/mpnc/mpncmodelkinetic.hh
index f57f4d89571681a16912c2a5b6f040b50814cd5f..398029e46784c3330d8f9a01d92c03ea263e66bc 100644
--- a/dumux/implicit/mpnc/mpncmodelkinetic.hh
+++ b/dumux/implicit/mpnc/mpncmodelkinetic.hh
@@ -238,26 +238,26 @@ public:
 //                    }
 //                }
 //
-//				// interfacial area check (interfacial area between fluid as well as solid phases)
-//				for(int phaseIdxI=0; phaseIdxI<numPhases+1; phaseIdxI++){
-//					const Scalar eps = 1e-6 ;
-//					for (int phaseIdxII=0; phaseIdxII< numPhases+1; ++ phaseIdxII){
-//						if (phaseIdxI == phaseIdxII)
-//							continue;
-//						assert(numEnergyEqs == 3) ; // otherwise this ia call does not make sense
-//						const Scalar ia = elemVolVars[scvIdx].interfacialArea(phaseIdxI, phaseIdxII);
-//						if (not std::isfinite(ia) or ia < 0.-eps ) {
-//							message <<"\nUnphysical Value in interfacial area: \n";
-//							message << "\tia" <<FluidSystem::phaseName(phaseIdxI)
-//											 <<FluidSystem::phaseName(phaseIdxII)<<"="
-//									<< ia << "\n" ;
-//							message << "\t S[0]=" << fluidState.saturation(0);
-//							message << "\t S[1]=" << fluidState.saturation(1);
-//							message << "\t p[0]=" << fluidState.pressure(0);
-//							message << "\t p[1]=" << fluidState.pressure(1);
-//						}
-//					}
-//				}
+//              // interfacial area check (interfacial area between fluid as well as solid phases)
+//              for(int phaseIdxI=0; phaseIdxI<numPhases+1; phaseIdxI++){
+//                  const Scalar eps = 1e-6 ;
+//                  for (int phaseIdxII=0; phaseIdxII< numPhases+1; ++ phaseIdxII){
+//                      if (phaseIdxI == phaseIdxII)
+//                          continue;
+//                      assert(numEnergyEqs == 3) ; // otherwise this ia call does not make sense
+//                      const Scalar ia = elemVolVars[scvIdx].interfacialArea(phaseIdxI, phaseIdxII);
+//                      if (not std::isfinite(ia) or ia < 0.-eps ) {
+//                          message <<"\nUnphysical Value in interfacial area: \n";
+//                          message << "\tia" <<FluidSystem::phaseName(phaseIdxI)
+//                                           <<FluidSystem::phaseName(phaseIdxII)<<"="
+//                                  << ia << "\n" ;
+//                          message << "\t S[0]=" << fluidState.saturation(0);
+//                          message << "\t S[1]=" << fluidState.saturation(1);
+//                          message << "\t p[0]=" << fluidState.pressure(0);
+//                          message << "\t p[1]=" << fluidState.pressure(1);
+//                      }
+//                  }
+//              }
 //
 //                // General Check
 //                for(int phaseIdx=0; phaseIdx<numPhases; phaseIdx++){
@@ -270,17 +270,17 @@ public:
 //                    }
 //                }
 //
-//				// velocity Check
+//              // velocity Check
 //                const unsigned int globalVertexIdx = this->problem_().vertexMapper().map(element, scvIdx, dim);
-//				for(int phaseIdx=0; phaseIdx<numPhases; phaseIdx++){
-//					const Scalar eps = 1e-6 ;
-//					const Scalar velocityTest = volumeDarcyMagVelocity(phaseIdx, globalVertexIdx);
-//					if (not std::isfinite(velocityTest) ){
-//						message <<"\nUnphysical Value in Velocity: \n";
-//						message << "\tv" <<"_"<<FluidSystem::phaseName(phaseIdx)<<"=" << std::scientific
-//						<< velocityTest << std::fixed << "\n";
-//					}
-//				}
+//              for(int phaseIdx=0; phaseIdx<numPhases; phaseIdx++){
+//                  const Scalar eps = 1e-6 ;
+//                  const Scalar velocityTest = volumeDarcyMagVelocity(phaseIdx, globalVertexIdx);
+//                  if (not std::isfinite(velocityTest) ){
+//                      message <<"\nUnphysical Value in Velocity: \n";
+//                      message << "\tv" <<"_"<<FluidSystem::phaseName(phaseIdx)<<"=" << std::scientific
+//                      << velocityTest << std::fixed << "\n";
+//                  }
+//              }
 //
 //
 //                // Some check wrote into the error-message, add some additional information and throw
diff --git a/dumux/implicit/mpnc/mpncvolumevariables.hh b/dumux/implicit/mpnc/mpncvolumevariables.hh
index 61f138d95dec3cc8e0f10b5da65de26d7d5603ac..abaee9efa3dba7bc4cde979951f6c7e1e112a3e0 100644
--- a/dumux/implicit/mpnc/mpncvolumevariables.hh
+++ b/dumux/implicit/mpnc/mpncvolumevariables.hh
@@ -103,12 +103,12 @@ public:
     /*!
      * \copydoc ImplicitVolumeVariables::update
      *
-     * 		\param priVars The primary Variables
-     * 		\param problem The Problem
-     * 		\param element The finite element
-     * 		\param fvGeometry The finite-volume geometry in the fully implicit scheme
-     * 		\param scvIdx The index of the sub-control volume
-     * 		\param isOldSol Specifies whether this is the previous solution or the current one
+     *      \param priVars The primary Variables
+     *      \param problem The Problem
+     *      \param element The finite element
+     *      \param fvGeometry The finite-volume geometry in the fully implicit scheme
+     *      \param scvIdx The index of the sub-control volume
+     *      \param isOldSol Specifies whether this is the previous solution or the current one
      */
     void update(const PrimaryVariables &priVars,
                 const Problem &problem,
@@ -315,7 +315,7 @@ public:
     /*!
      * \brief Returns the value of the NCP-function for a phase.
      *
-     * 		\param phaseIdx The local index of the phases
+     *      \param phaseIdx The local index of the phases
      */
     Scalar phaseNcp(const unsigned int phaseIdx) const
     {
diff --git a/dumux/implicit/mpnc/mpncvolumevariablesia.hh b/dumux/implicit/mpnc/mpncvolumevariablesia.hh
index 3bfeb1f3fa5515bceb3553560819f223c55aa6b3..c7bd81ab1e792e82250b9353685464227c4747f1 100644
--- a/dumux/implicit/mpnc/mpncvolumevariablesia.hh
+++ b/dumux/implicit/mpnc/mpncvolumevariablesia.hh
@@ -64,14 +64,14 @@ public:
     /*!
      * \brief Updates the volume specific interfacial area [m^2 / m^3] between the phases.
      *
-     * 		\param volVars The volume variables
-     * 		\param fluidState Container for all the secondary variables concerning the fluids
-     * 		\param paramCache Container for cache parameters
-     * 		\param priVars The primary Variables
-     * 		\param problem The problem
-     * 		\param element The finite element
-     * 		\param fvGeometry The finite-volume geometry in the fully implicit scheme
-     * 		\param scvIdx The index of the sub-control volumete element
+     *      \param volVars The volume variables
+     *      \param fluidState Container for all the secondary variables concerning the fluids
+     *      \param paramCache Container for cache parameters
+     *      \param priVars The primary Variables
+     *      \param problem The problem
+     *      \param element The finite element
+     *      \param fvGeometry The finite-volume geometry in the fully implicit scheme
+     *      \param scvIdx The index of the sub-control volumete element
      *
      */
     void update(const VolumeVariables & volVars,
diff --git a/dumux/implicit/mpnc/mpncvolumevariablesiakinetic.hh b/dumux/implicit/mpnc/mpncvolumevariablesiakinetic.hh
index 621f0f04c8303338cb1402cf9ca527d0b0384722..7032a01e2786ca4e7188d095b02d783c522e3fff 100644
--- a/dumux/implicit/mpnc/mpncvolumevariablesiakinetic.hh
+++ b/dumux/implicit/mpnc/mpncvolumevariablesiakinetic.hh
@@ -389,10 +389,10 @@ public:
             const Scalar density              = fluidState.density(phaseIdx);
             const Scalar kinematicViscosity   = dynamicViscosity / density;
             const Scalar heatCapacity         = FluidSystem::heatCapacity(fluidState,
-								                                          paramCache,
-								                                          phaseIdx);
+                                                                          paramCache,
+                                                                          phaseIdx);
             const Scalar thermalConductivity  = FluidSystem::thermalConductivity(fluidState,
-									                                       paramCache,
+                                                                           paramCache,
                                                                            phaseIdx);
 
             const Scalar porosity = problem.spatialParams().porosity(element,
diff --git a/dumux/implicit/mpnc/velomodelnewtoncontroller.hh b/dumux/implicit/mpnc/velomodelnewtoncontroller.hh
index 4653d8a67616049e2923353ea451512a199f5a32..46e15c703c2c1335ad8a38de1189368edcb47019 100644
--- a/dumux/implicit/mpnc/velomodelnewtoncontroller.hh
+++ b/dumux/implicit/mpnc/velomodelnewtoncontroller.hh
@@ -58,8 +58,8 @@ public:
 
         // Averages the face velocities of a vertex. Implemented in the model.
         // The velocities are stored in the model.
-    	if(velocityAveragingInModel)
-    		this->problem_().model().calcVelocityAverage();
+        if(velocityAveragingInModel)
+            this->problem_().model().calcVelocityAverage();
     }
 
     void newtonUpdate(SolutionVector &uCurrentIter,
diff --git a/dumux/io/cpgridcreator.hh b/dumux/io/cpgridcreator.hh
index 2b1a7d1463f74efd131f832bafe894797a3e9faa..d05dc637834ceee4d319f3c72e4acb5532709a53 100644
--- a/dumux/io/cpgridcreator.hh
+++ b/dumux/io/cpgridcreator.hh
@@ -83,7 +83,7 @@ public:
 
     /*!
      * \brief Returns a reference to the input deck.
-     * 
+     *
      * The input deck can be used to read parameters like porosity/permeability.
      */
     static Opm::DeckConstPtr &deck()
diff --git a/dumux/io/vtkmultiwriter.hh b/dumux/io/vtkmultiwriter.hh
index 187fe38d94df08b1e4c2faee3de7378cec8e7ca6..8c576beb03b6312bc37075d31ad0c85db25002ef 100644
--- a/dumux/io/vtkmultiwriter.hh
+++ b/dumux/io/vtkmultiwriter.hh
@@ -372,7 +372,7 @@ private:
         if (commSize_ > 1) {
             std::ostringstream oss;
             oss << "s" << std::setw(4) << std::setfill('0') << commSize_
-     	        << "-p" << std::setw(4) << std::setfill('0') << rank
+                << "-p" << std::setw(4) << std::setfill('0') << rank
                 << "-" << simName_ << "-"
                 << std::setw(5) << curWriterNum_;
             return oss.str();
diff --git a/dumux/linear/pardisobackend.hh b/dumux/linear/pardisobackend.hh
index 9ee5d39d0e4541607e2c1539697c16e41ee648f0..a9379229c5d660c017ee337c6df955f9e1276d8d 100644
--- a/dumux/linear/pardisobackend.hh
+++ b/dumux/linear/pardisobackend.hh
@@ -384,7 +384,7 @@ public:
     {}
 
     /*! \brief Apply one step of the preconditioner to the system \f$ A(v)=d \f$.
-	 *
+     *
      * On entry \f$ v=0 \f$ and \f$ d=b-A(x) \f$ (although this might not be
      * computed in that way. On exit v contains the update, i.e
      * one step computes \f$ v = M^{-1} d \f$ where \f$ M \f$ is the
@@ -422,11 +422,11 @@ public:
     }
 
     /*! \brief Clean up.
-	 *
+     *
      * This method is called after the last apply call for the
      * linear system to be solved. Memory may be deallocated safely
      * here. x is the solution of the linear equation.
-	 *
+     *
      * \param b The right hand side of the equation.
     */
     virtual void post (Vector& b)
diff --git a/dumux/material/binarycoefficients/brine_air.hh b/dumux/material/binarycoefficients/brine_air.hh
index e39e27ff5ca6b3acbeeab9cf188a2cce005beb35..51b5ffec1a819a79e9bc91fe747efb350f633877 100644
--- a/dumux/material/binarycoefficients/brine_air.hh
+++ b/dumux/material/binarycoefficients/brine_air.hh
@@ -55,12 +55,12 @@ public:
     static Scalar gasDiffCoeff(Scalar temperature, Scalar pressure) {
         //Diffusion coefficient of water in the Air phase
         const Scalar Theta=1.8;
-	    const Scalar Daw=2.13e-5;  /* reference value */
-	    const Scalar pg0=1.e5;     /* reference pressure */
-	    const Scalar T0=273.15;    /* reference temperature */
-	    Scalar Dgaw;
-	    Dgaw=Daw*(pg0/pressure)*pow((temperature/T0),Theta);
-	    return Dgaw;
+        const Scalar Daw=2.13e-5;  /* reference value */
+        const Scalar pg0=1.e5;     /* reference pressure */
+        const Scalar T0=273.15;    /* reference temperature */
+        Scalar Dgaw;
+        Dgaw=Daw*(pg0/pressure)*pow((temperature/T0),Theta);
+        return Dgaw;
     }
     ;
 
@@ -76,7 +76,7 @@ public:
         * Himmelblau by the temperature.
         * \param temperature The temperature \f$\mathrm{[K]}\f$
         * \param pressure The pressure \f$\mathrm{[Pa]}\f$
-        * 
+        *
         * See:
         *
         * R. Reid et al.: "The properties of Gases and Liquids", 4th edition,
@@ -112,13 +112,13 @@ public:
      * \param xlNaCl the xlNaCl
      */
     static void calculateMoleFractions(const Scalar temperature,
-                                       const Scalar pg, 
+                                       const Scalar pg,
                                        const Scalar XlNaCl,
                                        const int knownPhaseIdx,
-                                       Scalar &xlAir, 
+                                       Scalar &xlAir,
                                        Scalar &ygH2O,
                                      Scalar &xlNaCl) {
-    	DUNE_THROW(Dune::InvalidStateException, "Function: " << "calculateMoleFractions" << " is invalid.");
+        DUNE_THROW(Dune::InvalidStateException, "Function: " << "calculateMoleFractions" << " is invalid.");
 //        Scalar A = computeA_(temperature, pg);
 //
 //        /* XlNaCl: conversion from mass fraction to mol fraction */
@@ -140,7 +140,7 @@ public:
 //        // with the mutual solubility function
 //        if (knownPhaseIdx == lPhaseIdx) {
 ////            ygH2O = A * (1 - xlAir - xlNaCl);
-//        	DUNE_THROW(Dune::InvalidStateException, "phase index: " << "lPhaseIdx" << " is invalid.");
+//          DUNE_THROW(Dune::InvalidStateException, "phase index: " << "lPhaseIdx" << " is invalid.");
 //
 //        }
 //
@@ -150,7 +150,7 @@ public:
 //        if (knownPhaseIdx == gPhaseIdx) {
 //            //y_H2o = fluidstate.
 ////            xlAir = 1 - xlNaCl - ygH2O / A;
-//        	DUNE_THROW(Dune::InvalidStateException, "phase index: " << "gPhaseIdx" << " is invalid.");
+//          DUNE_THROW(Dune::InvalidStateException, "phase index: " << "gPhaseIdx" << " is invalid.");
 //        }
     }
 
@@ -227,7 +227,7 @@ private:
      */
     static Scalar massTomoleFrac_(Scalar XlNaCl) {
 
-    	DUNE_THROW(Dune::InvalidStateException, "Function: " << "massTomoleFrac_" << " is invalid.");
+        DUNE_THROW(Dune::InvalidStateException, "Function: " << "massTomoleFrac_" << " is invalid.");
 
 //        const Scalar Mw = H2O::molarMass(); /* molecular weight of water [kg/mol] */
 //        const Scalar Ms = 58.8e-3; /* molecular weight of NaCl  [kg/mol] */
diff --git a/dumux/material/binarycoefficients/h2o_o2.hh b/dumux/material/binarycoefficients/h2o_o2.hh
index db77b85969725420e1a4982ee11007ce5bf2c475..36ba3fa5d6d8a40b6c09b64b2eab877ddd70c843 100644
--- a/dumux/material/binarycoefficients/h2o_o2.hh
+++ b/dumux/material/binarycoefficients/h2o_o2.hh
@@ -56,7 +56,7 @@ public:
 
         return henryIAPWS(E, F, G, H, temperature);
     };
-    
+
     /*!
      * \brief Binary diffusion coefficent \f$\mathrm{[m^2/s]}\f$ for molecular water and nitrogen.
      *
diff --git a/dumux/material/components/brinevarsalinity.hh b/dumux/material/components/brinevarsalinity.hh
old mode 100755
new mode 100644
index 8bbff134b067e24c6d42f6bd314f594a792f1516..5373d06a4f7f78960388883668085f9e05edbd49
--- a/dumux/material/components/brinevarsalinity.hh
+++ b/dumux/material/components/brinevarsalinity.hh
@@ -125,7 +125,7 @@ public:
      * Equations given in:    - Palliser & McKibbin 1997
      *                         - Michaelides 1981
      *                         - Daubert & Danner 1989
-     * 
+     *
      */
     static const Scalar liquidEnthalpy(Scalar T,
                                        Scalar p, Scalar salinity)
@@ -240,7 +240,7 @@ public:
 
         Scalar rhow = H2O::liquidDensity(temperature, pressure);
 
-        	Scalar	density =  rhow +
+            Scalar  density =  rhow +
             1000*salinity*(
                 0.668 +
                 0.44*salinity +
@@ -314,7 +314,7 @@ public:
      * \param temperature temperature of component in \f$\mathrm{[K]}\f$
      * \param pressure pressure of component in \f$\mathrm{[Pa]}\f$
      * \param salinity The mass fraction of salt
-     * 
+     *
      * Equation given in:    - Batzle & Wang (1992)
      *                         - cited by: Bachu & Adams (2002)
      *                           "Equations of State for basin geofluids"
@@ -330,7 +330,7 @@ public:
         Scalar A = (0.42*pow((pow(salinity, 0.8)-0.17), 2) + 0.045)*pow(T_C, 0.8);
         Scalar mu_brine = 0.1 + 0.333*salinity + (1.65+91.9*salinity*salinity*salinity)*exp(-A);
         assert(mu_brine > 0.0);
-        return mu_brine/1000.0; 
+        return mu_brine/1000.0;
     }
 };
 } // end namespace
diff --git a/dumux/material/components/nacl.hh b/dumux/material/components/nacl.hh
index a75df18f37c28dcecffb533eeac9f5365064333e..c10722f427375d41cddb17c71e5d2cab0b14f29d 100644
--- a/dumux/material/components/nacl.hh
+++ b/dumux/material/components/nacl.hh
@@ -46,24 +46,24 @@ public:
      * \brief A human readable name for the NaCl.
      */
     static const char *name()
-    { 
-        return "NaCl"; 
+    {
+        return "NaCl";
     }
 
     /*!
      * \brief The molar mass of NaCl in \f$\mathrm{[kg/mol]}\f$.
      */
     static Scalar molarMass()
-    { 
-        return 58.4428e-3 ; 
+    {
+        return 58.4428e-3 ;
     }
 
     /*!
      * \brief The diffusion Coefficient \f$\mathrm{[m^2/s]}\f$ of NaCl in water.
      */
     static Scalar liquidDiffCoeff(Scalar temperature, Scalar pressure)
-    { 
-        return 2e-9; 
+    {
+        return 2e-9;
     }
 
     /*!
@@ -78,4 +78,4 @@ public:
 } // end namespace
 
 #endif
- 
+
diff --git a/dumux/material/constraintsolvers/compositionfromfugacities2pncmin.hh b/dumux/material/constraintsolvers/compositionfromfugacities2pncmin.hh
index 08af7ecef4d291dbf58e269e1283b2575e67e867..ad1d3cba5cb4ad963db6dcc62a57c7da1e252047 100644
--- a/dumux/material/constraintsolvers/compositionfromfugacities2pncmin.hh
+++ b/dumux/material/constraintsolvers/compositionfromfugacities2pncmin.hh
@@ -40,9 +40,9 @@ namespace Dumux {
 template <class Scalar, class FluidSystem>
 class compositionFromFugacities2pncmin
 {
-    enum { 
+    enum {
             numComponents = FluidSystem::numComponents,
-            numMajorComponents = FluidSystem::numMajorComponents 
+            numMajorComponents = FluidSystem::numMajorComponents
          };
 
     typedef typename FluidSystem::ParameterCache ParameterCache;
@@ -55,8 +55,8 @@ public:
      * \brief Guess an initial value for the composition of the phase.
      * \param fluidState Thermodynamic state of the fluids
      * \param paramCache  Container for cache parameters
-     * \param phaseIdx The phase index     
-     * \param phasePresence The presence index of the reference phase 
+     * \param phaseIdx The phase index
+     * \param phasePresence The presence index of the reference phase
      * \param fugVec fugacity vector of the component
      */
     template <class FluidState>
@@ -70,7 +70,7 @@ public:
             return;
 
         // Pure component fugacities
-        for (int i = 0; i < numComponents; ++ i) 
+        for (int i = 0; i < numComponents; ++ i)
         {
             fluidState.setMoleFraction(phaseIdx,i, 1.0/numComponents);
         }
@@ -78,10 +78,10 @@ public:
 
     /*!
      * \brief Calculates the chemical equilibrium from the component
-     *        fugacities in a phase. This constraint solver is developed for drying scenarios where 
-     *        salt component is restricted to liquid phase and still for the sake for equilibrium 
-     *        calculation some residual salt must be considered in the gas phase. In such cases for 
-     *        existence of gas phase only, in the theoretical liquid phase, we set the mole fraction 
+     *        fugacities in a phase. This constraint solver is developed for drying scenarios where
+     *        salt component is restricted to liquid phase and still for the sake for equilibrium
+     *        calculation some residual salt must be considered in the gas phase. In such cases for
+     *        existence of gas phase only, in the theoretical liquid phase, we set the mole fraction
      *        of salt to  1e-10.
      * \param fluidState Thermodynamic state of the fluids
      * \param paramCache  Container for cache parameters
@@ -100,7 +100,7 @@ public:
     {
         // use a much more efficient method in case the phase is an
         // ideal mixture
-        if (FluidSystem::isIdealMixture(phaseIdx)) 
+        if (FluidSystem::isIdealMixture(phaseIdx))
         {
             solveIdealMix_(fluidState, paramCache, phaseIdx, phasePresence, targetFug);
             return;
@@ -108,7 +108,7 @@ public:
         else
             DUNE_THROW(NumericalProblem, "This constraint solver is not tested for non-ideal mixtures: Please refer computefromfugacities.hh for details" );
     }
-    
+
 protected:
     // update the phase composition in case the phase is an ideal
     // mixture, i.e. the component's fugacity coefficients are
@@ -120,7 +120,7 @@ protected:
                                int phasePresence,
                                const ComponentVector &fugacities)
     {
-        for (int i = 0; i < numComponents; ++ i) 
+        for (int i = 0; i < numComponents; ++ i)
         {
             Scalar phi = FluidSystem::fugacityCoefficient(fluidState,
                                                           paramCache,
diff --git a/dumux/material/constraintsolvers/computefromreferencephase2pnc.hh b/dumux/material/constraintsolvers/computefromreferencephase2pnc.hh
index 340ac18c8ef522bc4ae7c7ea8d6db6918040cb1a..ad11b5d320ba1d524a314b8fdc01765dd57ef664 100644
--- a/dumux/material/constraintsolvers/computefromreferencephase2pnc.hh
+++ b/dumux/material/constraintsolvers/computefromreferencephase2pnc.hh
@@ -137,7 +137,7 @@ public:
                                                            refPhaseIdx));
 
         // compute the fugacities of all components in the reference phase
-        for (int compIdx = 0; compIdx < numComponents; ++compIdx) 
+        for (int compIdx = 0; compIdx < numComponents; ++compIdx)
         {
             fluidState.setFugacityCoefficient(refPhaseIdx,
                                               compIdx,
@@ -149,7 +149,7 @@ public:
         }
 
         // compute all quantities for the non-reference phases
-        for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) 
+        for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
         {
             if (phaseIdx == refPhaseIdx)
                 continue; // reference phase is already calculated
diff --git a/dumux/material/constraintsolvers/computefromreferencephase2pncmin.hh b/dumux/material/constraintsolvers/computefromreferencephase2pncmin.hh
index 3b0f150f3beb34fce6901bd06a3bbe58496c2705..ffb6b6de61d1d046e22203148cffc0080887385a 100644
--- a/dumux/material/constraintsolvers/computefromreferencephase2pncmin.hh
+++ b/dumux/material/constraintsolvers/computefromreferencephase2pncmin.hh
@@ -101,7 +101,7 @@ public:
      * \param fluidState Thermodynamic state of the fluids
      * \param paramCache  Container for cache parameters
      * \param refPhaseIdx The phase index of the reference phase
-     * \param phasePresence The presence index of the reference phase 
+     * \param phasePresence The presence index of the reference phase
      * \param setViscosity Specify whether the dynamic viscosity of
      *                     each phase should also be set.
      * \param setEnthalpy Specify whether the specific
diff --git a/dumux/material/constraintsolvers/miscible2pnccomposition.hh b/dumux/material/constraintsolvers/miscible2pnccomposition.hh
index 0dbca54ee9998db6120b46e1cfcb7ddae114a9a2..7a8567068412c1af242d8474e53c3a35b2510b55 100644
--- a/dumux/material/constraintsolvers/miscible2pnccomposition.hh
+++ b/dumux/material/constraintsolvers/miscible2pnccomposition.hh
@@ -109,7 +109,7 @@ public:
         Dune::FieldVector<Scalar, numComponents-numMajorComponents> xKnown(0.0);
         for (int knownCompIdx = 0; knownCompIdx < numComponents-numMajorComponents; ++knownCompIdx)
         {
-        	xKnown[knownCompIdx] = fluidState.moleFraction(knownPhaseIdx, knownCompIdx + numMajorComponents);
+            xKnown[knownCompIdx] = fluidState.moleFraction(knownPhaseIdx, knownCompIdx + numMajorComponents);
         }
 
         // compute all fugacity coefficients
@@ -148,10 +148,10 @@ public:
         //Components, of which the molefractions are known, set to molefraction(knownCompIdx)=xKnown
         for(int knownCompIdx = 0; knownCompIdx < numComponents-numMajorComponents; ++knownCompIdx)
         {
-        	int rowIdx = numComponents + numPhases + knownCompIdx;
-        	int colIdx = knownPhaseIdx*numComponents + knownCompIdx + numMajorComponents;
-        	M[rowIdx][colIdx] = 1.0;
-        	b[rowIdx] = xKnown[knownCompIdx];
+            int rowIdx = numComponents + numPhases + knownCompIdx;
+            int colIdx = knownPhaseIdx*numComponents + knownCompIdx + numMajorComponents;
+            M[rowIdx][colIdx] = 1.0;
+            b[rowIdx] = xKnown[knownCompIdx];
         }
 
         // assemble the equations expressing the fact that the
@@ -179,25 +179,25 @@ public:
         //prevents matrix meeting dune's singularity criteria
         for (int compIdx = 0; compIdx < numComponents; compIdx++)
         {
-        	if (compIdx < numMajorComponents)
-        	{
-        		for (int colIdx = 0; colIdx < numPhases*numComponents; colIdx++)
-        		{
-        			//Multiply row of main component (Raoult's Law) with 10e-5 (order of magn. of pressure)
-        			M[compIdx][colIdx] *= 10e-5;
-        		}
-        	} else {
-        		for (int colIdx = 0; colIdx < numPhases*numComponents; colIdx++)
-        		{
-        			//Multiply row of sec. components (Henry's Law) with 10e-9 (order of magn. of Henry constant)
-        			M[compIdx][colIdx] *= 10e-9;
-        		}
-        	}
+            if (compIdx < numMajorComponents)
+            {
+                for (int colIdx = 0; colIdx < numPhases*numComponents; colIdx++)
+                {
+                    //Multiply row of main component (Raoult's Law) with 10e-5 (order of magn. of pressure)
+                    M[compIdx][colIdx] *= 10e-5;
+                }
+            } else {
+                for (int colIdx = 0; colIdx < numPhases*numComponents; colIdx++)
+                {
+                    //Multiply row of sec. components (Henry's Law) with 10e-9 (order of magn. of Henry constant)
+                    M[compIdx][colIdx] *= 10e-9;
+                }
+            }
         }
 
         // solve for all mole fractions
         M.solve(x, b);
-        
+
         // set all mole fractions and the the additional quantities in
         // the fluid state
         for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
diff --git a/dumux/material/fluidmatrixinteractions/2p/philtophoblaw.hh b/dumux/material/fluidmatrixinteractions/2p/philtophoblaw.hh
index cb6537419c978f0c69c68d2b94baf3d4db668d7c..1df1877531a01a0189b4e97d4e86c625dd064a46 100644
--- a/dumux/material/fluidmatrixinteractions/2p/philtophoblaw.hh
+++ b/dumux/material/fluidmatrixinteractions/2p/philtophoblaw.hh
@@ -22,7 +22,7 @@
  * \brief This material law takes a material law defined for effective
  *        saturations and converts it to a material law defined on
  *        absolute saturations. It is valid for hydrophobic materials and is
- *        called with the non-wetting phase saturation and then calls the 
+ *        called with the non-wetting phase saturation and then calls the
  *        material law defined for the wetting phase saturation.
  */
 #ifndef DUMUX_PHIL_TO_PHOB_LAW_HH
@@ -53,15 +53,15 @@ namespace Dumux
  *        saturation from the model. Here, the wetting saturation is calculated from it: \f$\mathrm{S_w = 1 - S_n}\f$.
  *        Then the effective wetting saturations are calculated and handed to the material law.
  *
- *        The following definition shows the pc-Sw-relationship for the case in which the x-axis 
+ *        The following definition shows the pc-Sw-relationship for the case in which the x-axis
  *        is \f$\mathrm{S_w (=S_{water})}\f$ and the y-axis is \f$\mathrm{p_c = p_{non_wetting}-p_{wetting}}\f$.
  *        \image html pc_Sw_1.png
  *
  *        But DumuX actually uses the following curve, because \f$\mathrm{p_c = p_{non-water}-p_{water}}\f$ for the y-axis is used.
  *        \image html pc_Sw_2.png
  *
- *        This approach makes sure that in the "material laws" only effective wetting saturations are considered, 
- *        which makes sense, as these laws only deal with effective saturations. This also allows for changing 
+ *        This approach makes sure that in the "material laws" only effective wetting saturations are considered,
+ *        which makes sense, as these laws only deal with effective saturations. This also allows for changing
  *        the calculation of the effective saturations easily, as this is subject of discussion / may be problem specific.
  *
  *        Additionally, handing over effective saturations to the "material laws" in stead of them calculating effective
@@ -74,8 +74,8 @@ namespace Dumux
  *        - the definition of the material law in the spatial parameters is not really intuitive, but using it is:
  *          Hand in values, get back values, do not deal with conversion.
  *
- *        CAREFULL: here w and n still stand for wetting and non-wetting. In the hydrophobic model w stands for 
- *        water (non-wetting) and n for non-water (wetting). Hence, Swr is the wetting phase (gas) residual 
+ *        CAREFULL: here w and n still stand for wetting and non-wetting. In the hydrophobic model w stands for
+ *        water (non-wetting) and n for non-water (wetting). Hence, Swr is the wetting phase (gas) residual
  *        saturation and needs to be set accordingly in the spatial parameters.
  */
 
@@ -92,13 +92,13 @@ public:
      * \brief The capillary pressure-saturation curve.
      *
      *
-     * \param sn Absolute saturation of the non-wetting phase \f$\mathrm{\overline{S}_n}\f$. It is converted to the 
-     *                  effective saturation of the wetting phase and then handed over to the material law actually 
+     * \param sn Absolute saturation of the non-wetting phase \f$\mathrm{\overline{S}_n}\f$. It is converted to the
+     *                  effective saturation of the wetting phase and then handed over to the material law actually
      *                  used for calculation.
      * \param params A container object that is populated with the appropriate coefficients for the respective law.
-     *                  Therefore, in the (problem specific) spatialParameters  first, the material law is chosen, 
+     *                  Therefore, in the (problem specific) spatialParameters  first, the material law is chosen,
      *                  and then the params container is constructed accordingly. Afterwards the values are set there, too.
-     * \return Capillary pressure \f$\mathrm{p_c}\f$ in \f$\mathrm{[Pa]}\f$  calculated by specific constitutive relation 
+     * \return Capillary pressure \f$\mathrm{p_c}\f$ in \f$\mathrm{[Pa]}\f$  calculated by specific constitutive relation
      *                  (EffLaw e.g. Brooks & Corey, van Genuchten, linear...)*-1 to account for hydrophobic material.
      *
      */
@@ -110,11 +110,11 @@ public:
     /*!
      * \brief The saturation-capillary pressure curve.
      *
-     * \param pc Capillary pressure \f$\mathrm{p_c}\f$ in \f$\mathrm{[Pa]}\f$ 
+     * \param pc Capillary pressure \f$\mathrm{p_c}\f$ in \f$\mathrm{[Pa]}\f$
      * \param params A container object that is populated with the appropriate coefficients for the respective law.
-     *                  Therefore, in the (problem specific) spatialParameters  first, the material law is chosen, 
+     *                  Therefore, in the (problem specific) spatialParameters  first, the material law is chosen,
      *                  and then the params container is constructed accordingly. Afterwards the values are set there, too.
-     *\return Absolute wetting phase saturation calculated as inverse of 
+     *\return Absolute wetting phase saturation calculated as inverse of
      *                  (EffLaw e.g. Brooks & Corey, van Genuchten, linear...) constitutive relation.
      *
      * \return The absolute saturation of the wetting phase \f$\mathrm{S_w}\f$
@@ -135,9 +135,9 @@ public:
      }\f$
      * \param sw Absolute saturation of the wetting phase \f$\mathrm{\overline{S}_w}\f$.
      * \param params A container object that is populated with the appropriate coefficients for the respective law.
-     *                  Therefore, in the (problem specific) spatialParameters  first, the material law is chosen, and 
+     *                  Therefore, in the (problem specific) spatialParameters  first, the material law is chosen, and
      *                  then the params container is constructed accordingly. Afterwards the values are set there, too.
-     * \return          Partial derivative of \f$\mathrm{p_c}\f$ w.r.t. effective saturation according to EffLaw 
+     * \return          Partial derivative of \f$\mathrm{p_c}\f$ w.r.t. effective saturation according to EffLaw
      *                  e.g. Brooks & Corey, van Genuchten, linear... .
      */
     static Scalar dpc_dsw(const Params &params, Scalar sw)
@@ -156,7 +156,7 @@ public:
     }\f$
      *
      *
-     * \param pc Capillary pressure \f$\mathrm{p_c}\f$ in \f$\mathrm{[Pa]}\f$ 
+     * \param pc Capillary pressure \f$\mathrm{p_c}\f$ in \f$\mathrm{[Pa]}\f$
      * \param params A container object that is populated with the appropriate coefficients for the respective law.
      *                  Therefore, in the (problem specific) spatialParameters  first, the material law is chosen, and then the params container
      *                  is constructed accordingly. Afterwards the values are set there, too.
diff --git a/dumux/material/fluidmatrixinteractions/2p/regularizedlinearmaterial.hh b/dumux/material/fluidmatrixinteractions/2p/regularizedlinearmaterial.hh
index d0b7ea48c037ac741d64a96ffbf72d5746379006..8c8741d9f88fb04b705e518873f6cbddc2eb5f48 100644
--- a/dumux/material/fluidmatrixinteractions/2p/regularizedlinearmaterial.hh
+++ b/dumux/material/fluidmatrixinteractions/2p/regularizedlinearmaterial.hh
@@ -91,7 +91,7 @@ public:
      S_w = 1 - \frac{p_C - p_{C,entry}}{p_{C,max} - p_{C,entry}}
      }\f$
      *
-     * \param pc Capillary pressure \f$\mathrm{p_C}\f$ in \f$\mathrm{[Pa]}\f$ 
+     * \param pc Capillary pressure \f$\mathrm{p_C}\f$ in \f$\mathrm{[Pa]}\f$
      * \param params A container object that is populated with the appropriate coefficients for the respective law.
      *                  Therefore, in the (problem specific) spatialParameters  first, the material law is chosen, and then the params container
      *                  is constructed accordingly. Afterwards the values are set there, too.
@@ -125,7 +125,7 @@ public:
      * \brief Returns the partial derivative of the effective
      *        saturation to the capillary pressure.
      *
-     * \param pc Capillary pressure \f$\mathrm{p_C}\f$ in \f$\mathrm{[Pa]}\f$ 
+     * \param pc Capillary pressure \f$\mathrm{p_C}\f$ in \f$\mathrm{[Pa]}\f$
      * \param params A container object that is populated with the appropriate coefficients for the respective law.
      *                  Therefore, in the (problem specific) spatialParameters  first, the material law is chosen, and then the params container
      *                  is constructed accordingly. Afterwards the values are set there, too.
diff --git a/dumux/material/fluidmatrixinteractions/2p/thermalconductivityjohansen.hh b/dumux/material/fluidmatrixinteractions/2p/thermalconductivityjohansen.hh
index e23fb173b5754912b71bc1a32e9c41f96105894f..99939985aa14afa45e735ead8587988133b2190a 100644
--- a/dumux/material/fluidmatrixinteractions/2p/thermalconductivityjohansen.hh
+++ b/dumux/material/fluidmatrixinteractions/2p/thermalconductivityjohansen.hh
@@ -109,11 +109,11 @@ public:
      * \brief Returns the effective thermal conductivity \f$\mathrm{[W/(m K)]}\f$ after Johansen (1975).
      *
      * \param Sw The saturation of the wetting phase
-     * \param lambdaW The thermal conductivity of the wetting phase in \f$\mathrm{[W/(m K)]}\f$ 
-     * \param lambdaN The thermal conductivity of the non-wetting phase in \f$\mathrm{[W/(m K)]}\f$ 
-     * \param lambdaSolid The thermal conductivity of the solid phase in \f$\mathrm{[W/(m K)]}\f$ 
+     * \param lambdaW The thermal conductivity of the wetting phase in \f$\mathrm{[W/(m K)]}\f$
+     * \param lambdaN The thermal conductivity of the non-wetting phase in \f$\mathrm{[W/(m K)]}\f$
+     * \param lambdaSolid The thermal conductivity of the solid phase in \f$\mathrm{[W/(m K)]}\f$
      * \param porosity The porosity
-     * \param rhoSolid The density of solid phase in \f$\mathrm{[kg/m^3]}\f$ 
+     * \param rhoSolid The density of solid phase in \f$\mathrm{[kg/m^3]}\f$
      *
      * \return Effective thermal conductivity \f$\mathrm{[W/(m K)]}\f$ after Johansen (1975)
      */
diff --git a/dumux/material/fluidmatrixinteractions/2p/thermalconductivitysimplefluidlumping.hh b/dumux/material/fluidmatrixinteractions/2p/thermalconductivitysimplefluidlumping.hh
index 02b20c3222e0ac719240b3733993adc6391564f3..140e1b358798631ffaa2a5b7d85bd7d80e94d375 100644
--- a/dumux/material/fluidmatrixinteractions/2p/thermalconductivitysimplefluidlumping.hh
+++ b/dumux/material/fluidmatrixinteractions/2p/thermalconductivitysimplefluidlumping.hh
@@ -74,11 +74,11 @@ public:
      * \brief Returns the effective thermal conductivity \f$\mathrm{[W/(m K)]}\f$.
      *
      * \param sw The saturation of the wetting phase
-     * \param lambdaW The thermal conductivity of the wetting phase in \f$\mathrm{[W/(m K)]}\f$ 
-     * \param lambdaN The thermal conductivity of the non-wetting phase in \f$\mathrm{[W/(m K)]}\f$ 
-     * \param lambdaSolid The thermal conductivity of the solid phase in \f$\mathrm{[W/(m K)]}\f$ 
+     * \param lambdaW The thermal conductivity of the wetting phase in \f$\mathrm{[W/(m K)]}\f$
+     * \param lambdaN The thermal conductivity of the non-wetting phase in \f$\mathrm{[W/(m K)]}\f$
+     * \param lambdaSolid The thermal conductivity of the solid phase in \f$\mathrm{[W/(m K)]}\f$
      * \param porosity The porosity
-     * \param rhoSolid The density of the solid phase in \f$\mathrm{[kg/m^3]}\f$ 
+     * \param rhoSolid The density of the solid phase in \f$\mathrm{[kg/m^3]}\f$
      *
      * \return Effective thermal conductivity of the fluid phases
      */
diff --git a/dumux/material/fluidmatrixinteractions/2p/thermalconductivitysomerton.hh b/dumux/material/fluidmatrixinteractions/2p/thermalconductivitysomerton.hh
index 95713a8d801d9c352e885c78719b3243ca7e58ce..85a328be2e8aef3f8c69f59b6d9d500abab4067c 100644
--- a/dumux/material/fluidmatrixinteractions/2p/thermalconductivitysomerton.hh
+++ b/dumux/material/fluidmatrixinteractions/2p/thermalconductivitysomerton.hh
@@ -103,11 +103,11 @@ public:
      * \brief effective thermal conductivity \f$\mathrm{[W/(m K)]}\f$ after Somerton (1974)
      *
      * \param sw The saturation of the wetting phase
-     * \param lambdaW The thermal conductivity of the wetting phase in \f$\mathrm{[W/(m K)]}\f$ 
-     * \param lambdaN The thermal conductivity of the non-wetting phase in \f$\mathrm{[W/(m K)]}\f$ 
-     * \param lambdaSolid The thermal conductivity of the solid phase in \f$\mathrm{[W/(m K)]}\f$ 
+     * \param lambdaW The thermal conductivity of the wetting phase in \f$\mathrm{[W/(m K)]}\f$
+     * \param lambdaN The thermal conductivity of the non-wetting phase in \f$\mathrm{[W/(m K)]}\f$
+     * \param lambdaSolid The thermal conductivity of the solid phase in \f$\mathrm{[W/(m K)]}\f$
      * \param porosity The porosity
-     * \param rhoSolid The density of solid phase in \f$\mathrm{[kg/m^3]}\f$ 
+     * \param rhoSolid The density of solid phase in \f$\mathrm{[kg/m^3]}\f$
      *
      * \return effective thermal conductivity \f$\mathrm{[W/(m K)]}\f$ after Somerton (1974)
      */
diff --git a/dumux/material/fluidmatrixinteractions/2pia/efftoabslawia.hh b/dumux/material/fluidmatrixinteractions/2pia/efftoabslawia.hh
index a48dbe65daafa00334c78e516f9c83fd98faef6f..ba97f918fa64642e41f0de2252b9bc68df3139f3 100644
--- a/dumux/material/fluidmatrixinteractions/2pia/efftoabslawia.hh
+++ b/dumux/material/fluidmatrixinteractions/2pia/efftoabslawia.hh
@@ -71,7 +71,7 @@ public:
     typedef typename MaterialParams::Scalar Scalar;
 
     /*!
-     * \brief The interfacial area relation 
+     * \brief The interfacial area relation
      * \param sw Absolute saturation of the wetting phase \f$\mathrm{{S}_w}\f$.
      * \param iaParams parameter container for the interfacial area
      * \param pc Capillary pressure in \f$\mathrm{[Pa]}\f$
diff --git a/dumux/material/fluidmatrixinteractions/3p/parkervangen3p.hh b/dumux/material/fluidmatrixinteractions/3p/parkervangen3p.hh
index 8eadca1bdd672c8a5cb1f0a3166d9ce41e43ef3e..e07bd8e1fbd05a076589a83c8deaa6d6a6a253c9 100644
--- a/dumux/material/fluidmatrixinteractions/3p/parkervangen3p.hh
+++ b/dumux/material/fluidmatrixinteractions/3p/parkervangen3p.hh
@@ -62,7 +62,7 @@ public:
         DUNE_THROW(Dune::NotImplemented, "Capillary pressures for three phases is not so simple! Use pcgn, pcnw, and pcgw");
     }
    /*!
-     * \brief The capillary pressure-saturation curve copied from MUFTE/pml/constrel3p3cni.c 
+     * \brief The capillary pressure-saturation curve copied from MUFTE/pml/constrel3p3cni.c
      * \param params Array of parameters
      * \param sw wetting phase saturation or sum of wetting phase saturations
      *
@@ -112,7 +112,7 @@ public:
         }
     }
   /*!
-     * \brief The capillary pressure-saturation curve copied from MUFTE/pml/constrel3p3cni.c 
+     * \brief The capillary pressure-saturation curve copied from MUFTE/pml/constrel3p3cni.c
      * \param params Array of parameters
      * \param sw wetting phase saturation or sum of wetting phase saturations
      */
@@ -161,7 +161,7 @@ public:
         }
     }
     /*!
-     * \brief The capillary pressure-saturation curve copied from MUFTE/pml/constrel3p3cni.c 
+     * \brief The capillary pressure-saturation curve copied from MUFTE/pml/constrel3p3cni.c
      * \param params Array of parameters
      * \param St sum of wetting (liquid) phase saturations
      */
@@ -209,7 +209,7 @@ public:
         }
     }
  /*!
-     * \brief The capillary pressure-saturation curve copied from MUFTE/pml/constrel3p3cni.c 
+     * \brief The capillary pressure-saturation curve copied from MUFTE/pml/constrel3p3cni.c
      * \param params Array of parameters
      * \param sn Non-wetting liquid saturation
      */
diff --git a/dumux/material/fluidmatrixinteractions/3p/parkervangen3pparams.hh b/dumux/material/fluidmatrixinteractions/3p/parkervangen3pparams.hh
index 22f6c71200f928f75242b383d786914e6ee4ce65..ed542fc0add07558e18dcd46d835166088542c2e 100644
--- a/dumux/material/fluidmatrixinteractions/3p/parkervangen3pparams.hh
+++ b/dumux/material/fluidmatrixinteractions/3p/parkervangen3pparams.hh
@@ -181,7 +181,7 @@ public:
      */
     void setSgr(Scalar input)
     { sgr_ = input; }
-     
+
     /*!
      * \brief Set the residual gas saturation.
      */
diff --git a/dumux/material/fluidstates/immisciblefluidstate.hh b/dumux/material/fluidstates/immisciblefluidstate.hh
index 8e6e853d4198c5ca01855b6609de1dd77e5f6856..8f02ae5395d2b0074d9f17c4e4e67ff8d305b61d 100644
--- a/dumux/material/fluidstates/immisciblefluidstate.hh
+++ b/dumux/material/fluidstates/immisciblefluidstate.hh
@@ -130,7 +130,7 @@ public:
 
     /*!
      * \brief The molar volume of a fluid phase \f$\mathrm{[m^3/mol]}\f$
-     * 
+     *
      */
     Scalar molarVolume(int phaseIdx) const
     { return 1/molarDensity(phaseIdx); }
@@ -189,7 +189,7 @@ public:
     { return temperature_; }
 
     /*!
-     * \brief The fugacity of a component 
+     * \brief The fugacity of a component
      *
      * This assumes chemical equilibrium.
      */
diff --git a/dumux/material/fluidsystems/2pimmisciblefluidsystem.hh b/dumux/material/fluidsystems/2pimmisciblefluidsystem.hh
index 5742f4339991583d91fb45ecc8cbdce331e15757..8c4666afc0535be030b5c521f5a66403fd334879 100644
--- a/dumux/material/fluidsystems/2pimmisciblefluidsystem.hh
+++ b/dumux/material/fluidsystems/2pimmisciblefluidsystem.hh
@@ -333,7 +333,7 @@ public:
      * \param fluidState The fluid state of the two-phase model
      * \param phaseIdx Index of the fluid phase
      * \param compIdx index of the component
-     * 
+     *
      * Molecular diffusion of a compoent \f$\mathrm{\kappa}\f$ is caused by a
      * gradient of the chemical potential and follows the law
      *
diff --git a/dumux/material/fluidsystems/basefluidsystem.hh b/dumux/material/fluidsystems/basefluidsystem.hh
index 430899a9ba619d095aef31adf7c8f62ab1b5c3a8..0d9f42397c7ccbad1c2911b3c07b66e9eed02dba 100644
--- a/dumux/material/fluidsystems/basefluidsystem.hh
+++ b/dumux/material/fluidsystems/basefluidsystem.hh
@@ -41,7 +41,7 @@ public:
 
     /*!
      * \brief Calculate the density \f$\mathrm{[kg/m^3]}\f$ of a fluid phase
-     * \param fluidState The fluid state 
+     * \param fluidState The fluid state
      * \param phaseIdx Index of the fluid phase
      * \param paramCache mutable parameters
      */
@@ -56,11 +56,11 @@ public:
     /*!
      * \brief Calculate the fugacity coefficient \f$\mathrm{[Pa]}\f$ of an individual
      *        component in a fluid phase
-     * \param fluidState The fluid state 
+     * \param fluidState The fluid state
      * \param phaseIdx Index of the fluid phase
      * \param paramCache mutable parameters
      * \param compIdx Index of the component
-     * 
+     *
      * The fugacity coefficient \f$\mathrm{\phi_\kappa}\f$ is connected to the
      * fugacity \f$\mathrm{f_\kappa}\f$ and the component's molarity
      * \f$\mathrm{x_\kappa}\f$ by means of the relation
@@ -78,7 +78,7 @@ public:
 
     /*!
      * \brief Calculate the dynamic viscosity of a fluid phase \f$\mathrm{[Pa*s]}\f$
-     * \param fluidState The fluid state 
+     * \param fluidState The fluid state
      * \param phaseIdx Index of the fluid phase
      * \param paramCache mutable parameters
      */
@@ -93,7 +93,7 @@ public:
     /*!
      * \brief Calculate the binary molecular diffusion coefficient for
      *        a component in a fluid phase \f$\mathrm{[mol^2 * s / (kg*m^3)]}\f$
-     * \param fluidState The fluid state 
+     * \param fluidState The fluid state
      * \param phaseIdx Index of the fluid phase
      * \param paramCache mutable parameters
      * \param compIdx Index of the component
@@ -125,11 +125,11 @@ public:
      * \brief Given a phase's composition, temperature and pressure,
      *        return the binary diffusion coefficient \f$\mathrm{[m^2/s]}\f$ for components
      *        \f$\mathrm{i}\f$ and \f$\mathrm{j}\f$ in this phase.
-     * \param fluidState The fluid state 
+     * \param fluidState The fluid state
      * \param phaseIdx Index of the fluid phase
      * \param paramCache mutable parameters
      * \param compIIdx Index of the component i
-     * \param compJIdx Index of the component j 
+     * \param compJIdx Index of the component j
      */
     template <class FluidState>
     static Scalar binaryDiffusionCoefficient(const FluidState &fluidState,
@@ -145,10 +145,10 @@ public:
     /*!
      * \brief Given a phase's composition, temperature, pressure and
      *        density, calculate its specific enthalpy \f$\mathrm{[J/kg]}\f$.
-     * \param fluidState The fluid state 
+     * \param fluidState The fluid state
      * \param phaseIdx Index of the fluid phase
      * \param paramCache mutable parameters
-     * 
+     *
      *  \todo This fluid system neglects the contribution of
      *        gas-molecules in the liquid phase. This contribution is
      *        probably not big. Somebody would have to find out the
@@ -164,10 +164,10 @@ public:
 
     /*!
      * \brief Thermal conductivity of a fluid phase \f$\mathrm{[W/(m K)]}\f$.
-     * \param fluidState The fluid state 
+     * \param fluidState The fluid state
      * \param phaseIdx Index of the fluid phase
      * \param paramCache mutable parameters
-     * 
+     *
      * Use the conductivity of air and water as a first approximation.
      * Source:
      * http://en.wikipedia.org/wiki/List_of_thermal_conductivities
diff --git a/dumux/material/fluidsystems/brineairfluidsystem.hh b/dumux/material/fluidsystems/brineairfluidsystem.hh
index 406f79881343d7068772a7db3b77b50cb7f07c8f..b0bd731941a65edf060916c42bce6db124b5f5c6 100644
--- a/dumux/material/fluidsystems/brineairfluidsystem.hh
+++ b/dumux/material/fluidsystems/brineairfluidsystem.hh
@@ -20,7 +20,7 @@
  * \file
  *
  * \brief A fluid system with a liquid and a gaseous phase and \f$H_2O\f$, \f$Air\f$ and \f$S\f$ (dissolved minerals) as components.
- * 
+ *
  */
 #ifndef DUMUX_BRINE_AIR_SYSTEM_HH
 #define DUMUX_BRINE_AIR_SYSTEM_HH
@@ -45,9 +45,9 @@
 #include <dumux/common/basicproperties.hh>
 #endif
 
-namespace Dumux 
+namespace Dumux
 {
-namespace FluidSystems 
+namespace FluidSystems
 {
 /*!
  * \ingroup Fluidsystems
@@ -79,7 +79,7 @@ namespace FluidSystems
  * An adapter class using Dumux::FluidSystem<TypeTag> is also provided
  * at the end of this file.
  */
-template <class Scalar, 
+template <class Scalar,
           class H2Otype = Dumux::TabulatedComponent<Scalar, Dumux::H2O<Scalar>>,
           bool useComplexRelations=true>
 class BrineAir
@@ -89,7 +89,7 @@ class BrineAir
     typedef BaseFluidSystem <Scalar, ThisType> Base;
 
     typedef Dumux::IdealGas<Scalar> IdealGas;
-    
+
 public:
 
     typedef H2Otype H2O;
@@ -98,10 +98,10 @@ public:
     typedef Dumux::BinaryCoeff::Brine_Air<Scalar, Air> Brine_Air;
     typedef Dumux::BrineVarSalinity<Scalar, H2Otype> Brine;
     typedef Dumux::NaCl<Scalar> NaCl;
-     
+
     // the type of parameter cache objects. this fluid system does not
     typedef Dumux::NullParameterCache ParameterCache;
-    
+
     /****************************************
      * Fluid phase related static parameters
      ****************************************/
@@ -127,12 +127,12 @@ public:
      */
     static const char *phaseName(int phaseIdx)
     {
-    	switch (phaseIdx) {
-		case wPhaseIdx: return "liquid";
-		case nPhaseIdx: return "gas";
-		case sPhaseIdx: return "NaCl";
-		};
-		DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx);
+        switch (phaseIdx) {
+        case wPhaseIdx: return "liquid";
+        case nPhaseIdx: return "gas";
+        case sPhaseIdx: return "NaCl";
+        };
+        DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx);
     }
 
     /*!
@@ -165,8 +165,8 @@ public:
     {
         assert(0 <= phaseIdx && phaseIdx < numPhases);
         // we assume Henry's and Rault's laws for the water phase and
-		// and no interaction between gas molecules of different
-		// components, so all phases are ideal mixtures!
+        // and no interaction between gas molecules of different
+        // components, so all phases are ideal mixtures!
 
         return true;
     }
@@ -182,12 +182,12 @@ public:
      */
     static bool isCompressible(int phaseIdx)
     {
-    	assert(0 <= phaseIdx && phaseIdx < numPhases);
-		// ideal gases are always compressible
-		if (phaseIdx == nPhaseIdx)
-			return true;
-		// the water component decides for the liquid phase...
-		return H2O::liquidIsCompressible();
+        assert(0 <= phaseIdx && phaseIdx < numPhases);
+        // ideal gases are always compressible
+        if (phaseIdx == nPhaseIdx)
+            return true;
+        // the water component decides for the liquid phase...
+        return H2O::liquidIsCompressible();
     }
 
     /*!
@@ -223,13 +223,13 @@ public:
      */
     static const char *componentName(int compIdx)
     {
-    	switch (compIdx)
-		{
-		case H2OIdx: return H2O::name();
-		case AirIdx: return Air::name();
-		case NaClIdx:return "NaCl";
-		};
-		DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << compIdx);
+        switch (compIdx)
+        {
+        case H2OIdx: return H2O::name();
+        case AirIdx: return Air::name();
+        case NaClIdx:return "NaCl";
+        };
+        DUNE_THROW(Dune::InvalidStateException, "Invalid component index " << compIdx);
     }
 
     /*!
@@ -273,7 +273,7 @@ public:
      */
      static Scalar saltSpecificHeatCapacity(int phaseIdx)//Specific heat capacity per unit mole of solid salt phase (J/Kkg)
     {
-	return 36.79/molarMass(phaseIdx);
+    return 36.79/molarMass(phaseIdx);
     }
     /*!
      * \brief Return the molar density of the precipitate \f$\mathrm{[mol/m^3]}\f$.
@@ -282,7 +282,7 @@ public:
      */
     static Scalar precipitateMolarDensity(int phaseIdx)//Density of solid salt phase (mol/m3)
      {
-     	return precipitateDensity(phaseIdx)/molarMass(phaseIdx);
+        return precipitateDensity(phaseIdx)/molarMass(phaseIdx);
      }
 
     /****************************************
@@ -421,7 +421,7 @@ public:
      * \param  fluidState The fluid state
      * \param phaseIdx Index of the phase
      * \param compIdx Index of the component
-     * 
+     *
      * The fugacity coefficient \f$\mathrm{\phi^\kappa_\alpha}\f$ of
      * component \f$\mathrm{\kappa}\f$ in phase \f$\mathrm{\alpha}\f$ is connected to
      * the fugacity \f$\mathrm{f^\kappa_\alpha}\f$ and the component's mole
@@ -451,7 +451,7 @@ public:
         assert(p > 0);
 
         if (phaseIdx == gPhaseIdx)
-        	return 1.0;
+            return 1.0;
 
         else if (phaseIdx == lPhaseIdx)
         {
@@ -506,11 +506,11 @@ public:
             assert(compJIdx == AirIdx || compJIdx == NaClIdx);
             Scalar result = 0.0;
             if(compJIdx == AirIdx)
-            	result = Brine_Air::liquidDiffCoeff(temperature, pressure);
+                result = Brine_Air::liquidDiffCoeff(temperature, pressure);
             else if (compJIdx == NaClIdx)
-            	result = 0.12e-9; //http://webserver.dmt.upm.es/~isidoro/dat1/Mass%20diffusivity%20data.htm
+                result = 0.12e-9; //http://webserver.dmt.upm.es/~isidoro/dat1/Mass%20diffusivity%20data.htm
             else
-            	DUNE_THROW(Dune::NotImplemented, "Binary difussion coefficient : Incorrect compIdx");
+                DUNE_THROW(Dune::NotImplemented, "Binary difussion coefficient : Incorrect compIdx");
             Valgrind::CheckDefined(result);
             return result;
         }
@@ -555,7 +555,7 @@ public:
      */
     using Base::enthalpy;
     template <class FluidState>
-    static Scalar enthalpy(const FluidState &fluidState, 
+    static Scalar enthalpy(const FluidState &fluidState,
                            int phaseIdx)
     {
         assert(0 <= phaseIdx && phaseIdx < numPhases);
@@ -617,8 +617,8 @@ public:
         }
         DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx);
     }
-    
-    
+
+
     /*!
      * \brief Thermal conductivity of a fluid phase \f$\mathrm{[W/(m K)]}\f$.
      *
@@ -685,7 +685,7 @@ private:
                               Scalar pg,
                               Scalar xgH2O)
     {
-    	Scalar pH2O = xgH2O*pg; //Dalton' Law
+        Scalar pH2O = xgH2O*pg; //Dalton' Law
         Scalar pAir = pg - pH2O;
         Scalar gasDensityAir = Air::gasDensity(T, pAir);
         Scalar gasDensityH2O = H2O::gasDensity(T, pH2O);
@@ -794,11 +794,11 @@ private:
 //     Scalar R = 8.314;//[j/K mol] universal gas constant
 //     Scalar pi = (R * T * std::log(1- x))/vw;
 //     if (x > 0.26) // here we have hard coaded the solubility limit for NaCl
-//   	pi = (R * T * std::log(0.74))/vw;
+//      pi = (R * T * std::log(0.74))/vw;
 //     p_s = p_0 * std::exp((pi*vw)/(R*T));// Kelvin's law for reduction in saturation vapor pressure due to osmotic potential
 // #endif
       return p_s;
-    }   
+    }
 };
 
 } // end namespace
diff --git a/dumux/material/fluidsystems/brineco2fluidsystem.hh b/dumux/material/fluidsystems/brineco2fluidsystem.hh
index 201b52564add4884c1dd3cbc70d9ac3b3ed30618..ebb1bc7f52c04b8b0213595991207d72fc098c93 100644
--- a/dumux/material/fluidsystems/brineco2fluidsystem.hh
+++ b/dumux/material/fluidsystems/brineco2fluidsystem.hh
@@ -388,7 +388,7 @@ public:
     /*!
      * \brief Returns the equilibrium concentration of the dissolved component
      *        in a phase.
-     * \param fluidState An arbitrary fluid state 
+     * \param fluidState An arbitrary fluid state
      * \param paramCache Parameter cache
      * \param phaseIdx The index of the fluid phase to consider
      */
@@ -467,7 +467,7 @@ public:
      * \param fluidState An arbitrary fluid state
      * \param phaseIdx The index of the fluid phase to consider
      * \param compIIdx Index of the component i
-     * \param compJIdx Index of the component j 
+     * \param compJIdx Index of the component j
      */
     using Base::binaryDiffusionCoefficient;
     template <class FluidState>
diff --git a/dumux/material/fluidsystems/gasphase.hh b/dumux/material/fluidsystems/gasphase.hh
index c8b80c792a32c33a8e61501beb9346316f1833ad..5ad6d25462c0a4f93943c1324264c06727426054 100644
--- a/dumux/material/fluidsystems/gasphase.hh
+++ b/dumux/material/fluidsystems/gasphase.hh
@@ -67,19 +67,19 @@ public:
     {  return Component::molarMass(); }
 
     /*!
-     * \brief Returns the critical temperature in \f$\mathrm{[K]}\f$ of the component 
+     * \brief Returns the critical temperature in \f$\mathrm{[K]}\f$ of the component
      */
     static Scalar criticalTemperature()
     {  return Component::criticalTemperature(); }
 
     /*!
-     * \brief Returns the critical pressure in \f$\mathrm{[Pa]}\f$ of the component 
+     * \brief Returns the critical pressure in \f$\mathrm{[Pa]}\f$ of the component
      */
     static Scalar criticalPressure()
     {  return Component::criticalPressure(); }
 
     /*!
-     * \brief Returns the temperature in \f$\mathrm{[K]}\f$ at the component's triple point. 
+     * \brief Returns the temperature in \f$\mathrm{[K]}\f$ at the component's triple point.
      */
     static Scalar tripleTemperature()
     {  return Component::tripleTemperature(); }
diff --git a/dumux/material/fluidsystems/h2oairfluidsystem.hh b/dumux/material/fluidsystems/h2oairfluidsystem.hh
index 1f69742379651b1b01ac16cfc301d2560223a24d..4f5e4c9973604d44e80f8a440168c14901bdeedf 100644
--- a/dumux/material/fluidsystems/h2oairfluidsystem.hh
+++ b/dumux/material/fluidsystems/h2oairfluidsystem.hh
@@ -675,7 +675,7 @@ public:
      * \param fluidState An arbitrary fluid state
      * \param phaseIdx The index of the fluid phase to consider
      * \param componentIdx The index of the component to consider
-     * 
+     *
      */
     template <class FluidState>
     static Scalar componentEnthalpy(const FluidState &fluidState,
diff --git a/dumux/material/fluidsystems/h2oairmesitylenefluidsystem.hh b/dumux/material/fluidsystems/h2oairmesitylenefluidsystem.hh
index fcd753fba3409c5198e0ef20c0dc72c3b9e438ca..028d2f55c1611f21501a721c9c12d73d4a892bd9 100644
--- a/dumux/material/fluidsystems/h2oairmesitylenefluidsystem.hh
+++ b/dumux/material/fluidsystems/h2oairmesitylenefluidsystem.hh
@@ -225,7 +225,7 @@ public:
     /*!
      * \brief Given all mole fractions in a phase, return the phase
      *        density \f$\mathrm{[kg/m^3]}\f$.
-     * \param fluidState The fluid state 
+     * \param fluidState The fluid state
      * \param phaseIdx The index of the phase
      */
     using Base::density;
@@ -271,8 +271,8 @@ public:
 
     /*!
      * \brief Return the viscosity of a phase \f$\mathrm{[Pa s]}\f$.
-     * \param fluidState The fluid state 
-     * \param phaseIdx The index of the phase 
+     * \param fluidState The fluid state
+     * \param phaseIdx The index of the phase
      */
     using Base::viscosity;
     template <class FluidState>
@@ -343,7 +343,7 @@ public:
     /*!
      * \brief Given all mole fractions, return the diffusion
      *        coefficent in \f$\mathrm{[m^2/s]}\f$ of a component in a phase.
-     * \param fluidState The fluid state 
+     * \param fluidState The fluid state
      * \param phaseIdx The index of the phase
      * \param compIdx The index of the component
      */
@@ -418,7 +418,7 @@ public:
      * (fugacity equals the partial pressure of the component in the gas phase
      * respectively in the liquid phases it is the inverse of the
      * Henry coefficients scaled by pressure
-     * \param fluidState The fluid state 
+     * \param fluidState The fluid state
      * \param phaseIdx The index of the phase
      * \param compIdx The index of the component
      */
@@ -468,7 +468,7 @@ public:
     /*!
      * \brief Given all mole fractions in a phase, return the specific
      *        phase enthalpy \f$\mathrm{[J/kg]}\f$.
-     * \param fluidState The fluid state 
+     * \param fluidState The fluid state
      * \param phaseIdx The index of the phase
      */
     /*!
@@ -505,7 +505,7 @@ public:
     }
  /*!
      * \brief Return the heat capacity in \f$\mathrm{[J/(kg K)]}\f$.
-     * \param fluidState The fluid state 
+     * \param fluidState The fluid state
      * \param phaseIdx The index of the phase
      */
     using Base::heatCapacity;
@@ -515,10 +515,10 @@ public:
     {
         DUNE_THROW(Dune::NotImplemented, "FluidSystems::H2OAirMesitylene::heatCapacity()");
     }
-    
+
  /*!
      * \brief Return the thermal conductivity \f$\mathrm{[W/(m K)]}\f$.
-     * \param fluidState The fluid state 
+     * \param fluidState The fluid state
      * \param phaseIdx The index of the phase
      */
     using Base::thermalConductivity;
diff --git a/dumux/material/fluidsystems/h2oairxylenefluidsystem.hh b/dumux/material/fluidsystems/h2oairxylenefluidsystem.hh
index 40023d8fa40b0c82cd1d791c01c480ef169fbce1..6faf83710bbc9dc5a9aa5934c5f82d516d3f9441 100644
--- a/dumux/material/fluidsystems/h2oairxylenefluidsystem.hh
+++ b/dumux/material/fluidsystems/h2oairxylenefluidsystem.hh
@@ -223,8 +223,8 @@ public:
     /*!
      * \brief Given all mole fractions in a phase, return the phase
      *        density \f$\mathrm{[kg/m^3]}\f$.
-     * \param fluidState The fluid state 
-     * \param phaseIdx The index of the phase to consider 
+     * \param fluidState The fluid state
+     * \param phaseIdx The index of the phase to consider
      */
     using Base::density;
     template <class FluidState>
@@ -269,8 +269,8 @@ public:
 
     /*!
      * \brief Return the viscosity of a phase \f$\mathrm{[Pa s]}\f$.
-     * \param fluidState The fluid state 
-     * \param phaseIdx The index of the phase to consider 
+     * \param fluidState The fluid state
+     * \param phaseIdx The index of the phase to consider
      */
     using Base::viscosity;
     template <class FluidState>
@@ -344,8 +344,8 @@ public:
     /*!
      * \brief Given all mole fractions, return the diffusion
      *        coefficent \f$\mathrm{[m^2/s]}\f$ of a component in a phase.
-     * \param fluidState The fluid state 
-     * \param phaseIdx The index of the phase to consider 
+     * \param fluidState The fluid state
+     * \param phaseIdx The index of the phase to consider
      * \param compIdx The index of the component to consider
      */
     using Base::diffusionCoefficient;
@@ -413,8 +413,8 @@ public:
     /*!
      * \brief Returns the fugacity coefficient \f$\mathrm{[-]}\f$ of a component in a
      *        phase.
-     * \param fluidState The fluid state 
-     * \param phaseIdx The index of the phase to consider 
+     * \param fluidState The fluid state
+     * \param phaseIdx The index of the phase to consider
      * \param compIdx The index of the component to consider
      *
      * In this case, things are actually pretty simple. We have an ideal
@@ -469,8 +469,8 @@ public:
     /*!
      * \brief Given all mole fractions in a phase, return the specific
      *        phase enthalpy \f$\mathrm{[J/kg]}\f$.
-     * \param fluidState The fluid state 
-     * \param phaseIdx The index of the phase to consider 
+     * \param fluidState The fluid state
+     * \param phaseIdx The index of the phase to consider
      */
     /*!
      *  \todo This system neglects the contribution of gas-molecules in the liquid phase.
@@ -504,10 +504,10 @@ public:
         }
         DUNE_THROW(Dune::InvalidStateException, "Invalid phase index " << phaseIdx);
     }
-    
+
  /*!
      * \brief Return the heat capacity in \f$\mathrm{[J/(kg K)]}\f$.
-     * \param fluidState The fluid state 
+     * \param fluidState The fluid state
      * \param phaseIdx The index of the phase
      */
     using Base::heatCapacity;
@@ -517,10 +517,10 @@ public:
     {
         DUNE_THROW(Dune::NotImplemented, "FluidSystems::H2OAirXylene::heatCapacity()");
     }
-    
+
  /*!
      * \brief Return the thermal conductivity \f$\mathrm{[W/(m K)]}\f$.
-     * \param fluidState The fluid state 
+     * \param fluidState The fluid state
      * \param phaseIdx The index of the phase
      */
     using Base::thermalConductivity;
diff --git a/dumux/material/fluidsystems/h2on2fluidsystemkinetic.hh b/dumux/material/fluidsystems/h2on2fluidsystemkinetic.hh
index 48d548f8e16b3ddadf4adde856a50c959576db6f..01bbb59309736bff62923a98d74053815ea87313 100644
--- a/dumux/material/fluidsystems/h2on2fluidsystemkinetic.hh
+++ b/dumux/material/fluidsystems/h2on2fluidsystemkinetic.hh
@@ -254,13 +254,13 @@ public:
      *        phases was needed.
      *
      *        In the case of Henry, Raoult this would be
-     * 
-     *        ----- n-Comp w-Comp 
-     * 
+     *
+     *        ----- n-Comp w-Comp
+     *
      *        nPhase pn \f$\mathrm{x_n^n}\f$ pn \f$\mathrm{x_n^w}\f$
-     *      
+     *
      *        wPhase pv \f$\mathrm{w_w^w}\f$ H \f$\mathrm{x_w^n}\f$
-     *       
+     *
      *
      *        Plus additional relations for additional components.
      *
diff --git a/dumux/material/fluidsystems/spe5parametercache.hh b/dumux/material/fluidsystems/spe5parametercache.hh
index 159c726cb28636150d0927b288e73e4857bb1e70..16fab91feb1513bea4d5eebddc03cdd4aa645946 100644
--- a/dumux/material/fluidsystems/spe5parametercache.hh
+++ b/dumux/material/fluidsystems/spe5parametercache.hh
@@ -102,7 +102,7 @@ public:
      * concentration changed between two update*() calls. If more than
      * one concentration changed, call updatePhaseComposition() of
      * updatePhase()!
- 
+
      * \param fs An arbitrary fluid state
      * \param compIdx The index of the component to consider
      * \param phaseIdx The index of the fluid phase to consider
diff --git a/dumux/material/idealgas.hh b/dumux/material/idealgas.hh
index 51a4a0ed58c2be938b5d81fa3d125316172e2198..927b78f35559e0881d8690317881feab3f9e7dd1 100644
--- a/dumux/material/idealgas.hh
+++ b/dumux/material/idealgas.hh
@@ -42,9 +42,9 @@ public:
     /*!
      * \brief The density of the gas in \f$\mathrm{[kg/m^3]}\f$, depending on
      *        pressure, temperature and average molar mass of the gas.
-     * \param avgMolarMass The average molar mass of the gas 
-     * \param temperature The temperature of the gas 
-     * \param pressure The pressure of the gas 
+     * \param avgMolarMass The average molar mass of the gas
+     * \param temperature The temperature of the gas
+     * \param pressure The pressure of the gas
      */
     static Scalar density(Scalar avgMolarMass,
                           Scalar temperature,
@@ -54,8 +54,8 @@ public:
     /*!
      * \brief The pressure of the gas in \f$\mathrm{[Pa]}\f$, depending on
      *        the molar density and temperature.
-     * \param temperature The temperature of the gas 
-     * \param rhoMolar The molar density of the gas  
+     * \param temperature The temperature of the gas
+     * \param rhoMolar The molar density of the gas
      */
     static Scalar pressure(Scalar temperature,
                            Scalar rhoMolar)
@@ -64,8 +64,8 @@ public:
     /*!
      * \brief The molar density of the gas \f$\mathrm{[mol/m^3]}\f$,
      *        depending on pressure and temperature.
-     * \param temperature The temperature of the gas 
-     * \param pressure The pressure of the gas 
+     * \param temperature The temperature of the gas
+     * \param pressure The pressure of the gas
      */
     static Scalar molarDensity(Scalar temperature,
                                 Scalar pressure)
diff --git a/dune.module b/dune.module
index d98e48c614dc183e81f4ec6f28913c3d13fd4c8e..f5b67e6c62c8de74efa3738341bfd91a107c1048 100644
--- a/dune.module
+++ b/dune.module
@@ -3,3 +3,4 @@ Version: 2.9-git
 Maintainer: dumux@listserv.uni-stuttgart.de
 Depends: dune-grid (>= 2.4) dune-localfunctions (>= 2.4) dune-istl (>= 2.4)
 Suggests: dune-alugrid (>=2.4) dune-pdelab (>=2.0) dune-multidomain dune-cornerpoint
+Whitespace-Hook: Yes
diff --git a/test/decoupled/2p/buckleyleverettanalyticsolution.hh b/test/decoupled/2p/buckleyleverettanalyticsolution.hh
index e65b4f5acaa7b54fcca9e029cecfd7cf0edb3f81..a7d9d2cbd04b510955e22172ff13ae602091b336 100644
--- a/test/decoupled/2p/buckleyleverettanalyticsolution.hh
+++ b/test/decoupled/2p/buckleyleverettanalyticsolution.hh
@@ -96,7 +96,7 @@ private:
 
     void initializeAnalytic()
     {
-    	int size = problem_.gridView().size(0);
+        int size = problem_.gridView().size(0);
         analyticSolution_.resize(size);
         analyticSolution_ = 0;
         errorGlobal_.resize(size);
@@ -112,7 +112,7 @@ private:
      */
     void prepareAnalytic()
     {
-    	const auto& dummyElement = *problem_.gridView().template begin<0>();
+        const auto& dummyElement = *problem_.gridView().template begin<0>();
         const MaterialLawParams& materialLawParams(problem_.spatialParams().materialLawParams(dummyElement));
 
         swr_ = materialLawParams.swr();
@@ -129,12 +129,12 @@ private:
 
         if (CheckMaterialLaw<Scalar, MaterialLaw>::isLinear() && viscosityW == viscosityNW)
         {
-        	std::pair<Scalar, Scalar> entry;
-        	entry.first = 1 - snr_;
+            std::pair<Scalar, Scalar> entry;
+            entry.first = 1 - snr_;
 
-        	entry.second = vTot_  / (porosity * (1 - swr_ - snr_));
+            entry.second = vTot_  / (porosity * (1 - swr_ - snr_));
 
-        	frontParams_.push_back(entry);
+            frontParams_.push_back(entry);
         }
         else
         {
@@ -152,7 +152,7 @@ private:
 
         while (tangentSlopeNew >= tangentSlopeOld && sw1 < (1.0 - snr_))
         {
-        	tangentSlopeOld = tangentSlopeNew;
+            tangentSlopeOld = tangentSlopeNew;
             sw1 += deltaS_;
             fw1 = MaterialLaw::krw(materialLawParams, sw1)/viscosityW;
             fw1 /= (fw1 + MaterialLaw::krn(materialLawParams, sw1)/viscosityNW);
@@ -167,23 +167,23 @@ private:
         fw2 /= (fw2 + MaterialLaw::krn(materialLawParams, sw2)/viscosityNW);
         while (sw1 <= (1.0 - snr_))
         {
-        	std::pair<Scalar, Scalar> entry;
-        	entry.first = sw1;
+            std::pair<Scalar, Scalar> entry;
+            entry.first = sw1;
 
-        	Scalar dfwdsw = (fw2 - fw0)/(sw2 - sw0);
+            Scalar dfwdsw = (fw2 - fw0)/(sw2 - sw0);
 
-        	entry.second = vTot_  / porosity * dfwdsw;
+            entry.second = vTot_  / porosity * dfwdsw;
 
-        	frontParams_.push_back(entry);
+            frontParams_.push_back(entry);
 
-        	sw0 = sw1;
-        	sw1 = sw2;
-        	fw0 = fw1;
-        	fw1 = fw2;
+            sw0 = sw1;
+            sw1 = sw2;
+            fw0 = fw1;
+            fw1 = fw2;
 
-        	sw2 += deltaS_;
-        	fw2 = MaterialLaw::krw(materialLawParams, sw2)/viscosityW;
-        	fw2 /= (fw2 + MaterialLaw::krn(materialLawParams, sw2)/viscosityNW);
+            sw2 += deltaS_;
+            fw2 = MaterialLaw::krw(materialLawParams, sw2)/viscosityW;
+            fw2 /= (fw2 + MaterialLaw::krn(materialLawParams, sw2)/viscosityNW);
         }
         }
 
@@ -217,12 +217,12 @@ private:
 
         if (globalVolume > 0.0 && errorNorm > 0.0)
         {
-        	errorNorm = std::sqrt(errorNorm)/globalVolume;
-        	errorGlobal_ = errorNorm;
+            errorNorm = std::sqrt(errorNorm)/globalVolume;
+            errorGlobal_ = errorNorm;
         }
         else
         {
-        	errorGlobal_ = 0;
+            errorGlobal_ = 0;
         }
 
         return;
@@ -230,11 +230,11 @@ private:
 
     void updateExSol()
     {
-    	Scalar time = problem_.timeManager().time() + problem_.timeManager().timeStepSize();
+        Scalar time = problem_.timeManager().time() + problem_.timeManager().timeStepSize();
 
         // position of the fluid front
 
-    	Scalar xMax =  frontParams_[0].second * time;
+        Scalar xMax =  frontParams_[0].second * time;
 
         // iterate over vertices and get analytic saturation solution
         for (const auto& element : Dune::elements(problem_.gridView()))
@@ -245,21 +245,21 @@ private:
             int index = problem_.variables().index(element);
 
             if (globalPos[0] > xMax)
-            	analyticSolution_[index] = swr_;
+                analyticSolution_[index] = swr_;
             else
             {
-            	int size = frontParams_.size();
-            	if (size == 1)
-            	{
-            		analyticSolution_[index] = frontParams_[0].first;
-            	}
-            	else
-            	{
+                int size = frontParams_.size();
+                if (size == 1)
+                {
+                    analyticSolution_[index] = frontParams_[0].first;
+                }
+                else
+                {
                     // find x_f next to global coordinate of the vertex
                     for (int i = 0; i < size-1; i++)
                     {
-                    	Scalar x = frontParams_[i].second * time;
-                    	Scalar xMinus = frontParams_[i+1].second * time;
+                        Scalar x = frontParams_[i].second * time;
+                        Scalar xMinus = frontParams_[i+1].second * time;
                         if (globalPos[0] <= x && globalPos[0] > xMinus)
                         {
                             analyticSolution_[index] = frontParams_[i].first
@@ -269,7 +269,7 @@ private:
                         }
                     }
 
-            	}
+                }
             }
         }
 
@@ -296,7 +296,7 @@ public:
     template<class MultiWriter>
     void addOutputVtkFields(MultiWriter &writer)
     {
-    	int size = problem_.gridView().size(0);
+        int size = problem_.gridView().size(0);
         BlockVector *analyticSolution = writer.allocateManagedBuffer (size);
         BlockVector *errorGlobal = writer.allocateManagedBuffer (size);
         BlockVector *errorLocal = writer.allocateManagedBuffer (size);
@@ -319,7 +319,7 @@ public:
 
     void initialize(Scalar vTot)
     {
-    	vTot_ = vTot;
+        vTot_ = vTot;
         initializeAnalytic();
         prepareAnalytic();
     }
diff --git a/test/decoupled/2p/mcwhorteranalyticsolution.hh b/test/decoupled/2p/mcwhorteranalyticsolution.hh
index 6a5523bcc321ac0bfc46102a0090f0df230ecd94..4e564e305b28fbbd1a424e77ca72883661392797 100644
--- a/test/decoupled/2p/mcwhorteranalyticsolution.hh
+++ b/test/decoupled/2p/mcwhorteranalyticsolution.hh
@@ -78,7 +78,7 @@ private:
 
     void initializeAnalytic()
     {
-    	int size = problem_.gridView().size(0);
+        int size = problem_.gridView().size(0);
         analyticSolution_.resize(size);
         analyticSolution_=0;
         errorGlobal_.resize(size);
@@ -116,12 +116,12 @@ private:
 
         if (globalVolume > 0.0 && errorNorm > 0.0)
         {
-        	errorNorm = std::sqrt(errorNorm)/globalVolume;
-        	errorGlobal_ = errorNorm;
+            errorNorm = std::sqrt(errorNorm)/globalVolume;
+            errorGlobal_ = errorNorm;
         }
         else
         {
-        	errorGlobal_ = 0;
+            errorGlobal_ = 0;
         }
 
         return;
@@ -129,7 +129,7 @@ private:
 
     void prepareAnalytic()
     {
-    	const auto& dummyElement = *problem_.gridView().template begin<0>();
+        const auto& dummyElement = *problem_.gridView().template begin<0>();
         const MaterialLawParams& materialLawParams(problem_.spatialParams().materialLawParams(dummyElement));
 
         swr_ = materialLawParams.swr();
@@ -358,7 +358,7 @@ public:
      template<class MultiWriter>
      void addOutputVtkFields(MultiWriter &writer)
      {
-    	 int size = problem_.gridView().size(0);
+         int size = problem_.gridView().size(0);
          BlockVector *analyticSolution = writer.allocateManagedBuffer (size);
          BlockVector *errorGlobal = writer.allocateManagedBuffer (size);
          BlockVector *errorLocal = writer.allocateManagedBuffer (size);
diff --git a/test/decoupled/2p/test_mpfa2pproblem.hh b/test/decoupled/2p/test_mpfa2pproblem.hh
index 326a876bba4fe1c53341e69cc8eecdd8ffecd7a7..24670aabf9cb47b83eed8cd6ad8231020369b8af 100644
--- a/test/decoupled/2p/test_mpfa2pproblem.hh
+++ b/test/decoupled/2p/test_mpfa2pproblem.hh
@@ -208,7 +208,7 @@ ParentType(timeManager, gridView)
     int refinementFactor = 0;
     if (ParameterTree::tree().hasKey("Grid.RefinementFactor"))
     {
-    	refinementFactor = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Grid, RefinementFactor);
+        refinementFactor = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Grid, RefinementFactor);
     }
 
     this->grid().globalRefine(refinementFactor);
@@ -349,7 +349,7 @@ void boundaryTypesAtPos(BoundaryTypes &bcTypes, const GlobalPosition& globalPos)
 #elif  PROBLEM == 0
     if (globalPos[0] < eps_)
     {
-    	bcTypes.setAllDirichlet();
+        bcTypes.setAllDirichlet();
     }
     else if (globalPos[0] > this->bBoxMax()[0] - eps_)
     {
@@ -358,16 +358,16 @@ void boundaryTypesAtPos(BoundaryTypes &bcTypes, const GlobalPosition& globalPos)
     }
     else
     {
-    	bcTypes.setAllNeumann();
+        bcTypes.setAllNeumann();
     }
 #else
     if (globalPos[0] < eps_)
     {
-    	bcTypes.setAllDirichlet();
+        bcTypes.setAllDirichlet();
     }
     else
     {
-    	bcTypes.setAllNeumann();
+        bcTypes.setAllNeumann();
     }
 #endif
 }
@@ -414,7 +414,7 @@ void neumannAtPos(PrimaryVariables &values, const GlobalPosition& globalPos) con
 #elif PROBLEM == 0
     if (globalPos[0] > this->bBoxMax()[0] - eps_)
     {
-    	values[nPhaseIdx] = 3e-3;
+        values[nPhaseIdx] = 3e-3;
     }
 #endif
 }
diff --git a/test/decoupled/2p/test_mpfa2pspatialparams.hh b/test/decoupled/2p/test_mpfa2pspatialparams.hh
index 9ba6d4d982e03e49b8738501c684a7e6ddcc83a2..6323c4c42e80a402aa566acc5c46daf0b8a946db 100644
--- a/test/decoupled/2p/test_mpfa2pspatialparams.hh
+++ b/test/decoupled/2p/test_mpfa2pspatialparams.hh
@@ -94,7 +94,7 @@ public:
     Scalar porosity(const Element& element) const
     {
 #if PROBLEM == 0
-    	return 0.2;
+        return 0.2;
 #elif PROBLEM == 1
         return 0.3;
 #else
diff --git a/test/decoupled/2p2c/test_adaptive2p2c2dproblem.hh b/test/decoupled/2p2c/test_adaptive2p2c2dproblem.hh
index 5ea824b0e786eb4afc128fd56dd399f7cd28a89d..e696309d6c9f0f023d6a4da6953f4447b03ddae8 100644
--- a/test/decoupled/2p2c/test_adaptive2p2c2dproblem.hh
+++ b/test/decoupled/2p2c/test_adaptive2p2c2dproblem.hh
@@ -92,7 +92,7 @@ SET_INT_PROP(Adaptive2p2c2d, PressureFormulation, GET_PROP_TYPE(TypeTag, Indices
 /*!
  * \ingroup Adaptive2p2c
  * \ingroup IMPETtests
- * 
+ *
  * \brief test problem for the grid-adaptive sequential 2p2c model
  *
  * The domain is box shaped (2D). All sides are closed (Neumann 0 boundary)
diff --git a/test/decoupled/2p2c/test_adaptive2p2c3dproblem.hh b/test/decoupled/2p2c/test_adaptive2p2c3dproblem.hh
index 690e749235b8f2797599fbe19400bda69bcb8a0d..d67cb02a399c55c034913f42e38e34342a277ac8 100644
--- a/test/decoupled/2p2c/test_adaptive2p2c3dproblem.hh
+++ b/test/decoupled/2p2c/test_adaptive2p2c3dproblem.hh
@@ -92,7 +92,7 @@ SET_INT_PROP(Adaptive2p2c3d, PressureFormulation, GET_PROP_TYPE(TypeTag, Indices
 /*!
  * \ingroup Adaptive2p2c
  * \ingroup IMPETtests
- * 
+ *
  * \brief test problem for the grid-adaptive sequential 2p2c model
  *
  * The domain is box shaped (2D). All sides are closed (Neumann 0 boundary)
@@ -149,7 +149,7 @@ Adaptive2p2c3d(TimeManager &timeManager, const GridView& gridView) :
     GridCreator::grid().globalRefine(GET_PARAM_FROM_GROUP(TypeTag, int, GridAdapt, MaxLevel));
     this->setGrid(GridCreator::grid());
 
-    //Process parameter file 
+    //Process parameter file
     //Simulation Control
     const int outputInterval = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, int, Problem, OutputInterval);
     this->setOutputInterval(outputInterval);
@@ -188,7 +188,7 @@ Scalar temperatureAtPos(const GlobalPosition& globalPos) const
  */
 Scalar referencePressureAtPos(const GlobalPosition& globalPos) const
 {
-	return 1e6;
+    return 1e6;
 }
 
 /*!
diff --git a/test/freeflow/stokes/stokestestproblem.hh b/test/freeflow/stokes/stokestestproblem.hh
index 627e5f300cb1bb2667faa0e47607d5172c0e0906..f2b69f7d299c4be4cb660a6d1c30d7e45b9654d0 100644
--- a/test/freeflow/stokes/stokestestproblem.hh
+++ b/test/freeflow/stokes/stokestestproblem.hh
@@ -85,7 +85,7 @@ class StokesTestProblem : public StokesProblem<TypeTag>
 
         // Number of equations and grid dimension
         dim = GridView::dimension,
-	dimWorld = GridView::dimensionworld,
+    dimWorld = GridView::dimensionworld,
 
         // copy some indices for convenience
         massBalanceIdx = Indices::massBalanceIdx, //!< Index of the mass balance
diff --git a/test/implicit/2pminc/2pminctestspatialparams.hh b/test/implicit/2pminc/2pminctestspatialparams.hh
index 0fd33282558e863773aa6be113bec370d52c061d..19e4969d76fb9a1ac58de74ea32910f7005dace0 100644
--- a/test/implicit/2pminc/2pminctestspatialparams.hh
+++ b/test/implicit/2pminc/2pminctestspatialparams.hh
@@ -138,7 +138,7 @@ public:
     {
         if (nC == 0)
             return KFracture_;
-        
+
         return KMatrix_[nC];
     }
 
@@ -156,7 +156,7 @@ public:
     {
         if (nC == 0)
             return porosityFracture_;
-        
+
         return porosityMatrix_;
     }
 
diff --git a/test/implicit/2pminc/CMakeLists.txt b/test/implicit/2pminc/CMakeLists.txt
index e836a8b1691c6f2813706755c5297497269b1ed8..f3bbd15fa70dce483533e87639a0332572c51e77 100644
--- a/test/implicit/2pminc/CMakeLists.txt
+++ b/test/implicit/2pminc/CMakeLists.txt
@@ -7,7 +7,7 @@ add_dumux_test(test_box2pmincvol test_box2pminc test_box2pminc.cc
                ${CMAKE_CURRENT_BINARY_DIR}/box2pmincvol-00027.vtu
                --command "${CMAKE_CURRENT_BINARY_DIR}/test_box2pminc -ParameterFile ${CMAKE_CURRENT_SOURCE_DIR}/test_box2pmincvol.input"
                )
-  
+
 add_dumux_test(test_box2pmincdist test_box2pminc test_box2pminc.cc
                python ${CMAKE_SOURCE_DIR}/bin/runtest.py
                --script fuzzy
diff --git a/test/implicit/2pncmin/dissolutionproblem.hh b/test/implicit/2pncmin/dissolutionproblem.hh
index 5394006e9dc91a7fc2759e3e2dd6333be800f1a4..46ec8478bc01602a68f33d5ea11647b5a16f6a82 100644
--- a/test/implicit/2pncmin/dissolutionproblem.hh
+++ b/test/implicit/2pncmin/dissolutionproblem.hh
@@ -32,7 +32,7 @@
 
 namespace Dumux
 {
-    
+
 template <class TypeTag>
 class DissolutionProblem;
 
@@ -78,7 +78,7 @@ SET_INT_PROP(DissolutionProblem, Formulation, TwoPNCFormulation::pgSl);
  * The injected water phase migrates downwards due to increase in density as the precipitated salt dissolves.
  *
  * The model uses mole fractions of dissolved components and volume fractions of precipitated salt as primary variables. Make sure that the according units are used in the problem setup.
- * 
+ *
  * This problem uses the \ref TwoPNCMinModel.
  *
  * To run the simulation execute the following line in shell:
@@ -115,7 +115,7 @@ class DissolutionProblem : public ImplicitPorousMediaProblem<TypeTag>
             conti0EqIdx = Indices::conti0EqIdx,
             contiTotalMassIdx = conti0EqIdx + FluidSystem::AirIdx,
             precipNaClEqIdx = Indices::conti0EqIdx + FluidSystem::numComponents,
-            contiWEqIdx       =	conti0EqIdx + FluidSystem::H2OIdx,
+            contiWEqIdx       = conti0EqIdx + FluidSystem::H2OIdx,
 
             // Phase State
             wPhaseOnly = Indices::wPhaseOnly,
@@ -141,11 +141,11 @@ class DissolutionProblem : public ImplicitPorousMediaProblem<TypeTag>
     typedef Dune::FieldVector<Scalar, dimWorld> GlobalPosition;
 
 public:
-    DissolutionProblem(TimeManager &timeManager, 
+    DissolutionProblem(TimeManager &timeManager,
                        const GridView &gridView)
         : ParentType(timeManager, GridCreator::grid().leafGridView())
     {
-    
+
         outerSalinity_          = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, OuterSalinity);
         temperature_            = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, Temperature);
         reservoirPressure_      = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, ReservoirPressure);
@@ -159,7 +159,7 @@ public:
         innerPressure_          = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, InnerPressure);
         outerPressure_          = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, OuterPressure);
         reservoirSaturation_    = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, Problem, reservoirSaturation);
-        
+
         nTemperature_           = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, int, FluidSystem, NTemperature);
         nPressure_              = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, int, FluidSystem, NPressure);
         pressureLow_            = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, FluidSystem, PressureLow);
@@ -167,7 +167,7 @@ public:
         temperatureLow_         = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, FluidSystem, TemperatureLow);
         temperatureHigh_        = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, Scalar, FluidSystem, TemperatureHigh);
         name_                   = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, std::string, Problem, OutputName);
-        freqMassOutput_ 	    = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, int, Output, FreqMassOutput);
+        freqMassOutput_         = GET_RUNTIME_PARAM_FROM_GROUP(TypeTag, int, Output, FreqMassOutput);
         storageLastTimestep_    = Scalar(0);
         lastMassOutputTime_     = Scalar(0);
 
@@ -291,13 +291,13 @@ public:
      * influx.
      */
     void neumann(PrimaryVariables &values,
-				  const Element &element,
-				  const FVElementGeometry &fvGeometry,
-				  const Intersection &is,
-				  int scvIdx,
-				  int boundaryFaceIdx) const
+                  const Element &element,
+                  const FVElementGeometry &fvGeometry,
+                  const Intersection &is,
+                  int scvIdx,
+                  int boundaryFaceIdx) const
     {
-    	values = 0.0;
+        values = 0.0;
     }
 
     /*!
@@ -354,7 +354,7 @@ public:
         Scalar precipSalt = volVars.porosity() * volVars.molarDensity(wPhaseIdx)
                                             * volVars.saturation(wPhaseIdx)
                                             * pow(std::abs(moleFracNaCl_lPhase - moleFracNaCl_Max_lPhase), 1.0);
-                                            
+
         if (moleFracNaCl_lPhase < moleFracNaCl_Max_lPhase)
             precipSalt *= -1;
 
@@ -399,19 +399,19 @@ public:
 private:
 
     /*!
-	 * \brief Returns the molality of NaCl (mol NaCl / kg water) for a given mole fraction
-	 *
-	 * \param XlNaCl the XlNaCl [kg NaCl / kg solution]
-	 */
+     * \brief Returns the molality of NaCl (mol NaCl / kg water) for a given mole fraction
+     *
+     * \param XlNaCl the XlNaCl [kg NaCl / kg solution]
+     */
     static Scalar massTomoleFrac_(Scalar XlNaCl)
     {
-	   const Scalar Mw = 18.015e-3; /* molecular weight of water [kg/mol] */
-	   const Scalar Ms = 58.44e-3; /* molecular weight of NaCl  [kg/mol] */
+       const Scalar Mw = 18.015e-3; /* molecular weight of water [kg/mol] */
+       const Scalar Ms = 58.44e-3; /* molecular weight of NaCl  [kg/mol] */
 
-	   const Scalar X_NaCl = XlNaCl;
-	   /* XlNaCl: conversion from mass fraction to mol fraction */
-	   const Scalar xlNaCl = -Mw * X_NaCl / ((Ms - Mw) * X_NaCl - Ms);
-	   return xlNaCl;
+       const Scalar X_NaCl = XlNaCl;
+       /* XlNaCl: conversion from mass fraction to mol fraction */
+       const Scalar xlNaCl = -Mw * X_NaCl / ((Ms - Mw) * X_NaCl - Ms);
+       return xlNaCl;
     }
 
     int nTemperature_;
diff --git a/test/implicit/2pncmin/dissolutionspatialparams.hh b/test/implicit/2pncmin/dissolutionspatialparams.hh
index 7a35bd0646db335984b59c87363896ca6ee0df6d..dc40d8de56ad87e6349ae60fb817b9a177732334 100644
--- a/test/implicit/2pncmin/dissolutionspatialparams.hh
+++ b/test/implicit/2pncmin/dissolutionspatialparams.hh
@@ -98,7 +98,7 @@ class DissolutionSpatialparams : public ImplicitSpatialParams<TypeTag>
 public:
     DissolutionSpatialparams(const GridView &gridView)
         : ParentType(gridView),
-          K_(0)    
+          K_(0)
     {
         //set main diagonal entries of the permeability tensor to a value
         //setting to one value means: isotropic, homogeneous
@@ -136,11 +136,11 @@ public:
      *  could be defined, where globalPos is the vector including the global coordinates
      *  of the finite volume.
      */
-    const Dune::FieldMatrix<Scalar, dim, dim> &intrinsicPermeability(const Element &element, 
+    const Dune::FieldMatrix<Scalar, dim, dim> &intrinsicPermeability(const Element &element,
                                                                      const FVElementGeometry &fvGeometry,
                                                                      const int scvIdx) const
-    { 
-        return K_; 
+    {
+        return K_;
     }
 
     /*!
@@ -157,7 +157,7 @@ public:
      {
         return 1e-5;
      }
-     
+
     /*!
      * \brief Define the minimum porosity \f$[-]\f$ after clogging caused by mineralization
      *
diff --git a/test/implicit/3p/3pnispatialparams.hh b/test/implicit/3p/3pnispatialparams.hh
index 3bae1688fad9e657b0f744554b458e65b9b8780b..50fd59b8cdaa06c65a799d3425e54846a83d604d 100644
--- a/test/implicit/3p/3pnispatialparams.hh
+++ b/test/implicit/3p/3pnispatialparams.hh
@@ -155,7 +155,7 @@ public:
                                                const FVElementGeometry &fvGeometry,
                                                const int scvIdx) const
     {
-		return materialParams_;
+        return materialParams_;
     }
 
     /*!
diff --git a/test/implicit/3p3c/infiltrationspatialparameters.hh b/test/implicit/3p3c/infiltrationspatialparameters.hh
index 5b01955393f1a70ec6be2b981733953e3d0572f5..8743e13e2c1032fb3182bfedb6a88e5b96409fb8 100644
--- a/test/implicit/3p3c/infiltrationspatialparameters.hh
+++ b/test/implicit/3p3c/infiltrationspatialparameters.hh
@@ -65,7 +65,7 @@ class InfiltrationSpatialParams : public ImplicitSpatialParams<TypeTag>
     typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
     typedef typename Grid::ctype CoordScalar;
     enum {
-	dim = GridView::dimension,
+    dim = GridView::dimension,
         dimWorld=GridView::dimensionworld
     };
 
diff --git a/test/implicit/mpnc/combustionproblem1c.hh b/test/implicit/mpnc/combustionproblem1c.hh
index 2f6083efc32345f15b621b4d89973372baed3fe9..2bd31dc99b4f34f7f806ec78cae6f4fe87f4c0dc 100644
--- a/test/implicit/mpnc/combustionproblem1c.hh
+++ b/test/implicit/mpnc/combustionproblem1c.hh
@@ -65,7 +65,7 @@ class CombustionProblemOneComponent;
 
 namespace Properties {
 NEW_TYPE_TAG(CombustionProblemOneComponent,
-		INHERITS_FROM(BoxMPNCKinetic, CombustionSpatialParams));
+        INHERITS_FROM(BoxMPNCKinetic, CombustionSpatialParams));
 
 // Set the grid type
 SET_TYPE_PROP(CombustionProblemOneComponent, Grid, Dune::OneDGrid);
@@ -76,25 +76,25 @@ SET_TYPE_PROP(CombustionProblemOneComponent, LinearSolver, SuperLUBackend<TypeTa
 
 // Set the problem property
 SET_TYPE_PROP(CombustionProblemOneComponent,
-				Problem,
-				Dumux::CombustionProblemOneComponent<TTAG(CombustionProblemOneComponent)>);
+                Problem,
+                Dumux::CombustionProblemOneComponent<TTAG(CombustionProblemOneComponent)>);
 
 // Set fluid configuration
 SET_PROP(CombustionProblemOneComponent, FluidSystem){
 private:
-	typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
+    typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
 public:
-	typedef Dumux::FluidSystems::PureWaterSimpleFluidSystem<Scalar, /*useComplexRelations=*/false> type;
+    typedef Dumux::FluidSystems::PureWaterSimpleFluidSystem<Scalar, /*useComplexRelations=*/false> type;
 };
 
 // Set the newton controller
 SET_TYPE_PROP(CombustionProblemOneComponent, NewtonController,
-		Dumux::VeloModelNewtonController<TypeTag>);
+        Dumux::VeloModelNewtonController<TypeTag>);
 
 //! Set the default pressure formulation: either pw first or pn first
 SET_INT_PROP(CombustionProblemOneComponent,
-		PressureFormulation,
-		MpNcPressureFormulation::mostWettingFirst);
+        PressureFormulation,
+        MpNcPressureFormulation::mostWettingFirst);
 
 // Set the type used for scalar values
 SET_TYPE_PROP(CombustionProblemOneComponent, Scalar, double );
@@ -134,11 +134,11 @@ SET_BOOL_PROP(CombustionProblemOneComponent, EnableKinetic, false);
 SET_PROP(CombustionProblemOneComponent, FluidState)
 {
 private:
-	typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
+    typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
 private:
-	typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
+    typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
 public:
-	typedef Dumux::CompositionalFluidState<Scalar, FluidSystem> type;
+    typedef Dumux::CompositionalFluidState<Scalar, FluidSystem> type;
 };
 
 SET_BOOL_PROP(CombustionProblemOneComponent, UseMaxwellDiffusion, false);
@@ -160,602 +160,602 @@ SET_BOOL_PROP(CombustionProblemOneComponent, VelocityAveragingInModel, true);
  * \ingroup MpNcBoxproblems
  *
  * \brief Problem where water is injected from the left hand side into a porous media filled domain,
- *  	  which is supplied with energy from the right hand side to evaporate the water.
+ *        which is supplied with energy from the right hand side to evaporate the water.
  */
 template<class TypeTag>
 class CombustionProblemOneComponent: public ImplicitPorousMediaProblem<TypeTag> {
 
-	typedef CombustionProblemOneComponent<TypeTag> ThisType;
-	typedef ImplicitPorousMediaProblem<TypeTag> ParentType;
-	typedef typename GET_PROP_TYPE(TypeTag, Indices)Indices;
-	typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
-	typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
-	typedef typename GET_PROP_TYPE(TypeTag, PrimaryVariables) PrimaryVariables;
-	typedef typename GET_PROP_TYPE(TypeTag, BoundaryTypes) BoundaryTypes;
-	typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
-	/*!
-	 * \brief The fluid state which is used by the volume variables to
-	 *        store the thermodynamic state. This should be chosen
-	 *        appropriately for the model ((non-)isothermal, equilibrium, ...).
-	 *        This can be done in the problem.
-	 */
-	typedef typename GET_PROP_TYPE(TypeTag, FluidState) FluidState;
-	typedef typename FluidSystem::ParameterCache ParameterCache;
-	enum {dim = GridView::dimension}; // Grid and world dimension
-	enum {dimWorld = GridView::dimensionworld};
-	enum {numPhases = GET_PROP_VALUE(TypeTag, NumPhases)};
-	enum {numComponents = GET_PROP_VALUE(TypeTag, NumComponents)};
-	enum {s0EqIdx = Indices::s0Idx};
-	enum {p0EqIdx = Indices::p0Idx};
-	enum {conti00EqIdx = Indices::conti0EqIdx};
-	enum {temperature0Idx = Indices::temperatureIdx};
-	enum {energyEq0Idx = Indices::energyEqIdx};
-	enum {numEnergyEqs = Indices::numPrimaryEnergyVars};
-	enum {energyEqSolidIdx = energyEq0Idx + numEnergyEqs - 1};
-	enum {wPhaseIdx = FluidSystem::wPhaseIdx};
-	enum {nPhaseIdx = FluidSystem::nPhaseIdx};
-	enum {sPhaseIdx = FluidSystem::sPhaseIdx};
-	enum {wCompIdx = FluidSystem::H2OIdx};
-	enum {nCompIdx = FluidSystem::N2Idx};
-	enum {enableKinetic = GET_PROP_VALUE(TypeTag, EnableKinetic)};
-	enum {enableEnergy = GET_PROP_VALUE(TypeTag, EnableEnergy)};
-	enum {numEnergyEquations = GET_PROP_VALUE(TypeTag, NumEnergyEquations)};
-	enum {velocityOutput = GET_PROP_VALUE(TypeTag, VtkAddVelocities)};
-	enum {velocityAveragingInModel = GET_PROP_VALUE(TypeTag, VelocityAveragingInModel)};
-
-	// formulations
-	enum {
-		pressureFormulation = GET_PROP_VALUE(TypeTag, PressureFormulation),
-		mostWettingFirst = MpNcPressureFormulation::mostWettingFirst,
-		leastWettingFirst = MpNcPressureFormulation::leastWettingFirst
-	};
-
-	typedef typename GridView::template Codim<0>::Entity Element;
-	typedef typename GridView::template Codim<dim>::Entity Vertex;
-	typedef typename GridView::Intersection Intersection;
-	typedef typename GET_PROP_TYPE(TypeTag, FVElementGeometry) FVElementGeometry;
-	typedef typename GET_PROP_TYPE(TypeTag, MaterialLawParams) MaterialLawParams;
-	typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
-	typedef typename GET_PROP_TYPE(TypeTag, ElementVolumeVariables) ElementVolumeVariables;
-	typedef typename GET_PROP_TYPE(TypeTag, VolumeVariables) VolumeVariables;
-	typedef typename GET_PROP_TYPE(TypeTag, TimeManager) TimeManager;
-	typedef typename GET_PROP_TYPE(TypeTag, ElementBoundaryTypes) ElementBoundaryTypes;
-	typedef typename GET_PROP_TYPE(TypeTag, FluxVariables) FluxVariables;
-
-	typedef Dune::FieldVector<typename GridView::Grid::ctype, dim> DimVector;
-	typedef Dune::FieldVector<typename GridView::Grid::ctype, dimWorld> GlobalPosition;
-
-	typedef std::vector<Dune::FieldVector<Scalar, 1> > ScalarBuffer;
-	typedef std::array<ScalarBuffer, numPhases> PhaseBuffer;
-	typedef Dune::FieldVector<Scalar, dimWorld> VelocityVector;
-	typedef Dune::BlockVector<VelocityVector> VelocityField;
-	typedef std::array<VelocityField, numPhases> PhaseVelocityField;
+    typedef CombustionProblemOneComponent<TypeTag> ThisType;
+    typedef ImplicitPorousMediaProblem<TypeTag> ParentType;
+    typedef typename GET_PROP_TYPE(TypeTag, Indices)Indices;
+    typedef typename GET_PROP_TYPE(TypeTag, GridView) GridView;
+    typedef typename GET_PROP_TYPE(TypeTag, Scalar) Scalar;
+    typedef typename GET_PROP_TYPE(TypeTag, PrimaryVariables) PrimaryVariables;
+    typedef typename GET_PROP_TYPE(TypeTag, BoundaryTypes) BoundaryTypes;
+    typedef typename GET_PROP_TYPE(TypeTag, FluidSystem) FluidSystem;
+    /*!
+     * \brief The fluid state which is used by the volume variables to
+     *        store the thermodynamic state. This should be chosen
+     *        appropriately for the model ((non-)isothermal, equilibrium, ...).
+     *        This can be done in the problem.
+     */
+    typedef typename GET_PROP_TYPE(TypeTag, FluidState) FluidState;
+    typedef typename FluidSystem::ParameterCache ParameterCache;
+    enum {dim = GridView::dimension}; // Grid and world dimension
+    enum {dimWorld = GridView::dimensionworld};
+    enum {numPhases = GET_PROP_VALUE(TypeTag, NumPhases)};
+    enum {numComponents = GET_PROP_VALUE(TypeTag, NumComponents)};
+    enum {s0EqIdx = Indices::s0Idx};
+    enum {p0EqIdx = Indices::p0Idx};
+    enum {conti00EqIdx = Indices::conti0EqIdx};
+    enum {temperature0Idx = Indices::temperatureIdx};
+    enum {energyEq0Idx = Indices::energyEqIdx};
+    enum {numEnergyEqs = Indices::numPrimaryEnergyVars};
+    enum {energyEqSolidIdx = energyEq0Idx + numEnergyEqs - 1};
+    enum {wPhaseIdx = FluidSystem::wPhaseIdx};
+    enum {nPhaseIdx = FluidSystem::nPhaseIdx};
+    enum {sPhaseIdx = FluidSystem::sPhaseIdx};
+    enum {wCompIdx = FluidSystem::H2OIdx};
+    enum {nCompIdx = FluidSystem::N2Idx};
+    enum {enableKinetic = GET_PROP_VALUE(TypeTag, EnableKinetic)};
+    enum {enableEnergy = GET_PROP_VALUE(TypeTag, EnableEnergy)};
+    enum {numEnergyEquations = GET_PROP_VALUE(TypeTag, NumEnergyEquations)};
+    enum {velocityOutput = GET_PROP_VALUE(TypeTag, VtkAddVelocities)};
+    enum {velocityAveragingInModel = GET_PROP_VALUE(TypeTag, VelocityAveragingInModel)};
+
+    // formulations
+    enum {
+        pressureFormulation = GET_PROP_VALUE(TypeTag, PressureFormulation),
+        mostWettingFirst = MpNcPressureFormulation::mostWettingFirst,
+        leastWettingFirst = MpNcPressureFormulation::leastWettingFirst
+    };
+
+    typedef typename GridView::template Codim<0>::Entity Element;
+    typedef typename GridView::template Codim<dim>::Entity Vertex;
+    typedef typename GridView::Intersection Intersection;
+    typedef typename GET_PROP_TYPE(TypeTag, FVElementGeometry) FVElementGeometry;
+    typedef typename GET_PROP_TYPE(TypeTag, MaterialLawParams) MaterialLawParams;
+    typedef typename GET_PROP_TYPE(TypeTag, MaterialLaw) MaterialLaw;
+    typedef typename GET_PROP_TYPE(TypeTag, ElementVolumeVariables) ElementVolumeVariables;
+    typedef typename GET_PROP_TYPE(TypeTag, VolumeVariables) VolumeVariables;
+    typedef typename GET_PROP_TYPE(TypeTag, TimeManager) TimeManager;
+    typedef typename GET_PROP_TYPE(TypeTag, ElementBoundaryTypes) ElementBoundaryTypes;
+    typedef typename GET_PROP_TYPE(TypeTag, FluxVariables) FluxVariables;
+
+    typedef Dune::FieldVector<typename GridView::Grid::ctype, dim> DimVector;
+    typedef Dune::FieldVector<typename GridView::Grid::ctype, dimWorld> GlobalPosition;
+
+    typedef std::vector<Dune::FieldVector<Scalar, 1> > ScalarBuffer;
+    typedef std::array<ScalarBuffer, numPhases> PhaseBuffer;
+    typedef Dune::FieldVector<Scalar, dimWorld> VelocityVector;
+    typedef Dune::BlockVector<VelocityVector> VelocityField;
+    typedef std::array<VelocityField, numPhases> PhaseVelocityField;
 
 public:
-	CombustionProblemOneComponent(TimeManager &timeManager,
-			const GridView &gridView)
-	: ParentType(timeManager, gridView)
-	{
-	}
-
-	void init()
-	{
-			eps_ = 1e-6;
-			outputName_ = GET_RUNTIME_PARAM(TypeTag, std::string, Constants.outputName);
-			nRestart_ = GET_RUNTIME_PARAM(TypeTag, int, Constants.nRestart);
-			TInitial_ = GET_RUNTIME_PARAM(TypeTag, Scalar, InitialConditions.TInitial);
-			TRight_ = GET_RUNTIME_PARAM(TypeTag, Scalar, InitialConditions.TRight);
-			pnInitial_ = GET_RUNTIME_PARAM(TypeTag, Scalar,InitialConditions.pnInitial);
-			TBoundary_ = GET_RUNTIME_PARAM(TypeTag, Scalar,BoundaryConditions.TBoundary);
-			SwBoundary_ = GET_RUNTIME_PARAM(TypeTag, Scalar,BoundaryConditions.SwBoundary);
-			SwOneComponentSys_ = GET_RUNTIME_PARAM(TypeTag, Scalar,BoundaryConditions.SwOneComponentSys);
-			massFluxInjectedPhase_ = GET_RUNTIME_PARAM(TypeTag, Scalar,BoundaryConditions.massFluxInjectedPhase);
-			heatFluxFromRight_ = GET_RUNTIME_PARAM(TypeTag, Scalar,BoundaryConditions.heatFluxFromRight);
-			coldTime_ = GET_RUNTIME_PARAM(TypeTag, Scalar,BoundaryConditions.coldTime);
-
-		this->spatialParams().setInputInitialize();
-
-		ParentType::init();
-		this->model().initVelocityStuff();
-	}
-
-	/*!
-	 * \brief Returns the temperature \f$\mathrm{[K]}\f$ for an isothermal problem.
-	 *
-	 * This is not specific to the discretization.
-	 */
-	Scalar temperature() const
-	{	return TInitial_;}
-
-	/*!
-	 * \name Problem Params
-	 */
-	// \{
-	/*!
-	 * \brief The problem name.
-	 *
-	 * This is used as a prefix for files generated by the simulation.
-	 */
-	const std::string name() const
-	{	return outputName_;}
-
-	/*!
-	 * \brief Returns the temperature within the domain.
-	 *
-	 * \param element The finite element
-	 * \param fvGeometry The finite volume geometry of the element
-	 * \param scvIdx The local index of the sub-control volume
-	 *
-	 */
-	Scalar boxTemperature(const Element &element,
-			const FVElementGeometry &fvGeometry,
-			const unsigned int scvIdx) const
-	{	return TInitial_;}
-
-	// \}
-
-	/*!
-	 * \brief Called directly after the time integration.
-	 */
-	void postTimeStep()
-	{
-//    	if(this->timeManager().willBeFinished()){
-//			// write plot over Line to text files
-//			// each output gets it's writer object in order to allow for different header files
-//			PlotOverLine1D<TypeTag> writer;
+    CombustionProblemOneComponent(TimeManager &timeManager,
+            const GridView &gridView)
+    : ParentType(timeManager, gridView)
+    {
+    }
+
+    void init()
+    {
+            eps_ = 1e-6;
+            outputName_ = GET_RUNTIME_PARAM(TypeTag, std::string, Constants.outputName);
+            nRestart_ = GET_RUNTIME_PARAM(TypeTag, int, Constants.nRestart);
+            TInitial_ = GET_RUNTIME_PARAM(TypeTag, Scalar, InitialConditions.TInitial);
+            TRight_ = GET_RUNTIME_PARAM(TypeTag, Scalar, InitialConditions.TRight);
+            pnInitial_ = GET_RUNTIME_PARAM(TypeTag, Scalar,InitialConditions.pnInitial);
+            TBoundary_ = GET_RUNTIME_PARAM(TypeTag, Scalar,BoundaryConditions.TBoundary);
+            SwBoundary_ = GET_RUNTIME_PARAM(TypeTag, Scalar,BoundaryConditions.SwBoundary);
+            SwOneComponentSys_ = GET_RUNTIME_PARAM(TypeTag, Scalar,BoundaryConditions.SwOneComponentSys);
+            massFluxInjectedPhase_ = GET_RUNTIME_PARAM(TypeTag, Scalar,BoundaryConditions.massFluxInjectedPhase);
+            heatFluxFromRight_ = GET_RUNTIME_PARAM(TypeTag, Scalar,BoundaryConditions.heatFluxFromRight);
+            coldTime_ = GET_RUNTIME_PARAM(TypeTag, Scalar,BoundaryConditions.coldTime);
+
+        this->spatialParams().setInputInitialize();
+
+        ParentType::init();
+        this->model().initVelocityStuff();
+    }
+
+    /*!
+     * \brief Returns the temperature \f$\mathrm{[K]}\f$ for an isothermal problem.
+     *
+     * This is not specific to the discretization.
+     */
+    Scalar temperature() const
+    {   return TInitial_;}
+
+    /*!
+     * \name Problem Params
+     */
+    // \{
+    /*!
+     * \brief The problem name.
+     *
+     * This is used as a prefix for files generated by the simulation.
+     */
+    const std::string name() const
+    {   return outputName_;}
+
+    /*!
+     * \brief Returns the temperature within the domain.
+     *
+     * \param element The finite element
+     * \param fvGeometry The finite volume geometry of the element
+     * \param scvIdx The local index of the sub-control volume
+     *
+     */
+    Scalar boxTemperature(const Element &element,
+            const FVElementGeometry &fvGeometry,
+            const unsigned int scvIdx) const
+    {   return TInitial_;}
+
+    // \}
+
+    /*!
+     * \brief Called directly after the time integration.
+     */
+    void postTimeStep()
+    {
+//      if(this->timeManager().willBeFinished()){
+//          // write plot over Line to text files
+//          // each output gets it's writer object in order to allow for different header files
+//          PlotOverLine1D<TypeTag> writer;
 //
-//			// writer points: in the porous medium, not the outflow
-//			GlobalPosition pointOne ;
-//			GlobalPosition pointTwo ;
+//          // writer points: in the porous medium, not the outflow
+//          GlobalPosition pointOne ;
+//          GlobalPosition pointTwo ;
 //
-//			pointOne[0] = this->bBoxMin()[0] ;
-//			pointTwo[0] = (this->spatialParams().lengthPM()) ;
+//          pointOne[0] = this->bBoxMin()[0] ;
+//          pointTwo[0] = (this->spatialParams().lengthPM()) ;
 //
-//			writer.write(*this,
-//						  pointOne,
-//						  pointTwo,
-//						  "plotOverLine");
-//    	}
-
-		// Calculate storage terms of the individual phases
-		for (int phaseIdx = 0; phaseIdx < numEnergyEqs; ++phaseIdx) {
-			PrimaryVariables phaseStorage;
-			/* how much of the conserved quantity is in the system (in units of the balance eq.)
-			 * summed up storage term over the whole system
-			 */
-			this->model().globalPhaseStorage(phaseStorage, phaseIdx);
-
-			if (this->gridView().comm().rank() == 0) {
-				std::cout
-				<<"Storage in  "
-				<< FluidSystem::phaseName(phaseIdx)
-				<< "Phase: ["
-				<< phaseStorage
-				<< "]"
-				<< "\n";
-			}
-		}
-
-		// Calculate total storage terms
-		PrimaryVariables storage;
-		this->model().globalStorage(storage);
-
-		// Write mass balance information for rank 0
-		if (this->gridView().comm().rank() == 0) {
-			std::cout
-			<<"Storage total: [" << storage << "]"
-			<< "\n";
-		}
-	}
-
-	/*!
-	 * \brief Write a restart file?
-	 *yes of course:
-	 * - every nRestart_'th steps
-	 */
-	bool shouldWriteRestartFile() const
-	{
-		return this->timeManager().timeStepIndex() > 0 and
-		(this->timeManager().timeStepIndex() % nRestart_ == 0);
-	}
-
-	/*!
-	 * \brief Write a output file?
-	 */
-	bool shouldWriteOutput() const
-	{	return true;}
-
-	/*!
-	 * \brief Evaluate the source term for all phases within a given
-	 *        sub-control-volume.
-	 *
-	 * For this method, the \a values parameter stores the rate mass
-	 * of a component is generated or annihilate per volume
-	 * unit. Positive values mean that mass is created, negative ones
-	 * mean that it vanishes.
-	 *
-	 * \param priVars Stores the Dirichlet values for the conservation equations in
-	 * \param element The finite element
-	 * \param fvGeometry The finite volume geometry of the element
-	 * \param scvIdx The local index of the sub-control volume
-	 */
-	void source(PrimaryVariables & priVars,
-			const Element & element,
-			const FVElementGeometry & fvGeometry,
-			const unsigned int scvIdx) const
-	{
-		priVars = Scalar(0.0);
-		const GlobalPosition & globalPos = fvGeometry.subContVol[scvIdx].global;
-
-		const Scalar volume = fvGeometry.subContVol[scvIdx].volume;
-		const unsigned int numScv = fvGeometry.numScv; // box: numSCV, cc:1
-
-		if ((this->timeManager().time()) > coldTime_ )
-		{
-			if (onRightBoundaryPorousMedium_(globalPos))
-			{
-				if(numEnergyEquations>1) // for  2 and 3 energy equations
-				{
-					if(numEnergyEquations==2) {
-						// Testing the location of a vertex, but function is called for each associated scv. Compensate for that
-						priVars[energyEqSolidIdx] = heatFluxFromRight_ / volume / (Scalar)numScv;
-					}
-					else
-					// Multiply by volume, because afterwards we divide by it -> make independet of grid resolution
-					priVars[energyEq0Idx + sPhaseIdx] = heatFluxFromRight_ / volume / (Scalar)numScv;
-				}
-				else if (enableEnergy) {
-					// Multiply by volume, because afterwards we divide by it -> make independet of grid resolution
-					priVars[energyEq0Idx] = heatFluxFromRight_ / volume / (Scalar)numScv;
-				}
-			}
-		}
-	}
-
-	/*!
-	 * \brief Specifies which kind of boundary condition should be
-	 *        used for which equation on a given boundary segment.
-	 *
-	 * \param bTypes The bounentraldary types for the conservation equations
-	 * \param globalPos The global position
-
-	 */
-	void boundaryTypesAtPos(BoundaryTypes &bTypes,
+//          writer.write(*this,
+//                        pointOne,
+//                        pointTwo,
+//                        "plotOverLine");
+//      }
+
+        // Calculate storage terms of the individual phases
+        for (int phaseIdx = 0; phaseIdx < numEnergyEqs; ++phaseIdx) {
+            PrimaryVariables phaseStorage;
+            /* how much of the conserved quantity is in the system (in units of the balance eq.)
+             * summed up storage term over the whole system
+             */
+            this->model().globalPhaseStorage(phaseStorage, phaseIdx);
+
+            if (this->gridView().comm().rank() == 0) {
+                std::cout
+                <<"Storage in  "
+                << FluidSystem::phaseName(phaseIdx)
+                << "Phase: ["
+                << phaseStorage
+                << "]"
+                << "\n";
+            }
+        }
+
+        // Calculate total storage terms
+        PrimaryVariables storage;
+        this->model().globalStorage(storage);
+
+        // Write mass balance information for rank 0
+        if (this->gridView().comm().rank() == 0) {
+            std::cout
+            <<"Storage total: [" << storage << "]"
+            << "\n";
+        }
+    }
+
+    /*!
+     * \brief Write a restart file?
+     *yes of course:
+     * - every nRestart_'th steps
+     */
+    bool shouldWriteRestartFile() const
+    {
+        return this->timeManager().timeStepIndex() > 0 and
+        (this->timeManager().timeStepIndex() % nRestart_ == 0);
+    }
+
+    /*!
+     * \brief Write a output file?
+     */
+    bool shouldWriteOutput() const
+    {   return true;}
+
+    /*!
+     * \brief Evaluate the source term for all phases within a given
+     *        sub-control-volume.
+     *
+     * For this method, the \a values parameter stores the rate mass
+     * of a component is generated or annihilate per volume
+     * unit. Positive values mean that mass is created, negative ones
+     * mean that it vanishes.
+     *
+     * \param priVars Stores the Dirichlet values for the conservation equations in
+     * \param element The finite element
+     * \param fvGeometry The finite volume geometry of the element
+     * \param scvIdx The local index of the sub-control volume
+     */
+    void source(PrimaryVariables & priVars,
+            const Element & element,
+            const FVElementGeometry & fvGeometry,
+            const unsigned int scvIdx) const
+    {
+        priVars = Scalar(0.0);
+        const GlobalPosition & globalPos = fvGeometry.subContVol[scvIdx].global;
+
+        const Scalar volume = fvGeometry.subContVol[scvIdx].volume;
+        const unsigned int numScv = fvGeometry.numScv; // box: numSCV, cc:1
+
+        if ((this->timeManager().time()) > coldTime_ )
+        {
+            if (onRightBoundaryPorousMedium_(globalPos))
+            {
+                if(numEnergyEquations>1) // for  2 and 3 energy equations
+                {
+                    if(numEnergyEquations==2) {
+                        // Testing the location of a vertex, but function is called for each associated scv. Compensate for that
+                        priVars[energyEqSolidIdx] = heatFluxFromRight_ / volume / (Scalar)numScv;
+                    }
+                    else
+                    // Multiply by volume, because afterwards we divide by it -> make independet of grid resolution
+                    priVars[energyEq0Idx + sPhaseIdx] = heatFluxFromRight_ / volume / (Scalar)numScv;
+                }
+                else if (enableEnergy) {
+                    // Multiply by volume, because afterwards we divide by it -> make independet of grid resolution
+                    priVars[energyEq0Idx] = heatFluxFromRight_ / volume / (Scalar)numScv;
+                }
+            }
+        }
+    }
+
+    /*!
+     * \brief Specifies which kind of boundary condition should be
+     *        used for which equation on a given boundary segment.
+     *
+     * \param bTypes The bounentraldary types for the conservation equations
+     * \param globalPos The global position
+
+     */
+    void boundaryTypesAtPos(BoundaryTypes &bTypes,
                             const GlobalPosition &globalPos) const
-	{
-		// Default: Neumann
-		bTypes.setAllNeumann();
-
-		if(onRightBoundary_(globalPos) ) {
-			bTypes.setAllDirichlet();
-		}
-	}
-
-	/*!
-	 * \brief Evaluate the boundary conditions for a dirichlet
-	 *        boundary segment.
-	 *
-	 * \param priVars Stores the Dirichlet values for the conservation equations in
-	 *               \f$ [ \textnormal{unit of primary variable} ] \f$
-	 * \param globalPos The global position
-	 *
-	 * For this method, the \a values parameter stores primary variables.
-	 */
-	void dirichletAtPos(PrimaryVariables &priVars,
+    {
+        // Default: Neumann
+        bTypes.setAllNeumann();
+
+        if(onRightBoundary_(globalPos) ) {
+            bTypes.setAllDirichlet();
+        }
+    }
+
+    /*!
+     * \brief Evaluate the boundary conditions for a dirichlet
+     *        boundary segment.
+     *
+     * \param priVars Stores the Dirichlet values for the conservation equations in
+     *               \f$ [ \textnormal{unit of primary variable} ] \f$
+     * \param globalPos The global position
+     *
+     * For this method, the \a values parameter stores primary variables.
+     */
+    void dirichletAtPos(PrimaryVariables &priVars,
                         const GlobalPosition &globalPos) const
-	{
-		initial_(priVars, globalPos);
-	}
-
-	/*!
-	 * \brief Evaluate the boundary conditions for a neumann
-	 *        boundary segment.
-	 *
-	 * For this method, the \a values parameter stores the mass flux
-	 * in normal direction of each component. Negative values mean
-	 * influx.
-	 *
-	 * \param priVars Stores the Dirichlet values for the conservation equations in
-	 *               \f$ [ \textnormal{unit of primary variable} ] \f$
-	 * \param element The finite element
-	 * \param fvGeometry The finite volume geometry of the element
-	 * \param intersection The intersection between element and boundary
-	 * \param scvIdx The local index of the sub-control volume
-	 * \param boundaryFaceIdx The index of the boundary face
-	 * \param elemVolVars The element Volume Variables
-	 */
-
-	void solDependentNeumann(PrimaryVariables &priVars,
-			const Element &element,
-			const FVElementGeometry &fvGeometry,
-			const Intersection &intersection,
-			const int scvIdx,
-			const int boundaryFaceIdx,
-			const ElementVolumeVariables &elemVolVars) const {
-
-		const GlobalPosition & globalPos = fvGeometry.subContVol[scvIdx].global;
-		const Scalar massFluxInjectedPhase = massFluxInjectedPhase_;
-
-		FluidState fluidState;
-
-		const Scalar pn = elemVolVars[scvIdx].fluidState().pressure(nPhaseIdx);
-		const Scalar pw = elemVolVars[scvIdx].fluidState().pressure(wPhaseIdx);
-
-		const Scalar Tn = elemVolVars[scvIdx].fluidState().temperature(nPhaseIdx);
-		const Scalar Tw = elemVolVars[scvIdx].fluidState().temperature(wPhaseIdx);
-
-		fluidState.setPressure(nPhaseIdx, pn);
-		fluidState.setPressure(wPhaseIdx, pw);
-
-		if(numEnergyEquations == 3 ) { // only in this case the fluidstate hase several temperatures
-			fluidState.setTemperature(nPhaseIdx, Tn );
-			fluidState.setTemperature(wPhaseIdx, Tw );
-		}
-		else
-		fluidState.setTemperature(TBoundary_ );
-
-		ParameterCache dummyCache;
-		// This solves the system of equations defining x=x(p,T)
+    {
+        initial_(priVars, globalPos);
+    }
+
+    /*!
+     * \brief Evaluate the boundary conditions for a neumann
+     *        boundary segment.
+     *
+     * For this method, the \a values parameter stores the mass flux
+     * in normal direction of each component. Negative values mean
+     * influx.
+     *
+     * \param priVars Stores the Dirichlet values for the conservation equations in
+     *               \f$ [ \textnormal{unit of primary variable} ] \f$
+     * \param element The finite element
+     * \param fvGeometry The finite volume geometry of the element
+     * \param intersection The intersection between element and boundary
+     * \param scvIdx The local index of the sub-control volume
+     * \param boundaryFaceIdx The index of the boundary face
+     * \param elemVolVars The element Volume Variables
+     */
+
+    void solDependentNeumann(PrimaryVariables &priVars,
+            const Element &element,
+            const FVElementGeometry &fvGeometry,
+            const Intersection &intersection,
+            const int scvIdx,
+            const int boundaryFaceIdx,
+            const ElementVolumeVariables &elemVolVars) const {
+
+        const GlobalPosition & globalPos = fvGeometry.subContVol[scvIdx].global;
+        const Scalar massFluxInjectedPhase = massFluxInjectedPhase_;
+
+        FluidState fluidState;
+
+        const Scalar pn = elemVolVars[scvIdx].fluidState().pressure(nPhaseIdx);
+        const Scalar pw = elemVolVars[scvIdx].fluidState().pressure(wPhaseIdx);
+
+        const Scalar Tn = elemVolVars[scvIdx].fluidState().temperature(nPhaseIdx);
+        const Scalar Tw = elemVolVars[scvIdx].fluidState().temperature(wPhaseIdx);
+
+        fluidState.setPressure(nPhaseIdx, pn);
+        fluidState.setPressure(wPhaseIdx, pw);
+
+        if(numEnergyEquations == 3 ) { // only in this case the fluidstate hase several temperatures
+            fluidState.setTemperature(nPhaseIdx, Tn );
+            fluidState.setTemperature(wPhaseIdx, Tw );
+        }
+        else
+        fluidState.setTemperature(TBoundary_ );
+
+        ParameterCache dummyCache;
+        // This solves the system of equations defining x=x(p,T)
 //           FluidSystem::calculateEquilibriumMoleFractions(fluidState, dummyCache) ;
 
-		fluidState.setMoleFraction(wPhaseIdx, nCompIdx, 0.0);
-		fluidState.setMoleFraction(wPhaseIdx, wCompIdx, 1.0);
-		// compute density of injection phase
-		const Scalar density = FluidSystem::density(fluidState,
-				dummyCache,
-				wPhaseIdx);
-		fluidState.setDensity(wPhaseIdx, density);
-
-		for(int phaseIdx=0; phaseIdx<numPhases; phaseIdx++) {
-			const Scalar h = FluidSystem::enthalpy(fluidState,
-					dummyCache,
-					phaseIdx);
-			fluidState.setEnthalpy(phaseIdx, h);
-		}
-
-		const Scalar molarFlux = massFluxInjectedPhase / fluidState.averageMolarMass(wPhaseIdx);
-
-		// Default Neumann noFlow
-		priVars = 0.0;
-
-		if (onLeftBoundary_(globalPos))
-		{
-			if (enableKinetic) {
-				priVars[conti00EqIdx + numComponents*wPhaseIdx+wCompIdx] = - molarFlux * fluidState.moleFraction(wPhaseIdx, wCompIdx);
-				priVars[conti00EqIdx + numComponents*wPhaseIdx+nCompIdx] = - molarFlux * fluidState.moleFraction(wPhaseIdx, nCompIdx);
-				if (enableEnergy) {
-					if (numEnergyEquations == 3)
-					priVars[energyEq0Idx + wPhaseIdx] = - massFluxInjectedPhase * fluidState.enthalpy(wPhaseIdx);
-					else if(numEnergyEquations == 2)
-					priVars[energyEq0Idx] = - massFluxInjectedPhase * fluidState.enthalpy(wPhaseIdx);
-					else if(numEnergyEquations == 1)
-					priVars[energyEq0Idx] = - massFluxInjectedPhase * fluidState.enthalpy(wPhaseIdx);
-					else DUNE_THROW(Dune::InvalidStateException, "Formulation: " << pressureFormulation << " is invalid.");
-				}
-			}
-			else {
-				priVars[conti00EqIdx + wCompIdx] = - molarFlux * fluidState.moleFraction(wPhaseIdx, wCompIdx);;
-				priVars[conti00EqIdx + nCompIdx] = - molarFlux * fluidState.moleFraction(wPhaseIdx, nCompIdx);;
-				if (enableEnergy) {
-					if (numEnergyEquations == 3)
-					priVars[energyEq0Idx + wPhaseIdx] = - massFluxInjectedPhase * fluidState.enthalpy(wPhaseIdx);
-					else if(numEnergyEquations == 2)
-					priVars[energyEq0Idx] = - massFluxInjectedPhase * fluidState.enthalpy(wPhaseIdx);
-					else if(numEnergyEquations == 1)
-					priVars[energyEq0Idx] = - massFluxInjectedPhase * fluidState.enthalpy(wPhaseIdx);
-					else DUNE_THROW(Dune::InvalidStateException, "Formulation: " << pressureFormulation << " is invalid.");
-				}
-			}
-		}
-	}
-
-	/*!
-	 * \brief Evaluate the initial value for a control volume.
-	 *
-	 * For this method, the \a values parameter stores primary
-	 * variables.
-	 *
-	 * \param values Stores the Neumann values for the conservation equations in
-	 *               \f$ [ \textnormal{unit of conserved quantity} / (m^(dim-1) \cdot s )] \f$
-	 * \param globalPos The global position
-	 */
-	void initialAtPos(PrimaryVariables &values,
+        fluidState.setMoleFraction(wPhaseIdx, nCompIdx, 0.0);
+        fluidState.setMoleFraction(wPhaseIdx, wCompIdx, 1.0);
+        // compute density of injection phase
+        const Scalar density = FluidSystem::density(fluidState,
+                dummyCache,
+                wPhaseIdx);
+        fluidState.setDensity(wPhaseIdx, density);
+
+        for(int phaseIdx=0; phaseIdx<numPhases; phaseIdx++) {
+            const Scalar h = FluidSystem::enthalpy(fluidState,
+                    dummyCache,
+                    phaseIdx);
+            fluidState.setEnthalpy(phaseIdx, h);
+        }
+
+        const Scalar molarFlux = massFluxInjectedPhase / fluidState.averageMolarMass(wPhaseIdx);
+
+        // Default Neumann noFlow
+        priVars = 0.0;
+
+        if (onLeftBoundary_(globalPos))
+        {
+            if (enableKinetic) {
+                priVars[conti00EqIdx + numComponents*wPhaseIdx+wCompIdx] = - molarFlux * fluidState.moleFraction(wPhaseIdx, wCompIdx);
+                priVars[conti00EqIdx + numComponents*wPhaseIdx+nCompIdx] = - molarFlux * fluidState.moleFraction(wPhaseIdx, nCompIdx);
+                if (enableEnergy) {
+                    if (numEnergyEquations == 3)
+                    priVars[energyEq0Idx + wPhaseIdx] = - massFluxInjectedPhase * fluidState.enthalpy(wPhaseIdx);
+                    else if(numEnergyEquations == 2)
+                    priVars[energyEq0Idx] = - massFluxInjectedPhase * fluidState.enthalpy(wPhaseIdx);
+                    else if(numEnergyEquations == 1)
+                    priVars[energyEq0Idx] = - massFluxInjectedPhase * fluidState.enthalpy(wPhaseIdx);
+                    else DUNE_THROW(Dune::InvalidStateException, "Formulation: " << pressureFormulation << " is invalid.");
+                }
+            }
+            else {
+                priVars[conti00EqIdx + wCompIdx] = - molarFlux * fluidState.moleFraction(wPhaseIdx, wCompIdx);;
+                priVars[conti00EqIdx + nCompIdx] = - molarFlux * fluidState.moleFraction(wPhaseIdx, nCompIdx);;
+                if (enableEnergy) {
+                    if (numEnergyEquations == 3)
+                    priVars[energyEq0Idx + wPhaseIdx] = - massFluxInjectedPhase * fluidState.enthalpy(wPhaseIdx);
+                    else if(numEnergyEquations == 2)
+                    priVars[energyEq0Idx] = - massFluxInjectedPhase * fluidState.enthalpy(wPhaseIdx);
+                    else if(numEnergyEquations == 1)
+                    priVars[energyEq0Idx] = - massFluxInjectedPhase * fluidState.enthalpy(wPhaseIdx);
+                    else DUNE_THROW(Dune::InvalidStateException, "Formulation: " << pressureFormulation << " is invalid.");
+                }
+            }
+        }
+    }
+
+    /*!
+     * \brief Evaluate the initial value for a control volume.
+     *
+     * For this method, the \a values parameter stores primary
+     * variables.
+     *
+     * \param values Stores the Neumann values for the conservation equations in
+     *               \f$ [ \textnormal{unit of conserved quantity} / (m^(dim-1) \cdot s )] \f$
+     * \param globalPos The global position
+     */
+    void initialAtPos(PrimaryVariables &values,
                       const GlobalPosition &globalPos) const
-	{
-		initial_(values, globalPos);
-	}
+    {
+        initial_(values, globalPos);
+    }
 
 private:
-	// the internal method for the initial condition
-	void initial_(PrimaryVariables &priVars,
-			const GlobalPosition &globalPos) const
-	{
-		const Scalar curPos = globalPos[0];
-		const Scalar slope = (SwBoundary_-SwOneComponentSys_) / (this->spatialParams().lengthPM());
-		Scalar S[numPhases];
-		const Scalar thisSaturation = SwOneComponentSys_ + curPos * slope;
-
-		S[wPhaseIdx] = SwBoundary_;
-		if (inPM_(globalPos) ) {
-			S[wPhaseIdx] = thisSaturation;
-		}
-
-		S[nPhaseIdx] = 1. - S[wPhaseIdx];
-
-		//////////////////////////////////////
-		// Set saturation
-		//////////////////////////////////////
-		for (int i = 0; i < numPhases - 1; ++i) {
-			priVars[s0EqIdx + i] = S[i];
-		}
-
-		FluidState fluidState;
-
-		Scalar thisTemperature = TInitial_;
-		if(onRightBoundary_(globalPos))
-		thisTemperature = TRight_;
-
-		for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
-			fluidState.setSaturation(phaseIdx, S[phaseIdx]);
-
-			if(numEnergyEquations > 2) { // only in this case the fluidstate has more than one temperature
-				fluidState.setTemperature(phaseIdx, thisTemperature );
-			}
-			else
-			fluidState.setTemperature(thisTemperature );
-		}
-
-		//////////////////////////////////////
-		// Set temperature
-		//////////////////////////////////////
-		if(enableEnergy) {
-			for (int energyEqIdx=0; energyEqIdx< numEnergyEqs; ++energyEqIdx)
-			priVars[energyEq0Idx + energyEqIdx] = thisTemperature;
-		}
-
-		const MaterialLawParams &materialParams =
-		this->spatialParams().materialLawParamsAtPos(globalPos);
-		Scalar capPress[numPhases];
-
-		//obtain pc according to saturation
-		MaterialLaw::capillaryPressures(capPress, materialParams, fluidState);
-
-		Scalar p[numPhases];
-
-		p[wPhaseIdx] = pnInitial_ - std::abs(capPress[wPhaseIdx]);
-		p[nPhaseIdx] = p[wPhaseIdx] + std::abs(capPress[wPhaseIdx]);
-
-		for (int phaseIdx=0; phaseIdx<numPhases; phaseIdx++)
-		fluidState.setPressure(phaseIdx, p[phaseIdx]);
-
-		//////////////////////////////////////
-		// Set pressure
-		//////////////////////////////////////
-		if(pressureFormulation == mostWettingFirst) {
-			// This means that the pressures are sorted from the most wetting to the least wetting-1 in the primary variables vector.
-			// For two phases this means that there is one pressure as primary variable: pw
-			priVars[p0EqIdx] = p[wPhaseIdx];
-		}
-		else if(pressureFormulation == leastWettingFirst) {
-			// This means that the pressures are sorted from the least wetting to the most wetting-1 in the primary variables vector.
-			// For two phases this means that there is one pressure as primary variable: pn
-			priVars[p0EqIdx] = p[nPhaseIdx];
-		}
-		else DUNE_THROW(Dune::InvalidStateException, "Formulation: " << pressureFormulation << " is invalid.");
-
-		fluidState.setMoleFraction(wPhaseIdx, wCompIdx, 1.0);
-		fluidState.setMoleFraction(wPhaseIdx, nCompIdx, 0.0);
-
-		fluidState.setMoleFraction(nPhaseIdx, wCompIdx, 1.0);
-		fluidState.setMoleFraction(nPhaseIdx, nCompIdx, 0.0);
-
-		/* Difference between kinetic and MPNC:
-		 * number of component related primVar and how they are calculated (mole fraction, fugacities, resp.)          */
-		if(enableKinetic) {
-			for(int phaseIdx=0; phaseIdx < numPhases; ++ phaseIdx) {
-				for(int compIdx=0; compIdx <numComponents; ++compIdx) {
-					int offset = compIdx + phaseIdx * numComponents;
-					priVars[conti00EqIdx + offset] = fluidState.moleFraction(phaseIdx,compIdx);
-				}
-			}
-		}
-		else { //enableKinetic
-			int refPhaseIdx;
-
-			// on right boundary: reference is gas
-			refPhaseIdx = nPhaseIdx;
-
-			if(inPM_(globalPos)) {
-				refPhaseIdx = wPhaseIdx;
-			}
-
-			// obtain fugacities
-			typedef Dumux::ComputeFromReferencePhase<Scalar, FluidSystem> ComputeFromReferencePhase;
-			ParameterCache paramCache;
-			ComputeFromReferencePhase::solve(fluidState,
-					paramCache,
-					refPhaseIdx,
-					/*setViscosity=*/false,
-					/*setEnthalpy=*/false);
-
-			//////////////////////////////////////
-			// Set fugacities
-			//////////////////////////////////////
-			for (int compIdx = 0; compIdx < numComponents; ++compIdx) {
-				priVars[conti00EqIdx + compIdx] = fluidState.fugacity(refPhaseIdx,compIdx);
-			}
-		} // end not enableKinetic
-	}
-
-	/*!
-	 * \brief Give back whether the tested position (input) is a specific region (left) in the domain
-	 */
-	bool onLeftBoundary_(const GlobalPosition & globalPos) const
-	{	return globalPos[0] < this->bBoxMin()[0] + eps_;}
-
-	/*!
-	 * \brief Give back whether the tested position (input) is a specific region (right) in the domain
-	 */
-	bool onRightBoundary_(const GlobalPosition & globalPos) const
-	{	return globalPos[0] > this->bBoxMax()[0] - eps_;}
-
-	/*!
-	 * \brief Give back whether the tested position (input) is a specific region (right) in the domain
-	 * 	\todo this needs to be more sophisticated in order to allow for meshes with nodes not directly on the boundary
-	 */
-	bool onRightBoundaryPorousMedium_(const GlobalPosition & globalPos) const
-	{	return ( std::fabs(globalPos[0] - (this->spatialParams().lengthPM())) < eps_ );}
-
-	/*!
-	 * \brief Give back whether the tested position (input) is a specific region (right) in the domain
-	 */
-	bool inPM_(const GlobalPosition & globalPos) const
-	{	return
-		not this->spatialParams().inOutFlow(globalPos);
-	}
-
-	/*!
-	 * \brief Give back whether the tested position (input) is a specific region (down, (gravityDir)) in the domain
-	 */
-	bool onLowerBoundary_(const GlobalPosition & globalPos) const
-	{	return globalPos[dimWorld-1] < this->bBoxMin()[dimWorld-1] + eps_;}
-
-	/*!
-	 * \brief Give back whether the tested position (input) is a specific region (up, (gravityDir)) in the domain
-	 */
-	bool onUpperBoundary_(const GlobalPosition & globalPos) const
-	{	return globalPos[dimWorld-1] > this->bBoxMax()[dimWorld-1] - eps_;}
+    // the internal method for the initial condition
+    void initial_(PrimaryVariables &priVars,
+            const GlobalPosition &globalPos) const
+    {
+        const Scalar curPos = globalPos[0];
+        const Scalar slope = (SwBoundary_-SwOneComponentSys_) / (this->spatialParams().lengthPM());
+        Scalar S[numPhases];
+        const Scalar thisSaturation = SwOneComponentSys_ + curPos * slope;
+
+        S[wPhaseIdx] = SwBoundary_;
+        if (inPM_(globalPos) ) {
+            S[wPhaseIdx] = thisSaturation;
+        }
+
+        S[nPhaseIdx] = 1. - S[wPhaseIdx];
+
+        //////////////////////////////////////
+        // Set saturation
+        //////////////////////////////////////
+        for (int i = 0; i < numPhases - 1; ++i) {
+            priVars[s0EqIdx + i] = S[i];
+        }
+
+        FluidState fluidState;
+
+        Scalar thisTemperature = TInitial_;
+        if(onRightBoundary_(globalPos))
+        thisTemperature = TRight_;
+
+        for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
+            fluidState.setSaturation(phaseIdx, S[phaseIdx]);
+
+            if(numEnergyEquations > 2) { // only in this case the fluidstate has more than one temperature
+                fluidState.setTemperature(phaseIdx, thisTemperature );
+            }
+            else
+            fluidState.setTemperature(thisTemperature );
+        }
+
+        //////////////////////////////////////
+        // Set temperature
+        //////////////////////////////////////
+        if(enableEnergy) {
+            for (int energyEqIdx=0; energyEqIdx< numEnergyEqs; ++energyEqIdx)
+            priVars[energyEq0Idx + energyEqIdx] = thisTemperature;
+        }
+
+        const MaterialLawParams &materialParams =
+        this->spatialParams().materialLawParamsAtPos(globalPos);
+        Scalar capPress[numPhases];
+
+        //obtain pc according to saturation
+        MaterialLaw::capillaryPressures(capPress, materialParams, fluidState);
+
+        Scalar p[numPhases];
+
+        p[wPhaseIdx] = pnInitial_ - std::abs(capPress[wPhaseIdx]);
+        p[nPhaseIdx] = p[wPhaseIdx] + std::abs(capPress[wPhaseIdx]);
+
+        for (int phaseIdx=0; phaseIdx<numPhases; phaseIdx++)
+        fluidState.setPressure(phaseIdx, p[phaseIdx]);
+
+        //////////////////////////////////////
+        // Set pressure
+        //////////////////////////////////////
+        if(pressureFormulation == mostWettingFirst) {
+            // This means that the pressures are sorted from the most wetting to the least wetting-1 in the primary variables vector.
+            // For two phases this means that there is one pressure as primary variable: pw
+            priVars[p0EqIdx] = p[wPhaseIdx];
+        }
+        else if(pressureFormulation == leastWettingFirst) {
+            // This means that the pressures are sorted from the least wetting to the most wetting-1 in the primary variables vector.
+            // For two phases this means that there is one pressure as primary variable: pn
+            priVars[p0EqIdx] = p[nPhaseIdx];
+        }
+        else DUNE_THROW(Dune::InvalidStateException, "Formulation: " << pressureFormulation << " is invalid.");
+
+        fluidState.setMoleFraction(wPhaseIdx, wCompIdx, 1.0);
+        fluidState.setMoleFraction(wPhaseIdx, nCompIdx, 0.0);
+
+        fluidState.setMoleFraction(nPhaseIdx, wCompIdx, 1.0);
+        fluidState.setMoleFraction(nPhaseIdx, nCompIdx, 0.0);
+
+        /* Difference between kinetic and MPNC:
+         * number of component related primVar and how they are calculated (mole fraction, fugacities, resp.)          */
+        if(enableKinetic) {
+            for(int phaseIdx=0; phaseIdx < numPhases; ++ phaseIdx) {
+                for(int compIdx=0; compIdx <numComponents; ++compIdx) {
+                    int offset = compIdx + phaseIdx * numComponents;
+                    priVars[conti00EqIdx + offset] = fluidState.moleFraction(phaseIdx,compIdx);
+                }
+            }
+        }
+        else { //enableKinetic
+            int refPhaseIdx;
+
+            // on right boundary: reference is gas
+            refPhaseIdx = nPhaseIdx;
+
+            if(inPM_(globalPos)) {
+                refPhaseIdx = wPhaseIdx;
+            }
+
+            // obtain fugacities
+            typedef Dumux::ComputeFromReferencePhase<Scalar, FluidSystem> ComputeFromReferencePhase;
+            ParameterCache paramCache;
+            ComputeFromReferencePhase::solve(fluidState,
+                    paramCache,
+                    refPhaseIdx,
+                    /*setViscosity=*/false,
+                    /*setEnthalpy=*/false);
+
+            //////////////////////////////////////
+            // Set fugacities
+            //////////////////////////////////////
+            for (int compIdx = 0; compIdx < numComponents; ++compIdx) {
+                priVars[conti00EqIdx + compIdx] = fluidState.fugacity(refPhaseIdx,compIdx);
+            }
+        } // end not enableKinetic
+    }
+
+    /*!
+     * \brief Give back whether the tested position (input) is a specific region (left) in the domain
+     */
+    bool onLeftBoundary_(const GlobalPosition & globalPos) const
+    {   return globalPos[0] < this->bBoxMin()[0] + eps_;}
+
+    /*!
+     * \brief Give back whether the tested position (input) is a specific region (right) in the domain
+     */
+    bool onRightBoundary_(const GlobalPosition & globalPos) const
+    {   return globalPos[0] > this->bBoxMax()[0] - eps_;}
+
+    /*!
+     * \brief Give back whether the tested position (input) is a specific region (right) in the domain
+     *  \todo this needs to be more sophisticated in order to allow for meshes with nodes not directly on the boundary
+     */
+    bool onRightBoundaryPorousMedium_(const GlobalPosition & globalPos) const
+    {   return ( std::fabs(globalPos[0] - (this->spatialParams().lengthPM())) < eps_ );}
+
+    /*!
+     * \brief Give back whether the tested position (input) is a specific region (right) in the domain
+     */
+    bool inPM_(const GlobalPosition & globalPos) const
+    {   return
+        not this->spatialParams().inOutFlow(globalPos);
+    }
+
+    /*!
+     * \brief Give back whether the tested position (input) is a specific region (down, (gravityDir)) in the domain
+     */
+    bool onLowerBoundary_(const GlobalPosition & globalPos) const
+    {   return globalPos[dimWorld-1] < this->bBoxMin()[dimWorld-1] + eps_;}
+
+    /*!
+     * \brief Give back whether the tested position (input) is a specific region (up, (gravityDir)) in the domain
+     */
+    bool onUpperBoundary_(const GlobalPosition & globalPos) const
+    {   return globalPos[dimWorld-1] > this->bBoxMax()[dimWorld-1] - eps_;}
 
 private:
-	Scalar eps_;
-	int nTemperature_;
-	int nPressure_;
-	std::string outputName_;
-	int nRestart_;
-	Scalar TInitial_;
-	Scalar TRight_;
-
-	Scalar pnInitial_;
-
-	Dune::ParameterTree inputParameters_;
-	Scalar x_[numPhases][numComponents];
-
-	Scalar TBoundary_;
-	Scalar SwBoundary_;
-	Scalar SwOneComponentSys_;
-
-	Scalar massFluxInjectedPhase_;
-	Scalar heatFluxFromRight_;
-	PhaseVelocityField volumeDarcyVelocity_;
-	PhaseBuffer volumeDarcyMagVelocity_;
-	ScalarBuffer boxSurface_;
-	Scalar lengthPM_;
-	Scalar coldTime_;
+    Scalar eps_;
+    int nTemperature_;
+    int nPressure_;
+    std::string outputName_;
+    int nRestart_;
+    Scalar TInitial_;
+    Scalar TRight_;
+
+    Scalar pnInitial_;
+
+    Dune::ParameterTree inputParameters_;
+    Scalar x_[numPhases][numComponents];
+
+    Scalar TBoundary_;
+    Scalar SwBoundary_;
+    Scalar SwOneComponentSys_;
+
+    Scalar massFluxInjectedPhase_;
+    Scalar heatFluxFromRight_;
+    PhaseVelocityField volumeDarcyVelocity_;
+    PhaseBuffer volumeDarcyMagVelocity_;
+    ScalarBuffer boxSurface_;
+    Scalar lengthPM_;
+    Scalar coldTime_;
 
 public:
 
-	Dune::ParameterTree getInputParameters() const
-	{	return inputParameters_;}
+    Dune::ParameterTree getInputParameters() const
+    {   return inputParameters_;}
 };
 
 }
diff --git a/test/implicit/mpnc/combustionspatialparams.hh b/test/implicit/mpnc/combustionspatialparams.hh
index e37c81ac97308a35a5c786edb7bbde73c2036002..c0f9a392c39a27951d51fe9383b43712ed85a300 100644
--- a/test/implicit/mpnc/combustionspatialparams.hh
+++ b/test/implicit/mpnc/combustionspatialparams.hh
@@ -91,9 +91,9 @@ class CombustionSpatialParams : public ImplicitSpatialParams<TypeTag>
     enum {wPhaseIdx = FluidSystem::wPhaseIdx};
     enum {nPhaseIdx = FluidSystem::nPhaseIdx};
     enum {sPhaseIdx = FluidSystem::sPhaseIdx};
-    enum {numEnergyEquations  	= GET_PROP_VALUE(TypeTag, NumEnergyEquations)};
+    enum {numEnergyEquations    = GET_PROP_VALUE(TypeTag, NumEnergyEquations)};
     enum {numPhases = GET_PROP_VALUE(TypeTag, NumPhases)};
-    enum {enableEnergy 			= GET_PROP_VALUE(TypeTag, EnableEnergy)};
+    enum {enableEnergy          = GET_PROP_VALUE(TypeTag, EnableEnergy)};
 
     typedef typename GET_PROP_TYPE(TypeTag, VolumeVariables) VolumeVariables;
     typedef typename GET_PROP_TYPE(TypeTag, ElementVolumeVariables) ElementVolumeVariables;
@@ -124,28 +124,28 @@ public:
                 // BEWARE! First the input values have to be set, then the material parameters can be set
 
                 // this is the parameter value from file part
-                porosity_                 		= GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.PorousMedium.porosity);
+                porosity_                       = GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.PorousMedium.porosity);
 
-                intrinsicPermeabilityOutFlow_ 	= GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.OutFlow.permeabilityOutFlow);
-                porosityOutFlow_ 				= GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.OutFlow.porosityOutFlow);
+                intrinsicPermeabilityOutFlow_   = GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.OutFlow.permeabilityOutFlow);
+                porosityOutFlow_                = GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.OutFlow.porosityOutFlow);
                 solidThermalConductivityOutflow_ = GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.OutFlow.soilThermalConductivityOutFlow);
 
-                solidDensity_                		= GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.soil.density);
-                solidThermalConductivity_    		= GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.soil.thermalConductivity);
-                solidHeatCapacity_               		= GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.soil.heatCapacity);
+                solidDensity_                       = GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.soil.density);
+                solidThermalConductivity_           = GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.soil.thermalConductivity);
+                solidHeatCapacity_                      = GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.soil.heatCapacity);
 
-                interfacialTension_               	= GET_RUNTIME_PARAM(TypeTag, Scalar, Constants.interfacialTension);
+                interfacialTension_                 = GET_RUNTIME_PARAM(TypeTag, Scalar, Constants.interfacialTension);
 
                 Swr_            = GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.soil.Swr);
                 Snr_            = GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.soil.Snr);
 
-                characteristicLength_	= GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.PorousMedium.meanPoreSize);
+                characteristicLength_   = GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.PorousMedium.meanPoreSize);
                 intrinsicPermeability_  =  (std::pow(characteristicLength_,2.0)  * std::pow(porosity_,3.0)) / (150.0 * std::pow((1.0-porosity_),2.0)); // 1.69e-10 ; //
 
                 factorEnergyTransfer_         = GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.PorousMedium.factorEnergyTransfer);
                 factorMassTransfer_           = GET_RUNTIME_PARAM(TypeTag, Scalar, SpatialParams.PorousMedium.factorMassTransfer);
 
-                lengthPM_      			= GET_RUNTIME_PARAM(TypeTag, Scalar,Grid.lengthPM);
+                lengthPM_               = GET_RUNTIME_PARAM(TypeTag, Scalar,Grid.lengthPM);
 
             // residual saturations
             materialParams_.setSwr(Swr_) ;
@@ -172,10 +172,10 @@ public:
                                  const unsigned int scvIdx) const
     {
         const  GlobalPosition & globalPos =  fvGeometry.subContVol[scvIdx].global ;
-    	if ( inOutFlow(globalPos) )
-    		return intrinsicPermeabilityOutFlow_ ;
-    	else
-    		return intrinsicPermeability_ ;
+        if ( inOutFlow(globalPos) )
+            return intrinsicPermeabilityOutFlow_ ;
+        else
+            return intrinsicPermeability_ ;
     }
 
     /*! \brief Return the porosity \f$[-]\f$ of the soil
@@ -190,10 +190,10 @@ public:
                     const unsigned int scvIdx) const
     {
         const  GlobalPosition & globalPos =  fvGeometry.subContVol[scvIdx].global ;
-    	if ( inOutFlow(globalPos) )
-    		return porosityOutFlow_ ;
-    	else
-    		return porosity_ ;
+        if ( inOutFlow(globalPos) )
+            return porosityOutFlow_ ;
+        else
+            return porosity_ ;
     }
 
     /*!
@@ -254,7 +254,7 @@ public:
      * \param globalPos The position in global coordinates.*/
     const Scalar factorEnergyTransferAtPos(const  GlobalPosition & globalPos) const
     {
-    	return factorEnergyTransfer_ ;
+        return factorEnergyTransfer_ ;
     }
 
 
@@ -274,7 +274,7 @@ public:
      * \param globalPos The position in global coordinates.*/
     const Scalar factorMassTransferAtPos(const  GlobalPosition & globalPos) const
     {
-    	return factorMassTransfer_ ;
+        return factorMassTransfer_ ;
     }
 
 
@@ -306,10 +306,10 @@ public:
                                     const unsigned int scvIdx)const
     {
         const  GlobalPosition & globalPos =  fvGeometry.subContVol[scvIdx].global ;
-    	if ( inOutFlow(globalPos) )
-    		return solidThermalConductivityOutflow_ ;
-    	else
-    		return solidThermalConductivity_ ;
+        if ( inOutFlow(globalPos) )
+            return solidThermalConductivityOutflow_ ;
+        else
+            return solidThermalConductivity_ ;
     } // conductivity of solid  [W / (m K ) ]
 
     /*!
@@ -323,7 +323,7 @@ public:
      */
     Scalar lengthPM() const
     {
-    	return lengthPM_ ;
+        return lengthPM_ ;
     }
 
     /*!
@@ -331,7 +331,7 @@ public:
      */
     Scalar interfacialTension() const
     {
-    	return interfacialTension_ ;
+        return interfacialTension_ ;
     }
 
 private:
diff --git a/test/implicit/mpnc/evaporationatmosphereproblem.hh b/test/implicit/mpnc/evaporationatmosphereproblem.hh
index 3706d14def88d54ede84663c765338d5d190b854..2287165e12be91cde42d1629550e3d6566410ca0 100644
--- a/test/implicit/mpnc/evaporationatmosphereproblem.hh
+++ b/test/implicit/mpnc/evaporationatmosphereproblem.hh
@@ -175,7 +175,7 @@ SET_BOOL_PROP(EvaporationAtmosphereProblem, VelocityAveragingInModel, true);
  * \ingroup MpNcBoxproblems
  *
  * \brief Problem that simulates the coupled heat and mass transfer processes resulting form the evaporation of liquid water from
- *  	  a porous medium sub-domain into a gas filled "quasi-freeflow" sub-domain.
+ *        a porous medium sub-domain into a gas filled "quasi-freeflow" sub-domain.
  */
 template <class TypeTag>
 class EvaporationAtmosphereProblem
@@ -247,7 +247,7 @@ public:
      * \param gridView The grid view
      */
     EvaporationAtmosphereProblem(TimeManager &timeManager,
-    		const GridView &gridView)
+            const GridView &gridView)
         : ParentType(timeManager, gridView), gnuplot_(false)
     {
         this->timeManager().startNextEpisode(24.* 3600.);
@@ -403,7 +403,7 @@ public:
      */
     bool shouldWriteRestartFile() const
     {
-    	return this->timeManager().timeStepIndex() > 0 and
+        return this->timeManager().timeStepIndex() > 0 and
             (this->timeManager().timeStepIndex() % nRestart_  == 0);
     }
 
@@ -468,7 +468,7 @@ public:
     void dirichletAtPos(PrimaryVariables &priVars,
                         const GlobalPosition &globalPos) const
     {
-    	initial_(priVars, globalPos);
+        initial_(priVars, globalPos);
     }
 
     /*!
@@ -596,7 +596,7 @@ public:
      * \param fvGeometry The finite volume geometry of the element
      * \param scvIdx The local index of the sub-control volume
      *
-     * 		Positive values mean that mass is created, negative ones mean that it vanishes.
+     *      Positive values mean that mass is created, negative ones mean that it vanishes.
      */
     void source(PrimaryVariables & priVars,
                 const Element & element,
@@ -680,8 +680,8 @@ private:
 
         // temperature
         if(enableEnergy or numEnergyEquations)
-        	for (int energyEqIdx=0; energyEqIdx< numEnergyEqs; ++energyEqIdx)
-        		priVars[energyEq0Idx + energyEqIdx] = T;
+            for (int energyEqIdx=0; energyEqIdx< numEnergyEqs; ++energyEqIdx)
+                priVars[energyEq0Idx + energyEqIdx] = T;
 
         for (int phaseIdx=0; phaseIdx<numPhases; phaseIdx++)
              equilibriumFluidState.setPressure(phaseIdx, p[phaseIdx]);
diff --git a/test/implicit/mpnc/evaporationatmospherespatialparams.hh b/test/implicit/mpnc/evaporationatmospherespatialparams.hh
index 2b18af70c0a0daa35f14483f9fd3cc766de5f926..3d8dbdd56e9a16ac4dfcefbab327091a66240361 100644
--- a/test/implicit/mpnc/evaporationatmospherespatialparams.hh
+++ b/test/implicit/mpnc/evaporationatmospherespatialparams.hh
@@ -586,7 +586,7 @@ private:
     Scalar heightDomain_ ;
 
     AwnSurfaceParams    aWettingNonWettingSurfaceParams_;
-    AnsSurfaceParams 	aNonWettingSolidSurfaceParams_ ;
+    AnsSurfaceParams    aNonWettingSolidSurfaceParams_ ;
 
     AwnSurfaceParams    aWettingNonWettingSurfaceParamsFreeFlow_;
     AnsSurfaceParams    aNonWettingSolidSurfaceParamsFreeFlow_ ;
diff --git a/test/implicit/richards/richardsanalyticalproblem.hh b/test/implicit/richards/richardsanalyticalproblem.hh
index 5f50ebd023969d926750a30cc192b820602b6d64..861f06629c7935c5dd421727600491a42c1c1813 100644
--- a/test/implicit/richards/richardsanalyticalproblem.hh
+++ b/test/implicit/richards/richardsanalyticalproblem.hh
@@ -361,7 +361,7 @@ public:
                 Dune::FieldVector<Scalar, 1> values(0.0);
                 analyticalSolution(values, time, globalPos);
                 // add contributino of current quadrature point
-                l2error += (numericalSolution - values[0]) * (numericalSolution - values[0]) * 
+                l2error += (numericalSolution - values[0]) * (numericalSolution - values[0]) *
                     qIt->weight() * geometry.integrationElement(qIt->position());
                 l2analytic += values[0] * values[0] *
                     qIt->weight() * geometry.integrationElement(qIt->position());
diff --git a/test/multidomain/2cnistokes2p2cni/2cnistokes2p2cnispatialparams.hh b/test/multidomain/2cnistokes2p2cni/2cnistokes2p2cnispatialparams.hh
index bf456adea6862d7a73f292159fa5c3a1361f8a02..f08bb88887b00c0c3ec1ea47e69f7701d2333d71 100644
--- a/test/multidomain/2cnistokes2p2cni/2cnistokes2p2cnispatialparams.hh
+++ b/test/multidomain/2cnistokes2p2cni/2cnistokes2p2cnispatialparams.hh
@@ -124,8 +124,8 @@ public:
      * \param scvIdx The local index of the sub-control volume
      */
     const Scalar intrinsicPermeability(const Element &element,
-                                 	   const FVElementGeometry &fvGeometry,
-                                 	   const int scvIdx) const
+                                       const FVElementGeometry &fvGeometry,
+                                       const int scvIdx) const
     {
         return permeability_;
     }
diff --git a/test/multidomain/2cnistokes2p2cni/2p2cnisubproblem.hh b/test/multidomain/2cnistokes2p2cni/2p2cnisubproblem.hh
index 4c766e64396a2710900413588437281b68bef592..b1bcaadbc3e890625cdffc080590f7b16ceb3edb 100644
--- a/test/multidomain/2cnistokes2p2cni/2p2cnisubproblem.hh
+++ b/test/multidomain/2cnistokes2p2cni/2p2cnisubproblem.hh
@@ -44,7 +44,7 @@ class TwoPTwoCNISubProblem;
 namespace Properties
 {
 NEW_TYPE_TAG(TwoPTwoCNISubProblem,
-	INHERITS_FROM(BoxTwoPTwoCNI, SubDomain, TwoCNIStokesTwoPTwoCNISpatialParams));
+    INHERITS_FROM(BoxTwoPTwoCNI, SubDomain, TwoCNIStokesTwoPTwoCNISpatialParams));
 
 // Set the problem property
 SET_TYPE_PROP(TwoPTwoCNISubProblem, Problem, TwoPTwoCNISubProblem<TTAG(TwoPTwoCNISubProblem)>);
@@ -115,7 +115,7 @@ class TwoPTwoCNISubProblem : public ImplicitPorousMediaProblem<TypeTag>
         energyEqIdx = Indices::energyEqIdx
     };
     enum { // the indices of the primary variables
-		pNIdx = Indices::pressureIdx,
+        pNIdx = Indices::pressureIdx,
         sWIdx = Indices::switchIdx,
         temperatureIdx = Indices::temperatureIdx
     };
@@ -123,7 +123,7 @@ class TwoPTwoCNISubProblem : public ImplicitPorousMediaProblem<TypeTag>
         wCompIdx = Indices::wCompIdx,
         nCompIdx = Indices::nCompIdx
     };
-	enum { // the indices for the phase presence
+    enum { // the indices for the phase presence
         wPhaseOnly = Indices::wPhaseOnly,
         nPhaseOnly = Indices::nPhaseOnly,
         bothPhases = Indices::bothPhases
@@ -463,8 +463,8 @@ private:
 
     bool onBoundary_(const GlobalPosition &globalPos) const
     {
-    	return (onLeftBoundary_(globalPos) || onRightBoundary_(globalPos)
-    			|| onLowerBoundary_(globalPos) || onUpperBoundary_(globalPos));
+        return (onLeftBoundary_(globalPos) || onRightBoundary_(globalPos)
+                || onLowerBoundary_(globalPos) || onUpperBoundary_(globalPos));
     }
 
     static constexpr Scalar eps_ = 1e-8;
diff --git a/test/multidomain/2cnistokes2p2cni/stokes2cnisubproblem.hh b/test/multidomain/2cnistokes2p2cni/stokes2cnisubproblem.hh
index 132aae89fb4f7af2083ef3c844f9f3cb00092f91..c90bd7f466ed629fef4d18d8204b68ad7fc57a63 100644
--- a/test/multidomain/2cnistokes2p2cni/stokes2cnisubproblem.hh
+++ b/test/multidomain/2cnistokes2p2cni/stokes2cnisubproblem.hh
@@ -42,7 +42,7 @@ class Stokes2cniSubProblem;
 namespace Properties
 {
 NEW_TYPE_TAG(Stokes2cniSubProblem,
-	INHERITS_FROM(BoxStokesncni, SubDomain, TwoCNIStokesTwoPTwoCNISpatialParams));
+    INHERITS_FROM(BoxStokesncni, SubDomain, TwoCNIStokesTwoPTwoCNISpatialParams));
 
 // Set the problem property
 SET_TYPE_PROP(Stokes2cniSubProblem, Problem, Dumux::Stokes2cniSubProblem<TypeTag>);
diff --git a/test/multidomain/2cstokes2p2c/2cstokes2p2cproblem.hh b/test/multidomain/2cstokes2p2c/2cstokes2p2cproblem.hh
index 6f4f6f137c44059e9976174702d73d0acc1abe9d..5264684b071c90bc4cb86ac63b791343dab42c5f 100644
--- a/test/multidomain/2cstokes2p2c/2cstokes2p2cproblem.hh
+++ b/test/multidomain/2cstokes2p2c/2cstokes2p2cproblem.hh
@@ -173,17 +173,17 @@ class TwoCStokesTwoPTwoCProblem : public MultiDomainProblem<TypeTag>
     typedef typename GET_PROP_TYPE(Stokes2cTypeTag, FluidSystem) FluidSystem;
 
     enum { numEq1 = GET_PROP_VALUE(Stokes2cTypeTag, NumEq) };
-	enum { // indices in the Stokes domain
+    enum { // indices in the Stokes domain
         momentumXIdx1 = Stokes2cIndices::momentumXIdx, //!< Index of the x-component of the momentum balance
         momentumYIdx1 = Stokes2cIndices::momentumYIdx, //!< Index of the y-component of the momentum balance
         momentumZIdx1 = Stokes2cIndices::momentumZIdx, //!< Index of the z-component of the momentum balance
         lastMomentumIdx1 = Stokes2cIndices::lastMomentumIdx, //!< Index of the last component of the momentum balance
         massBalanceIdx1 = Stokes2cIndices::massBalanceIdx, //!< Index of the mass balance
         transportEqIdx1 = Stokes2cIndices::transportEqIdx, //!< Index of the transport equation
-	};
+    };
 
     enum { numEq2 = GET_PROP_VALUE(TwoPTwoCTypeTag, NumEq) };
-	enum { // indices in the Darcy domain
+    enum { // indices in the Darcy domain
         massBalanceIdx2 = TwoPTwoCIndices::pressureIdx,
         switchIdx2 = TwoPTwoCIndices::switchIdx,
         contiTotalMassIdx2 = TwoPTwoCIndices::contiNEqIdx,
@@ -195,7 +195,7 @@ class TwoCStokesTwoPTwoCProblem : public MultiDomainProblem<TypeTag>
         wPhaseIdx2 = TwoPTwoCIndices::wPhaseIdx,
         nPhaseIdx2 = TwoPTwoCIndices::nPhaseIdx
     };
-	enum { phaseIdx = GET_PROP_VALUE(Stokes2cTypeTag, PhaseIdx) };
+    enum { phaseIdx = GET_PROP_VALUE(Stokes2cTypeTag, PhaseIdx) };
     enum { transportCompIdx1 = Stokes2cIndices::transportCompIdx };
 
 
@@ -339,7 +339,7 @@ public:
     /*!
      * \brief Calculates fluxes and coupling terms at the interface
      *        for the Stokes model.
-     * 
+     *
      * Flux output files are created and the summarized flux is
      * written to a file.
      */
@@ -420,12 +420,12 @@ public:
                 // compute summarized fluxes for output
                 if (shouldWriteVaporFlux())
                 {
-                	int boundaryFaceIdx =
-                    	fvGeometry1.boundaryFaceIndex(firstFaceIdx, nodeInFace);
+                    int boundaryFaceIdx =
+                        fvGeometry1.boundaryFaceIndex(firstFaceIdx, nodeInFace);
 
-                	const BoundaryVariables1 boundaryVars1(this->sdProblem1(),
-                                                       	   sdElement1,
-                                                       	   fvGeometry1,
+                    const BoundaryVariables1 boundaryVars1(this->sdProblem1(),
+                                                           sdElement1,
+                                                           fvGeometry1,
                                                            boundaryFaceIdx,
                                                            elemVolVarsCur1,
                                                            /*onBoundary=*/true);
@@ -433,7 +433,7 @@ public:
                     advectiveWaterVaporFlux += computeAdvectiveVaporFluxes1(elemVolVarsCur1, boundaryVars1, vertInElem1);
                     diffusiveWaterVaporFlux += computeDiffusiveVaporFluxes1(elemVolVarsCur1, boundaryVars1, vertInElem1);
                     totalWaterComponentFlux += firstVertexDefect[vertInElem1][transportEqIdx1];
-				}
+                }
             }
         } // end loop over element faces on interface
 
@@ -551,14 +551,14 @@ public:
                 }
                 if (shouldWriteVaporFlux())
                 {
-	                if (!existing) // add phase storage only once per vertex
-    	                sumWaterFluxInGasPhase +=
-        	                this->localResidual2().evalPhaseStorageDerivative(vertInElem2);
-
-            	    totalWaterComponentFlux += secondVertexDefect[vertInElem2][contiWEqIdx2];
-                	sumWaterFluxInGasPhase +=
-                    	this->localResidual2().elementFluxes(vertInElem2);
-          		}
+                    if (!existing) // add phase storage only once per vertex
+                        sumWaterFluxInGasPhase +=
+                            this->localResidual2().evalPhaseStorageDerivative(vertInElem2);
+
+                    totalWaterComponentFlux += secondVertexDefect[vertInElem2][contiWEqIdx2];
+                    sumWaterFluxInGasPhase +=
+                        this->localResidual2().elementFluxes(vertInElem2);
+                }
             }
         }
 
diff --git a/test/multidomain/2cstokes2p2c/2p2csubproblem.hh b/test/multidomain/2cstokes2p2c/2p2csubproblem.hh
index e894ae3804e668d65784b14b6b95b7dc3cb91d73..bc08e3649290c5c2aab9b3884fd74d1348375806 100644
--- a/test/multidomain/2cstokes2p2c/2p2csubproblem.hh
+++ b/test/multidomain/2cstokes2p2c/2p2csubproblem.hh
@@ -439,7 +439,7 @@ private:
 
     Scalar runUpDistanceX_;
     Scalar initializationTime_;
-    std::ofstream outfile;    
+    std::ofstream outfile;
 };
 } //end namespace Dumux
 
diff --git a/test/multidomain/2cstokes2p2c/stokes2csubproblem.hh b/test/multidomain/2cstokes2p2c/stokes2csubproblem.hh
index d2f86c71c31482f57bb8495154326d545d80f42d..e52e70d7435ae4b2bb4b83091c540634fe50c926 100644
--- a/test/multidomain/2cstokes2p2c/stokes2csubproblem.hh
+++ b/test/multidomain/2cstokes2p2c/stokes2csubproblem.hh
@@ -260,9 +260,9 @@ public:
         {
             if (time > initializationTime_)
                 values.setDirichlet(pressureIdx);
-            else 
-            	if (!onLowerBoundary_(globalPos) && !onUpperBoundary_(globalPos))
-                	values.setDirichlet(pressureIdx);
+            else
+                if (!onLowerBoundary_(globalPos) && !onUpperBoundary_(globalPos))
+                    values.setDirichlet(pressureIdx);
         }
     }
 
@@ -311,8 +311,8 @@ public:
         if (onLeftBoundary_(globalPos)
                 && globalPos[1] > bBoxMin_[1] && globalPos[1] < bBoxMax_[1])
         {
-        	// rho*v*X at inflow
-        	values[transportEqIdx] = -xVelocity * density * refMassfrac();
+            // rho*v*X at inflow
+            values[transportEqIdx] = -xVelocity * density * refMassfrac();
         }
     }
 
@@ -337,7 +337,7 @@ public:
      */
     Scalar permeability(const Element &element,
                         const FVElementGeometry &fvGeometry,
-                 		const int scvIdx) const
+                        const int scvIdx) const
     {
         return spatialParams_.intrinsicPermeability(element,
                                                     fvGeometry,
diff --git a/tutorial/solutions_coupled/ex2_tutorialspatialparams_coupled.hh b/tutorial/solutions_coupled/ex2_tutorialspatialparams_coupled.hh
index 25ff6af52dac42c0d7b17297dcf4429c88f7dffd..01b623b9b5be4d3747ae6bdb54a9a6722ce9aeeb 100644
--- a/tutorial/solutions_coupled/ex2_tutorialspatialparams_coupled.hh
+++ b/tutorial/solutions_coupled/ex2_tutorialspatialparams_coupled.hh
@@ -77,7 +77,7 @@ class Ex2TutorialSpatialParamsCoupled: public ImplicitSpatialParams<TypeTag> /*@
     enum
     {
         dim = Grid::dimension,
-	dimWorld = Grid::dimensionworld
+    dimWorld = Grid::dimensionworld
     };
 
     // Get object types for function arguments
diff --git a/tutorial/solutions_coupled/ex5_tutorialproblem_coupled.hh b/tutorial/solutions_coupled/ex5_tutorialproblem_coupled.hh
index b4f3f2fb3f71bc6c733a3a9d989968e93d695c14..d373fe4649c2c823453c1804c326ed618853787c 100644
--- a/tutorial/solutions_coupled/ex5_tutorialproblem_coupled.hh
+++ b/tutorial/solutions_coupled/ex5_tutorialproblem_coupled.hh
@@ -142,7 +142,7 @@ public:
     //! as a VTK file
     bool shouldWriteOutput() const /*@\label{tutorial-coupled:output}@*/
     {
-    	const Scalar time = this->timeManager().time() + this->timeManager().timeStepSize();
+        const Scalar time = this->timeManager().time() + this->timeManager().timeStepSize();
 
         return
             this->timeManager().timeStepIndex() > 0 &&
@@ -174,14 +174,14 @@ public:
         values[Indices::pwIdx] = 200.0e3; // 200 kPa = 2 bar
         values[Indices::snIdx] = 1.0; // 1 % oil saturation on left boundary
 
-    	const Scalar time = this->timeManager().time() + this->timeManager().timeStepSize();
-    	const Scalar time_end = this->timeManager().endTime();
-    	Scalar injection_begin = time_end/5.0;
-    	Scalar injection_end = 4.0/5.0*time_end;
+        const Scalar time = this->timeManager().time() + this->timeManager().timeStepSize();
+        const Scalar time_end = this->timeManager().endTime();
+        Scalar injection_begin = time_end/5.0;
+        Scalar injection_end = 4.0/5.0*time_end;
 
-    	if(injection_begin < time && time < injection_end){
-			values[Indices::snIdx] = 1.0 - sin(M_PI*(time-injection_begin)/(injection_end-injection_begin));
-    	}
+        if(injection_begin < time && time < injection_end){
+            values[Indices::snIdx] = 1.0 - sin(M_PI*(time-injection_begin)/(injection_end-injection_begin));
+        }
     }
 
     //! Evaluates the boundary conditions for a Neumann boundary
diff --git a/tutorial/solutions_decoupled/ex2tutorialspatialparams_decoupled.hh b/tutorial/solutions_decoupled/ex2tutorialspatialparams_decoupled.hh
index 9bd8e05b758bc8c98656dc2cba2f21f5cd91abcf..df97e95213eeb83a4acca568186171136ed9b0f0 100644
--- a/tutorial/solutions_decoupled/ex2tutorialspatialparams_decoupled.hh
+++ b/tutorial/solutions_decoupled/ex2tutorialspatialparams_decoupled.hh
@@ -90,10 +90,10 @@ public:
      */
     const FieldMatrix& intrinsicPermeabilityAtPos(const GlobalPosition& globalPos) const
     {
-    	if(((globalPos[0]>25)&&(globalPos[0]<75))&&((globalPos[1]>15)&&(globalPos[1]<35)))
+        if(((globalPos[0]>25)&&(globalPos[0]<75))&&((globalPos[1]>15)&&(globalPos[1]<35)))
             return KLense_;
-    	else
-    		return K_;
+        else
+            return K_;
     }
 
     /**
@@ -107,9 +107,9 @@ public:
      */
     double porosityAtPos(const GlobalPosition& globalPos) const
     {
-    	if(((globalPos[0]>25)&&(globalPos[0]<75))&&((globalPos[1]>15)&&(globalPos[1]<35)))
-    	return 0.15;
-    	else
+        if(((globalPos[0]>25)&&(globalPos[0]<75))&&((globalPos[1]>15)&&(globalPos[1]<35)))
+        return 0.15;
+        else
         return 0.2;
     }
 
@@ -124,10 +124,10 @@ public:
      */
     const MaterialLawParams& materialLawParamsAtPos(const GlobalPosition& globalPos) const
     {
-    	if(((globalPos[0]>25)&&(globalPos[0]<75))&&((globalPos[1]>15)&&(globalPos[1]<35)))
-    	return materialParamsLense_;
-    	else
-    	return materialParams_;
+        if(((globalPos[0]>25)&&(globalPos[0]<75))&&((globalPos[1]>15)&&(globalPos[1]<35)))
+        return materialParamsLense_;
+        else
+        return materialParams_;
     }
 
     //! Constructor
@@ -137,19 +137,19 @@ public:
         for (int i = 0; i < dim; i++)
                 K_[i][i] = 1e-7;
         for (int i = 0; i < dim; i++)
-        	    KLense_[i][i] = 1e-9;
+                KLense_[i][i] = 1e-9;
 
 
-    	//set residual saturations
+        //set residual saturations
             materialParams_.setSwr(0.0);                /*@\label{tutorial-coupled:setLawParams}@*/
             materialParams_.setSnr(0.0);
-    	    materialParamsLense_.setSwr(0.0);                /*@\label{tutorial-coupled:setLawParams}@*/
+            materialParamsLense_.setSwr(0.0);                /*@\label{tutorial-coupled:setLawParams}@*/
             materialParamsLense_.setSnr(0.0);
 
             //parameters of Brooks & Corey Law
-    	    materialParams_.setPe(1000.0);
+            materialParams_.setPe(1000.0);
             materialParams_.setLambda(1.8);
-    	    materialParamsLense_.setPe(1500.0);
+            materialParamsLense_.setPe(1500.0);
             materialParamsLense_.setLambda(2.0);
     }
 
diff --git a/tutorial/solutions_decoupled/ex5_tutorialproblem_decoupled.hh b/tutorial/solutions_decoupled/ex5_tutorialproblem_decoupled.hh
index 05d4080bc2ec1cfea24831b722c742ea68f7fea3..96b2d523dd7a94eb50ce74f130d468be1b975f91 100644
--- a/tutorial/solutions_decoupled/ex5_tutorialproblem_decoupled.hh
+++ b/tutorial/solutions_decoupled/ex5_tutorialproblem_decoupled.hh
@@ -242,17 +242,17 @@ public:
      */
     void dirichlet(PrimaryVariables &values, const Intersection& intersection) const /*@\label{tutorial-decoupled:dirichlet}@*/
     {
-		values[pwIdx] = 2e5;
-		values[swIdx] = 0.0;
+        values[pwIdx] = 2e5;
+        values[swIdx] = 0.0;
 
-    	const Scalar time = this->timeManager().time() + this->timeManager().timeStepSize();
-    	const Scalar time_end = this->timeManager().endTime();
-    	Scalar injection_begin = time_end/5.0;
-    	Scalar injection_end = 4.0/5.0*time_end;
+        const Scalar time = this->timeManager().time() + this->timeManager().timeStepSize();
+        const Scalar time_end = this->timeManager().endTime();
+        Scalar injection_begin = time_end/5.0;
+        Scalar injection_end = 4.0/5.0*time_end;
 
-    	if(injection_begin < time && time < injection_end){
-			values[swIdx] = sin(M_PI*(time-injection_begin)/(injection_end-injection_begin));
-    	}
+        if(injection_begin < time && time < injection_end){
+            values[swIdx] = sin(M_PI*(time-injection_begin)/(injection_end-injection_begin));
+        }
     }
     //! Value for neumann boundary condition \f$ [\frac{kg}{m^3 \cdot s}] \f$ at position globalPos.
     /*! In case of a neumann boundary condition, the flux of matter