diff --git a/exercises/exercise-coupling-ff-pm/README.md b/exercises/exercise-coupling-ff-pm/README.md
index 1faa342529f5e560b36115583f664ce207341197..50756d140e289f8860026a2d608fc64c001806c1 100644
--- a/exercises/exercise-coupling-ff-pm/README.md
+++ b/exercises/exercise-coupling-ff-pm/README.md
@@ -15,12 +15,12 @@ Both single-phase flow and two-phase flow will be considered in the porous domai
 
 There are three sub folders: `interface` (Exercise 1), `models` (Exercise 2) and `turbulence` (Exercise 3).
 
-The task-related files for the simualtion exercises are as follows:
+The task-related files for the simulation exercises are as follows:
 * Three __main files__ for the three sub-tasks :`interface/main.cc`, `models/main.cc`, `turbulence/main.cc`,
 * Three __free flow problem files__: `interface/freeflowsubproblem.hh`, `models/freeflowsubproblem.hh`, `turbulence/freeflowsubproblem.hh`
 * Three __porous medium flow problem files__: `interface/porousmediumsubproblem.hh`, `models/porousmediumsubproblem.hh`, `turbulence/porousmediumsubproblem.hh`
 * Three __properties files__: `interface/properties.hh`, `models/properties.hh`, `turbulence/properties.hh`
-* The __input files__: `interface/params.input`, `models/parmas.input`, `turbulence/params.input`,
+* The __input files__: `interface/params.input`, `models/params.input`, `turbulence/params.input`,
 * The __spatial parameters files__: `1pspatialparams.hh`, `2pspatialparams.hh`
 
 In the main file, `TypeTags` for both submodels are defined, `FreeflowTypeTag` and `PorousMediumTypeTag`. These `TypeTags` collect all of the properties associated with each subproblem.
@@ -29,8 +29,8 @@ Since we use a monolithic coupling scheme, there is only one `Assembler` and one
 
 The problem files very much look like the "regular", uncoupled problem files seen in previous exercises, with the exception that they hold a pointer to the `CouplingManager`.
 This allows them to evaluate the coupling conditions and to exchange information between the coupled models.
-The coupling conditions are realized technically in terms of boundary condition. For instance,
-in `interface/freeflowsubproblem.hh`, `couplingNeumann` boundary conditions are set, which means that the free flow models evaluates the
+The coupling conditions are realized technically in terms of boundary conditions. For instance,
+in `interface/freeflowsubproblem.hh`, `couplingNeumann` boundary conditions are set, which means that the free flow model evaluates the
 mass and momentum fluxes coming from the porous domain and uses these values as boundary conditions at the interface.
 
 Note that certain checks are performed when combining different models, e.g., the fluid system has to be the same for both domains
@@ -39,10 +39,10 @@ and the sub-control-volume faces at the interface have to match.
 We will use a staggered grid to discretize the free-flow domain and a cell-centered finite volume method for the porous medium domain.
 Keep in mind that the staggered grid implementation distinguishes between face variables (velocity components) and
 cell center variables (all other variables).
-For this reason either the `CouplingManager::stokesCellCenterIdx` or the `CouplingManager::stokesFaceIdx` index sets need to be used, while for the porous medium all variables can be accessed with `CouplingManager::darcyIdx`.
+For this reason one distinguishes between `CouplingManager::stokesCellCenterIdx` and `CouplingManager::stokesFaceIdx` indices (see `main.cc`), while for the porous medium all variables can be accessed with `CouplingManager::darcyIdx`.
 
 __Task__:
-Take a closer look at the listed files before moving to the next part of the exercise:
+Take a closer look at the above listed files before moving to the next part of the exercise.
 
 
 ### 1. Changing the interface
@@ -149,7 +149,7 @@ the analytical solution of $`v_x`$ on the free flow domain. Play around with the
 __Task C: Change the shape of interface__
 
 Now we want to include a non-flat interface between the two domains.
-We use [`dune-subgrid`](https://doi.org/10.1007/s00607-009-0067-2) to construct two grids for the two domains from one common host grid.
+We use [`dune-subgrid`](https://doi.org/10.1007/s00607-009-0067-2) to construct two grids for the two domains from one common host grid. Thus, for the following tasks subgrid needs to be correctly installed.
 Our hostgrid will be a Dune-Yasp grid (`Dune::YaspGrid<2, Dune::TensorProductCoordinates<double, dim> >`)
 and the bounds and resolution of the domain will be set in the `params.input`file under the group `[Grid]`.
 This hostgrid, along with `elementSelector` functions defining some spatial cut of this domain, are passed to the grid manager to create each subdomain grid.
@@ -173,14 +173,14 @@ auto elementSelectorPorousMedium = [&](const auto& element)
 };
 ```
 
-Make sure that you have uncommented the lines including the grid managers in both properties files
+Make sure that you have uncommented the line for including the grid manager in the properties file, i.e.
 ```cpp
 #include <dumux/io/grid/gridmanager_sub.hh>
 ```
 
 and do the changes in the respective lines for the `Grid` property.
 
-The problem should now compile. However, an error occurs due to the coupling conditions.
+The problem should now compile. However, a runtime error occurs due to the coupling conditions.
 So far, we assumed a flat interface, therefore the normal momentum coupling condition
 
  $`[\sigma \cdot \mathbf{n}]^{FF} = p^{PM}`$
@@ -193,7 +193,7 @@ values.setCouplingNeumann(Indices::momentumYBalanceIdx);
  ```cpp
 values.setCouplingNeumann(scvf.directionIndex());
  ```
-in `main.cc` in the subfolder `interface`.
+in `freeflowsubproblem.hh` in the subfolder `interface`.
 
 The same is true for the BJS condition, however, here we need to consider the tangential direction:
 ```cpp
@@ -232,9 +232,9 @@ In the first task, the porous-medium model will be changed from a 1p2c system to
 Although a 2p2c system is plugged in, we still want to simulate the same situation as before, i.e., air with water vapor in both domains.
 The following changes have to be made in the porous-medium model (`models/properties.hh`):
 * Include the 2pnc model: include the respective headers and inherit from the new model `TwoPNC`
-* Exchange the spatial parameters for the 1-phase system by those for a 2-phase system (hint: two occurrences).
-* Since two phases are involved now, we do not need to use the `OnePAdapter` anymore. Change to property of the `FluidSystem` such that `H2OAir` is used directly.
-  Afterwards, set the `transportCompIdx` to `Indices::switchIdx`.
+* Exchange the spatial parameters for the 1-phase system by those for a 2-phase system.
+* Since two phases are involved now, we do not need to use the `OnePAdapter` anymore. Change the property of the `FluidSystem` such that `H2OAir` is used directly.
+  Afterwards, set the `transportCompIdx` to `Indices::switchIdx` in `porousmediumsubproblem.hh`.
 
 One big difference between the 1p and 2p model is, that the primary variables for which
 the problem is solved, are not fixed.
@@ -255,7 +255,7 @@ In the case under investigation, we want to use the gas pressure and the liquid
 However, if only the gas phase is present, the liquid saturation is always zero.
 In this case, the chosen formulation will set the given value as the mole fraction of water vapor in the gas phase.
 * To tell to program which phases are present in which parts of the domain at the beginning of the simulation,
-  you have to call `values.setState(MY_PHASE_PRESENCE);` in `initialAtPos()`. `MY_PHASE_PRESENCE` should be replaced with the correct value for the case where only a gas-phase is present.
+  you have to call `values.setState(MY_PHASE_PRESENCE);` in `initialAtPos()`. `MY_PHASE_PRESENCE` should be replaced with the correct value for the case where only a gas-phase (second phase of the fluid system) is present.
   Have a look at the `indices.hh` in the `2pnc` model (as the 2p2c model is a special case of the 2pnc model) in the subfolder `porousmediumflow` in your DuMuX directory. (hint: the numbering of phase indices begins with 0, the numbering of the phase presence states begins with 1. Take a look at your formulation to find out which phase index to use for the gas phase.)
 
 __Task B: Add output__:
@@ -263,7 +263,7 @@ __Task B: Add output__:
 In the next step, we want to add some output to the simulation. The standard method for providing simulation output for visualization is via a `VtkWriter`. These tools take the grid geometries of each domain, as well as the solutions and write spatial output that one can view in visualization tools such as paraview.
 
 Although this Vtk output is very useful, some output is more suited for other forms of visualization.
-Two examples of data ouput formats are `.csv` files and `.json` files.
+Two examples of data output formats are `.csv` files and `.json` files.
 From these files, data can be visualized using many different tools.
 In this exercise, data exported to a `.csv` file will be plotted using gnuplot, and data exported to a `.json` output will be plotted using the matplotlib python library.
 
@@ -274,23 +274,21 @@ The total storage can be plotted over time, whereas the water vapor flux will va
 First, we evaluate the storage term of the water component.
 * Have a look at the function `evaluateWaterMassStorageTerm()` in the `porousmediumsubproblem.hh`.
   Then implement a calculation of the total water mass:
-  $`\sum_{\alpha \in \textrm{g,l}} \left( \phi S_\alpha X^\text{w}_\alpha \varrho_\alpha V_\textrm{scv} \right)`$.
-  Afterwards, adapt the method `startStorageEvaluation()` such that the variable `initialWaterContent_` is initialized correctly using the `evaluateWaterMassStorageTerm()` method and assign that value also to the variable `lastWaterMass_`.
+  $`M^\text{w}:=\sum_{\alpha \in \textrm{g,l}} \left( \phi S_\alpha X^\text{w}_\alpha \varrho_\alpha V_\textrm{scv} \right)`$.
+* Afterwards, adapt the method `startStorageEvaluation()` such that the variable `initialWaterContent_` is
+  initialized correctly using the `evaluateWaterMassStorageTerm()` method and assign that value also to the variable `lastWaterMass_`.
 
-In order to evaluate how the stored water mass evolves over time, we should make sure the following is implemented:
-* Calculate the initial water mass at the beginning of the simulation.
-* For each time step, add the water mass loss to the output.
-* Based on the water mass loss you can derive the evaporation rate. The evaporation in [mm/d] can be calculated by:
-  $`e = \frac{\textrm{d}\, M^\text{w}}{\textrm{d}\, t} / A_\textrm{interface}`$.
-  (all boundaries have Neumann no-flow conditions and no sources are present, meaning the water can only leave the system
-  via the porous-medium free-flow interface.)
+After these functions are correctly implemented, the evaporation in [mm/d] is calculated in the code as
+$`e = \frac{\textrm{d}\, M^\text{w}}{\textrm{d}\, t} / A_\textrm{interface}`$.
+(all boundaries have Neumann no-flow conditions and no sources are present, meaning the water can only leave the system
+via the porous-medium free-flow interface.)
 
 During the simulation a `.csv` file is created (enable by setting `[Problem.ExportStorage] = true` in the `params.input` file). In addition, a gnuplot script `StorageOverTime.gp` is written and the mass loss and cumulative mass are plotted over time (`StorageOverTime.png`), when plotting is enabled with `[Problem.PlotStorage] = true`.
 
 Second, we want to know the distribution of the water mass fluxes across the interface.
 * Have a look at the function `evaluateInterfaceFluxes()` in the porous medium problem.
-  Use the facilities therein to return the values of ...massCouplingCondition... from the `couplingManager`
-  for each coupling scvf. Multiply this with the relevant face sizes, extrusion factor, massfraction and seconds per day to get the [mm/d] units as seen in the storage evaluation.
+  Use the facilities therein to return the values of `massCouplingCondition()` from the `couplingManager`
+  for each coupling scvf. Multiply this with the relevant face areas, extrusion factor, mass fraction, and seconds per day to get the [mm/d] units as seen in the storage evaluation.
 * When the `[Problem.ExportFluxes] = true` parameter is enabled, simulation data will be exported for this simulation to a `.json` file.
   This file can flexibly handle data of different types, which in this case is helpful as we have both temporal and spatial data.
 * Use the python plotting script `plotFluxes.py` to visualize the flux distribution across the surface at different times.
@@ -311,7 +309,7 @@ interface from any region of the porous medium and of course this poses numerica
 Therefore, the capillary pressure-saturation relationship has to be regularized for low saturations.
 The regularization has a strong influence on how long liquid water is present at the interface, see
 [Mosthaf (2014)](http://dx.doi.org/10.18419/opus-519).
-* Take a look at how the regularization is set in the `2pspatialparams.hh` and see how adapating
+* Take a look at how the regularization is set in the `2pspatialparams.hh` and see how adapting
   the parameter for `pcLowSw` and `pcHighSw` affects the results.
 
 The porous-medium model can now reach a liquid saturation of zero. As already explained above,
@@ -355,14 +353,13 @@ and in problem file:
 The includes for the NavierStokesNC model and the NavierStokesProblem are no longer needed and can be removed.
 
 Make sure your free flow problem inherits from the correct parent type:
-* Change the entry in the `FreeflowModel` definition accordingly (multi-component non-isothermal K-Omega SST model, SSTNCNI) in the properties file,
+* Change the entry in the `FreeflowModel` definition accordingly (multi-component non-isothermal K-Omega SST model, `SSTNCNI`) in the properties file,
 * Adapt the inheritance of the problem class in problem file (Use `RANSProblem<TypeTag>` rather than `NavierStokesStaggeredProblem<TypeTag>`).
-  Take a look into the two headers included above to see how the correct TypeTag and the Problem class the inherit from are called.
 
-Here, the turbulent free-flow is wall-bounded, meaning shear stress and turbulence in the flow develop primarily due to the walls.
-The porous medium at the bottom and also the channel wall are examples of such walls.
+Here, the turbulent free flow is wall-bounded, meaning shear stress and turbulence in the flow develop primarily due to the walls.
+The porous medium at the bottom and also the channel wall (upper boundary) act here as such walls.
 Within the problem the location of boundaries representing walls need to be defined.
-To do this, add the following function to the correct locations within the `boundaryTypes` function:
+To do this, add the following function to the upper and lower boundaries within the `boundaryTypes` function in `freeflowsubproblem.hh`:
 ```cpp
   values.setWall();
 ```
@@ -385,11 +382,14 @@ freeflowProblem->updateDynamicWallProperties(freeflowSol);
 In addition to designating the locations of walls,
 additional boundary conditions and initial conditions need to be set for the two new primary variables $k$ and $\omega$.
 In the `boundaryTypes` function, set both variables on all walls to be dirichlet, except for the right boundary, which should have outflow conditions.
-For the initial conditions, reynolds number specific base conditions should be applied everywhere.
-In the problem constructor, uncomment the code calcualting these terms,
-then apply the `turbulentKineticEnergy_`and `dissipation_` variables to their primary variables in all locations,
-except for on the wall boundaries, where these values can be set to zero.
-In addition, dirichlet constraints for the dissipation or $\omega$ variable will be set for all wall adjacent cells.
+
+For the initial conditions, Reynolds number specific base conditions should be applied everywhere.
+In the problem constructor, uncomment the code calculating these terms,
+then apply the `turbulentKineticEnergy_`and `dissipation_` variables to their primary variables in all locations.
+Within the dirichlet function for cell faces (`dirichlet(element, scvf)`),
+we also need to specify that these variables should be fixed to 0 at the wall.
+
+In addition, dirichlet cell constraints for the dissipation or $\omega$ variable will be set for all wall adjacent cells.
 This is done in the `isDirichletCell` function, as well as the `dirichlet` function already, and requires no further changes.
 
 Compile and run your new coupled problem and take a look at the results in Paraview.
@@ -408,7 +408,7 @@ you can use symmetric boundary conditions at the top boundary of your free flow
 values.setAllSymmetry();
 ```
 
-In addition, you have to remove the condition `onUpperBoundary_(globalPos)` from the `initialAtPos(globalPos)` method.
+In addition, you have to remove the condition `onUpperBoundary_(globalPos)` from the `initialAtPos(globalPos)` method and the `dirichlet(element, scvf)` method.
 
 __Task C__:
 
diff --git a/exercises/exercise-coupling-ff-pm/interface/main.cc b/exercises/exercise-coupling-ff-pm/interface/main.cc
index c7cce2fd18504e4452c88439b0c127a93f88489c..b869b4b98656934306bcae76640ed65dbc61709c 100644
--- a/exercises/exercise-coupling-ff-pm/interface/main.cc
+++ b/exercises/exercise-coupling-ff-pm/interface/main.cc
@@ -86,7 +86,7 @@ int main(int argc, char** argv)
 
     // ******************** uncomment this section for the last exercise ****************** //
 
-//        // use dune-subgrid to create the individual grids
+//     // use dune-subgrid to create the individual grids
 //     static constexpr int dim = 2;
 //     using HostGrid = Dune::YaspGrid<2, Dune::TensorProductCoordinates<double, dim> >;
 //     using HostGridManager = Dumux::GridManager<HostGrid>;
@@ -180,7 +180,6 @@ int main(int argc, char** argv)
     const auto freeflowName = getParam<std::string>("Problem.Name") + "_" + freeflowProblem->name();
     const auto porousMediumName = getParam<std::string>("Problem.Name") + "_" + porousMediumProblem->name();
 
-
     StaggeredVtkOutputModule<FreeflowGridVariables, decltype(freeflowSol)> freeflowVtkWriter(*freeflowGridVariables, freeflowSol, freeflowName);
     GetPropType<FreeflowTypeTag, Properties::IOFields>::initOutputModule(freeflowVtkWriter);
 
diff --git a/exercises/exercise-coupling-ff-pm/models/porousmediumsubproblem.hh b/exercises/exercise-coupling-ff-pm/models/porousmediumsubproblem.hh
index b05d07338875d6084c8d70c257a3cc7866096719..08c97a2ff8b862cb91c7aba403def39fc1ca83cb 100644
--- a/exercises/exercise-coupling-ff-pm/models/porousmediumsubproblem.hh
+++ b/exercises/exercise-coupling-ff-pm/models/porousmediumsubproblem.hh
@@ -97,7 +97,7 @@ public:
         exportFluxes_ = getParamFromGroup<bool>(this->paramGroup(), "Problem.ExportFluxes", false);
         if (exportFluxes_)
         {
-            fluxFileName_ = "flux" + getParam<std::string>("Problem.Name");
+            fluxFileName_ = "flux_" + getParam<std::string>("Problem.Name");
             simulationKey_ = getParam<std::string>("Problem.Name", "case1");
             initializeFluxOutput();
         }
@@ -127,7 +127,6 @@ public:
     {
         // TODO: dumux-course-task 2.B
         // Initialize `initialWaterContent_` and assign that to `lastWaterMass_`.
-
     }
 
     void startFluxEvaluation()
diff --git a/exercises/exercise-coupling-ff-pm/models/properties.hh b/exercises/exercise-coupling-ff-pm/models/properties.hh
index d3313c732a9bfbcf0211f5e9f8582710c8a59061..8c3faf3b87508d901d8acd99a33427dc38440501 100644
--- a/exercises/exercise-coupling-ff-pm/models/properties.hh
+++ b/exercises/exercise-coupling-ff-pm/models/properties.hh
@@ -56,7 +56,7 @@ namespace Dumux::Properties {
 // Create new type tags
 namespace TTag {
 // TODO: dumux-course-task 2.A
-// Change to property of the `FluidSystem` such that `H2OAir` is used directly.
+// Change the inheritance such that the correct model is used.
 struct PorousMediumOnePNC { using InheritsFrom = std::tuple<OnePNC, CCTpfaModel>; };
 struct FreeflowNC { using InheritsFrom = std::tuple<NavierStokesNC, StaggeredFreeFlowModel>; };
 } // end namespace TTag
@@ -82,6 +82,8 @@ template<class TypeTag>
 struct Problem<TypeTag, TTag::FreeflowNC> { using type = Dumux::FreeFlowSubProblem<TypeTag> ; };
 
 // The fluid system
+// TODO: dumux-course-task 2.A
+// Change to property of the `FluidSystem` such that `H2OAir` is used directly.
 template<class TypeTag>
 struct FluidSystem<TypeTag, TTag::PorousMediumOnePNC>
 {
diff --git a/exercises/exercise-coupling-ff-pm/turbulence/freeflowsubproblem.hh b/exercises/exercise-coupling-ff-pm/turbulence/freeflowsubproblem.hh
index 5a68f51303acfc7ef060cce07447f35855f5f69a..6859b9c661a999743db4a8fc4245d994505e69ca 100644
--- a/exercises/exercise-coupling-ff-pm/turbulence/freeflowsubproblem.hh
+++ b/exercises/exercise-coupling-ff-pm/turbulence/freeflowsubproblem.hh
@@ -54,7 +54,6 @@ class FreeFlowSubProblem : public NavierStokesStaggeredProblem<TypeTag>
 // TODO: dumux-course-task 3.A
 // Change the boundary types to Dumux::RANSBoundaryTypes<ModelTraits, ModelTraits::numEq()>
     using BoundaryTypes = Dumux::NavierStokesBoundaryTypes<GetPropType<TypeTag, Properties::ModelTraits>::numEq()>;
-
     using FVGridGeometry = GetPropType<TypeTag, Properties::GridGeometry>;
     using FVElementGeometry = typename FVGridGeometry::LocalView;
     using SubControlVolume = typename FVElementGeometry::SubControlVolume;
@@ -92,20 +91,20 @@ public:
                                                                      getParamFromGroup<std::string>(this->paramGroup(),
                                                                      "Problem.InterfaceDiffusionCoefficientAvg"));
 
-// TODO:   // ******************** uncomment this section for the first task 3.A ****************** //
-//         FluidSystem::init();
-//         Dumux::TurbulenceProperties<Scalar, dimWorld, true> turbulenceProperties;
-//         FluidState fluidState;
-//         const auto phaseIdx = 0;
-//         fluidState.setPressure(phaseIdx, refPressure_);
-//         fluidState.setTemperature(this->spatialParams().temperatureAtPos({}));
-//         fluidState.setMassFraction(phaseIdx, phaseIdx, 1.0);
-//         Scalar density = FluidSystem::density(fluidState, phaseIdx);
-//         Scalar kinematicViscosity = FluidSystem::viscosity(fluidState, phaseIdx) / density;
-//         Scalar diameter = this->gridGeometry().bBoxMax()[1] - this->gridGeometry().bBoxMin()[1];
-//         turbulentKineticEnergy_ = turbulenceProperties.turbulentKineticEnergy(refVelocity_, diameter, kinematicViscosity);
-//         dissipation_ = turbulenceProperties.dissipationRate(refVelocity_, diameter, kinematicViscosity);
-//         // ************************************************************************************* //
+        // TODO:   // ******************** uncomment this section for the first task 3.A ****************** //
+        // FluidSystem::init();
+        // Dumux::TurbulenceProperties<Scalar, dimWorld, true> turbulenceProperties;
+        // FluidState fluidState;
+        // const auto phaseIdx = 0;
+        // fluidState.setPressure(phaseIdx, refPressure_);
+        // fluidState.setTemperature(this->spatialParams().temperatureAtPos({}));
+        // fluidState.setMassFraction(phaseIdx, phaseIdx, 1.0);
+        // Scalar density = FluidSystem::density(fluidState, phaseIdx);
+        // Scalar kinematicViscosity = FluidSystem::viscosity(fluidState, phaseIdx) / density;
+        // Scalar diameter = this->gridGeometry().bBoxMax()[1] - this->gridGeometry().bBoxMin()[1];
+        // turbulentKineticEnergy_ = turbulenceProperties.turbulentKineticEnergy(refVelocity_, diameter, kinematicViscosity);
+        // dissipation_ = turbulenceProperties.dissipationRate(refVelocity_, diameter, kinematicViscosity);
+        // // ************************************************************************************* //
     }
 
     /*!
@@ -130,9 +129,9 @@ public:
             values.setDirichlet(Indices::energyEqIdx);
         }
 
-// TODO: dumux-course-task 3.A
-// set wall conditions for the turbulence model at the wall (values.setWall()) at the upper and lower boundary
-// set boundary conditions for the turbulence model primary variables everywhere (outflow on right boundary, otherwise dirichlet)
+        // TODO: dumux-course-task 3.A
+        // set wall conditions for the turbulence model at the walls (values.setWall()) corresponding to the upper and lower boundary
+        // set boundary conditions for the turbulence model primary variables everywhere (outflow on right boundary, otherwise dirichlet)
         if (onLowerBoundary_(globalPos))
         {
             values.setDirichlet(Indices::velocityXIdx);
@@ -182,6 +181,18 @@ public:
     {
         const auto globalPos = scvf.ipGlobal();
         PrimaryVariables values(initialAtPos(globalPos));
+
+        // TODO: dumux-course-task 3.A
+        // Add dirichlet conditions setting TKE and Dissipation to zero on the upper and lower walls.
+        // TODO: dumux-course-task 3.B
+        // Remove the condition `onUpperBoundary_(globalPos)` here.
+        // if (onUpperBoundary_(globalPos) || onLowerBoundary_(globalPos))
+        // {
+        //     values[Indices::turbulentKineticEnergyIdx] = 0.0;
+        //     values[Indices::dissipationIdx] = 0.0;
+        // }
+    }
+
         return values;
     }
 
@@ -290,15 +301,8 @@ public:
 
         // TODO: dumux-course-task 3.A
         // Set initial conditions for the TKE and the Dissipation. Values calculated in the constructor
-//         values[Indices::turbulentKineticEnergyIdx] = TODO??;
-//         values[Indices::dissipationIdx] = TODO??;
-//         // TODO: dumux-course-task 3.B
-//         // Remove the condition `onUpperBoundary_(globalPos)` here.
-//         if(onUpperBoundary_(globalPos) || onLowerBoundary_(globalPos))
-//         {
-//             values[Indices::turbulentKineticEnergyIdx] = 0.0;
-//             values[Indices::dissipationIdx] = 0.0;
-//         }
+        // values[Indices::turbulentKineticEnergyIdx] = TODO??;
+        // values[Indices::dissipationIdx] = TODO??;
 
         // TODO: dumux-course-task 3.B
         // Remove the condition `onUpperBoundary_(globalPos)` here.
diff --git a/exercises/exercise-coupling-ff-pm/turbulence/properties.hh b/exercises/exercise-coupling-ff-pm/turbulence/properties.hh
index 14bc29ff7ff023b927eb5ad26b808e6d8fbf04be..decd6cc83d3f8781fe4007713369b7c938745fe3 100644
--- a/exercises/exercise-coupling-ff-pm/turbulence/properties.hh
+++ b/exercises/exercise-coupling-ff-pm/turbulence/properties.hh
@@ -18,13 +18,12 @@
  *****************************************************************************/
 /*!
  * \file
-
  * \brief The coupled exercise properties file or the turbulent case.
  */
 #ifndef DUMUX_EXERCISE_COUPLED_TURBULENCE_PROPERTIES_HH
 #define DUMUX_EXERCISE_COUPLED_TURBULENCE_PROPERTIES_HH
 
-// Both domain
+// Both domains
 #include <dune/grid/yaspgrid.hh>
 #include <dumux/multidomain/staggeredtraits.hh>
 #include <dumux/multidomain/boundary/stokesdarcy/couplingmanager.hh>
diff --git a/exercises/solution/exercise-coupling-ff-pm/1pspatialparams.hh b/exercises/solution/exercise-coupling-ff-pm/1pspatialparams.hh
index b8b230fb4a8daba93cb3e5824947703a84044c88..fc7f7093f84961aed8bbe19ed6170aebd9cbe74c 100644
--- a/exercises/solution/exercise-coupling-ff-pm/1pspatialparams.hh
+++ b/exercises/solution/exercise-coupling-ff-pm/1pspatialparams.hh
@@ -54,7 +54,7 @@ public:
         permeability_ = getParam<Scalar>("SpatialParams.Permeability");
         porosity_ = getParam<Scalar>("SpatialParams.Porosity");
         alphaBJ_ = getParam<Scalar>("SpatialParams.AlphaBeaversJoseph");
-        temperature_ = getParam<Scalar>("SpatialParams.PorousMediumTemperature");
+        temperature_ = getParam<Scalar>("SpatialParams.Temperature");
     }
 
     /*!
@@ -88,7 +88,6 @@ public:
     { return temperature_; }
 
 
-
 private:
     Scalar permeability_;
     Scalar porosity_;
diff --git a/exercises/solution/exercise-coupling-ff-pm/2pspatialparams.hh b/exercises/solution/exercise-coupling-ff-pm/2pspatialparams.hh
index 1b6cc619770f4eb52dc2f38e15e1e985b092cc53..f1a4938b7f4faf89d2efba8b670e6ee394877f0f 100644
--- a/exercises/solution/exercise-coupling-ff-pm/2pspatialparams.hh
+++ b/exercises/solution/exercise-coupling-ff-pm/2pspatialparams.hh
@@ -59,7 +59,7 @@ public:
         permeability_ = getParam<Scalar>("SpatialParams.Permeability");
         porosity_ = getParam<Scalar>("SpatialParams.Porosity");
         alphaBJ_ = getParam<Scalar>("SpatialParams.AlphaBeaversJoseph");
-        temperature_ = getParam<Scalar>("SpatialParams.PorousMediumTemperature");
+        temperature_ = getParam<Scalar>("SpatialParams.Temperature");
     }
 
     /*!
diff --git a/exercises/solution/exercise-coupling-ff-pm/interface/freeflowsubproblem.hh b/exercises/solution/exercise-coupling-ff-pm/interface/freeflowsubproblem.hh
index e083d38aece298468b187e93ee673462738703ae..25606fdb90b73880207b6e9d5b0ebbb2404f8951 100644
--- a/exercises/solution/exercise-coupling-ff-pm/interface/freeflowsubproblem.hh
+++ b/exercises/solution/exercise-coupling-ff-pm/interface/freeflowsubproblem.hh
@@ -20,18 +20,17 @@
  * \file
  * \brief The free flow sub problem
  */
-#ifndef DUMUX_FREEFLOW_INTERFACE_SUBPROBLEM_HH
-#define DUMUX_FREEFLOW_INTERFACE_SUBPROBLEM_HH
+#ifndef DUMUX_FREEFLOW_SUBPROBLEM_HH
+#define DUMUX_FREEFLOW_SUBPROBLEM_HH
 
+#include <dumux/freeflow/navierstokes/staggered/problem.hh>
+#include <dumux/freeflow/navierstokes/boundarytypes.hh>
 #include <dumux/common/properties.hh>
 #include <dumux/common/boundarytypes.hh>
-#include <dumux/common/timeloop.hh>
 #include <dumux/common/numeqvector.hh>
 
-#include <dumux/freeflow/navierstokes/staggered/problem.hh>
-#include <dumux/freeflow/navierstokes/boundarytypes.hh>
-
 namespace Dumux {
+
 /*!
  * \brief The free flow sub problem
  */
@@ -61,8 +60,11 @@ class FreeFlowSubProblem : public NavierStokesStaggeredProblem<TypeTag>
     using CouplingManager = GetPropType<TypeTag, Properties::CouplingManager>;
 
 public:
-    FreeFlowSubProblem(std::shared_ptr<const FVGridGeometry> fvGridGeometry, std::shared_ptr<CouplingManager> couplingManager)
-    : ParentType(fvGridGeometry, "Freeflow"), eps_(1e-6), couplingManager_(couplingManager)
+    FreeFlowSubProblem(std::shared_ptr<const FVGridGeometry> fvGridGeometry,
+                       std::shared_ptr<CouplingManager> couplingManager)
+    : ParentType(fvGridGeometry, "Freeflow"),
+    eps_(1e-6),
+    couplingManager_(couplingManager)
     {
         deltaP_ = getParamFromGroup<Scalar>(this->paramGroup(), "Problem.PressureDifference");
     }
@@ -88,6 +90,7 @@ public:
             values.setDirichlet(Indices::velocityYIdx);
         }
 
+        // left/right wall
         if (onRightBoundary_(globalPos) || (onLeftBoundary_(globalPos)))
         {
             values.setDirichlet(Indices::velocityXIdx);
@@ -111,6 +114,7 @@ public:
         }
 #endif
 
+        // coupling interface
         if(couplingManager().isCoupledEntity(CouplingManager::stokesIdx, scvf))
         {
             values.setCouplingNeumann(Indices::conti0EqIdx);
@@ -180,15 +184,6 @@ public:
         return values;
     }
 
-
-    //! Set the coupling manager
-    void setCouplingManager(std::shared_ptr<CouplingManager> cm)
-    { couplingManager_ = cm; }
-
-    //! Get the coupling manager
-    const CouplingManager& couplingManager() const
-    { return *couplingManager_; }
-
     /*!
      * \brief Return the sources within the domain.
      *
@@ -269,6 +264,14 @@ public:
         return analyticalVelocityX_;
     }
 
+    //! Set the coupling manager
+    void setCouplingManager(std::shared_ptr<CouplingManager> cm)
+    { couplingManager_ = cm; }
+
+    //! Get the coupling manager
+    const CouplingManager& couplingManager() const
+    { return *couplingManager_; }
+
     // \}
 
 private:
diff --git a/exercises/solution/exercise-coupling-ff-pm/interface/main.cc b/exercises/solution/exercise-coupling-ff-pm/interface/main.cc
index b4db7eb11730c8da40b60fc4f839535cde3baefe..94536d74b8dc55434d52e01ac9bc4165e54f6472 100644
--- a/exercises/solution/exercise-coupling-ff-pm/interface/main.cc
+++ b/exercises/solution/exercise-coupling-ff-pm/interface/main.cc
@@ -28,6 +28,7 @@
 #include <dumux/common/initialize.hh>
 #include <dumux/common/properties.hh>
 #include <dumux/common/parameters.hh>
+
 #include <dumux/common/partial.hh>
 
 #include <dumux/linear/istlsolvers.hh>
@@ -61,7 +62,7 @@ int main(int argc, char** argv)
 
     // Define the sub problem type tags
     using FreeflowTypeTag = Properties::TTag::FreeflowOneP;
-    using PorousMediumTypeTag = Properties::TTag::PorousMediumOneP;
+    using PorousMediumTypeTag = Properties::TTag::PorousMediumFlowOneP;
 
 #if EXNUMBER < 3
     // try to create a grid (from the given grid file or the input file)
@@ -167,7 +168,6 @@ int main(int argc, char** argv)
     auto porousMediumGridVariables = std::make_shared<PorousMediumGridVariables>(porousMediumProblem, porousMediumFvGridGeometry);
     porousMediumGridVariables->init(sol[porousMediumIdx]);
 
-
     // intialize the vtk output module
     const auto freeflowName = getParam<std::string>("Problem.Name") + "_" + freeflowProblem->name();
     const auto porousMediumName = getParam<std::string>("Problem.Name") + "_" + porousMediumProblem->name();
@@ -181,7 +181,11 @@ int main(int argc, char** argv)
 
     freeflowVtkWriter.write(0.0);
 
-    VtkOutputModule<PorousMediumGridVariables, GetPropType<PorousMediumTypeTag, Properties::SolutionVector>> porousMediumVtkWriter(*porousMediumGridVariables, sol[porousMediumIdx], porousMediumName);
+    using PorousMediumSolutionVector = GetPropType<PorousMediumTypeTag, Properties::SolutionVector>;
+    VtkOutputModule<PorousMediumGridVariables, PorousMediumSolutionVector> porousMediumVtkWriter(*porousMediumGridVariables,
+                                                                                                 sol[porousMediumIdx],
+                                                                                                 porousMediumName);
+
     using PorousMediumVelocityOutput = GetPropType<PorousMediumTypeTag, Properties::VelocityOutput>;
     porousMediumVtkWriter.addVelocityOutput(std::make_shared<PorousMediumVelocityOutput>(*porousMediumGridVariables));
     GetPropType<PorousMediumTypeTag, Properties::IOFields>::initOutputModule(porousMediumVtkWriter);
diff --git a/exercises/solution/exercise-coupling-ff-pm/interface/params.input b/exercises/solution/exercise-coupling-ff-pm/interface/params.input
index 2f1dc08b438a052ce1ef5ff5ab57b34612c2c991..eab4c116545ff15c89fd02e2ea78671d5d31b49a 100644
--- a/exercises/solution/exercise-coupling-ff-pm/interface/params.input
+++ b/exercises/solution/exercise-coupling-ff-pm/interface/params.input
@@ -26,18 +26,17 @@ Cells1 = 20
 Grading1 = 1
 
 [Freeflow.Problem]
-Name = stokes
+Name = freeflow
 PressureDifference = 1e-9
 EnableInertiaTerms = false
 
 [PorousMedium.Problem]
-Name = darcy
+Name = porousmedium
 
 [SpatialParams]
 Permeability = 1e-6 # m^2
 Porosity = 0.4
 AlphaBeaversJoseph = 1.0
-PorousMediumTemperature = 283.15
 Temperature = 283.15
 
 [Problem]
diff --git a/exercises/solution/exercise-coupling-ff-pm/interface/porousmediumsubproblem.hh b/exercises/solution/exercise-coupling-ff-pm/interface/porousmediumsubproblem.hh
index 85e1978c2d7fe20a8785d431c2d2cfce8ddafc31..e7aa0ba37d4515ab449830792403c4a609fb5df1 100644
--- a/exercises/solution/exercise-coupling-ff-pm/interface/porousmediumsubproblem.hh
+++ b/exercises/solution/exercise-coupling-ff-pm/interface/porousmediumsubproblem.hh
@@ -21,8 +21,8 @@
  *
  * \brief The porous medium flow sub problem
  */
-#ifndef DUMUX_POROUSMEDIUMFLOW_INTERFACE_SUBPROBLEM_HH
-#define DUMUX_POROUSMEDIUMFLOW_INTERFACE_SUBPROBLEM_HH
+#ifndef DUMUX_POROUSMEDIUMFLOW_SUBPROBLEM_HH
+#define DUMUX_POROUSMEDIUMFLOW_SUBPROBLEM_HH
 
 #include <dumux/porousmediumflow/problem.hh>
 #include <dumux/common/properties.hh>
@@ -60,8 +60,10 @@ class PorousMediumSubProblem : public PorousMediumFlowProblem<TypeTag>
 
 public:
     PorousMediumSubProblem(std::shared_ptr<const FVGridGeometry> fvGridGeometry,
-                   std::shared_ptr<CouplingManager> couplingManager)
-    : ParentType(fvGridGeometry, "PorousMedium"), eps_(1e-7), couplingManager_(couplingManager)
+                           std::shared_ptr<CouplingManager> couplingManager)
+    : ParentType(fvGridGeometry, "PorousMedium"),
+    eps_(1e-7),
+    couplingManager_(couplingManager)
     {}
 
     /*!
@@ -74,6 +76,8 @@ public:
     BoundaryTypes boundaryTypes(const Element &element, const SubControlVolumeFace &scvf) const
     {
         BoundaryTypes values;
+
+        // set Neumann BCs to all boundaries first
         values.setAllNeumann();
 
 #if EXNUMBER == 0 // flow from top to bottom
@@ -97,6 +101,7 @@ public:
      */
     PrimaryVariables dirichlet(const Element &element, const SubControlVolumeFace &scvf) const
     {
+        // set p = 0 at the bottom
         PrimaryVariables values(0.0);
         values = initial(element);
 
@@ -120,8 +125,10 @@ public:
                         const ElementFluxVariablesCache& elemFluxVarsCache,
                         const SubControlVolumeFace& scvf) const
     {
+        // no-flow everywhere ...
         NumEqVector values(0.0);
 
+        // ... except at the coupling interface
         if (couplingManager().isCoupledEntity(CouplingManager::darcyIdx, scvf))
             values[Indices::conti0EqIdx] = couplingManager().couplingData().massCouplingCondition(element, fvGeometry, elemVolVars, scvf);
 
@@ -153,7 +160,7 @@ public:
      * variables.
      */
     PrimaryVariables initial(const Element &element) const
-    { return PrimaryVariables(0.0);}
+    { return PrimaryVariables(0.0); }
 
     //! Set the coupling manager
     void setCouplingManager(std::shared_ptr<CouplingManager> cm)
diff --git a/exercises/solution/exercise-coupling-ff-pm/interface/properties.hh b/exercises/solution/exercise-coupling-ff-pm/interface/properties.hh
index f8a8587c2f9c1299ae570b97d53ed1755142a8e1..36da2f0dc8ca4e4dbb69228ddaf476148ed38654 100644
--- a/exercises/solution/exercise-coupling-ff-pm/interface/properties.hh
+++ b/exercises/solution/exercise-coupling-ff-pm/interface/properties.hh
@@ -18,13 +18,12 @@
  *****************************************************************************/
 /*!
  * \file
- *
  * \brief The coupled exercise properties file or the interface case.
  */
 #ifndef DUMUX_EXERCISE_COUPLED_INTERFACE_PROPERTIES_HH
 #define DUMUX_EXERCISE_COUPLED_INTERFACE_PROPERTIES_HH
 
-// Both domains
+// Both Domains
 #include <dune/grid/yaspgrid.hh>
 #include <dumux/multidomain/staggeredtraits.hh>
 #include <dumux/multidomain/boundary/stokesdarcy/couplingmanager.hh>
@@ -32,8 +31,15 @@
 #if EXNUMBER >= 3
 #include <dumux/io/grid/gridmanager_sub.hh>
 #endif
-#include <dumux/material/components/simpleh2o.hh>
+
 #include <dumux/material/fluidsystems/1pliquid.hh>
+#include <dumux/material/components/simpleh2o.hh>
+
+// Free-flow domain
+#include <dumux/discretization/staggered/freeflow/properties.hh>
+#include <dumux/freeflow/navierstokes/model.hh>
+
+#include "freeflowsubproblem.hh"
 
 // Porous medium flow domain
 #include <dumux/discretization/cctpfa.hh>
@@ -42,29 +48,23 @@
 #include "../1pspatialparams.hh"
 #include "porousmediumsubproblem.hh"
 
-// Free-flow domain
-#include <dumux/discretization/staggered/freeflow/properties.hh>
-#include <dumux/freeflow/navierstokes/model.hh>
-
-#include "freeflowsubproblem.hh"
-
 namespace Dumux::Properties {
 
 // Create new type tags
 namespace TTag {
 struct FreeflowOneP { using InheritsFrom = std::tuple<NavierStokes, StaggeredFreeFlowModel>; };
-struct PorousMediumOneP { using InheritsFrom = std::tuple<OneP, CCTpfaModel>; };
+struct PorousMediumFlowOneP { using InheritsFrom = std::tuple<OneP, CCTpfaModel>; };
 } // end namespace TTag
 
 // Set the coupling manager
 template<class TypeTag>
 struct CouplingManager<TypeTag, TTag::FreeflowOneP>
 {
-    using Traits = StaggeredMultiDomainTraits<TypeTag, TypeTag, Properties::TTag::PorousMediumOneP>;
+    using Traits = StaggeredMultiDomainTraits<TypeTag, TypeTag, Properties::TTag::PorousMediumFlowOneP>;
     using type = Dumux::StokesDarcyCouplingManager<Traits>;
 };
 template<class TypeTag>
-struct CouplingManager<TypeTag, TTag::PorousMediumOneP>
+struct CouplingManager<TypeTag, TTag::PorousMediumFlowOneP>
 {
     using Traits = StaggeredMultiDomainTraits<Properties::TTag::FreeflowOneP, Properties::TTag::FreeflowOneP, TypeTag>;
     using type = Dumux::StokesDarcyCouplingManager<Traits>;
@@ -72,13 +72,27 @@ struct CouplingManager<TypeTag, TTag::PorousMediumOneP>
 
 // Set the problem property
 template<class TypeTag>
-struct Problem<TypeTag, TTag::PorousMediumOneP> { using type = Dumux::PorousMediumSubProblem<TypeTag>; };
+struct Problem<TypeTag, TTag::PorousMediumFlowOneP> { using type = Dumux::PorousMediumSubProblem<TypeTag>; };
 template<class TypeTag>
 struct Problem<TypeTag, TTag::FreeflowOneP> { using type = Dumux::FreeFlowSubProblem<TypeTag> ; };
 
+// the fluid system
+template<class TypeTag>
+struct FluidSystem<TypeTag, TTag::PorousMediumFlowOneP>
+{
+    using Scalar = GetPropType<TypeTag, Properties::Scalar>;
+    using type = FluidSystems::OnePLiquid<Scalar, Dumux::Components::SimpleH2O<Scalar> > ;
+};
+template<class TypeTag>
+struct FluidSystem<TypeTag, TTag::FreeflowOneP>
+{
+    using Scalar = GetPropType<TypeTag, Properties::Scalar>;
+    using type = FluidSystems::OnePLiquid<Scalar, Dumux::Components::SimpleH2O<Scalar> > ;
+};
+
 // Set the grid type
 template<class TypeTag>
-struct Grid<TypeTag, TTag::PorousMediumOneP>
+struct Grid<TypeTag, TTag::PorousMediumFlowOneP>
 {
     static constexpr auto dim = 2;
     using Scalar = GetPropType<TypeTag, Properties::Scalar>;
@@ -106,22 +120,8 @@ struct Grid<TypeTag, TTag::FreeflowOneP>
 #endif
 };
 
-// the fluid system
-template<class TypeTag>
-struct FluidSystem<TypeTag, TTag::PorousMediumOneP>
-{
-    using Scalar = GetPropType<TypeTag, Properties::Scalar>;
-    using type = FluidSystems::OnePLiquid<Scalar, Dumux::Components::SimpleH2O<Scalar> > ;
-};
-template<class TypeTag>
-struct FluidSystem<TypeTag, TTag::FreeflowOneP>
-{
-    using Scalar = GetPropType<TypeTag, Properties::Scalar>;
-    using type = FluidSystems::OnePLiquid<Scalar, Dumux::Components::SimpleH2O<Scalar> > ;
-};
-
 template<class TypeTag>
-struct SpatialParams<TypeTag, TTag::PorousMediumOneP> {
+struct SpatialParams<TypeTag, TTag::PorousMediumFlowOneP> {
     using type = OnePSpatialParams<GetPropType<TypeTag, GridGeometry>, GetPropType<TypeTag, Scalar>>;
 };
 
@@ -132,6 +132,6 @@ struct EnableGridFluxVariablesCache<TypeTag, TTag::FreeflowOneP> { static conste
 template<class TypeTag>
 struct EnableGridVolumeVariablesCache<TypeTag, TTag::FreeflowOneP> { static constexpr bool value = true; };
 
-} // end namespace Dumux::Properties
+} //end namespace Dumux::Properties
 
 #endif
diff --git a/exercises/solution/exercise-coupling-ff-pm/models/freeflowsubproblem.hh b/exercises/solution/exercise-coupling-ff-pm/models/freeflowsubproblem.hh
index d8ae90cde7d62a9915c673a355dc7731e7301f31..227e276af325321fe70ddf435172519bf84893c0 100644
--- a/exercises/solution/exercise-coupling-ff-pm/models/freeflowsubproblem.hh
+++ b/exercises/solution/exercise-coupling-ff-pm/models/freeflowsubproblem.hh
@@ -21,8 +21,8 @@
  * \ingroup NavierStokesTests
  * \brief A simple Stokes test problem for the staggered grid (Navier-)Stokes model.
  */
-#ifndef DUMUX_FREEFLOW_MODEL_SUBPROBLEM_HH
-#define DUMUX_FREEFLOW_MODEL_SUBPROBLEM_HH
+#ifndef DUMUX_FREEFLOW1P2C_SUBPROBLEM_HH
+#define DUMUX_FREEFLOW1P2C_SUBPROBLEM_HH
 
 #include <dumux/common/properties.hh>
 #include <dumux/common/boundarytypes.hh>
@@ -68,12 +68,12 @@ class FreeFlowSubProblem : public NavierStokesStaggeredProblem<TypeTag>
 
     static constexpr bool useMoles = GetPropType<TypeTag, Properties::ModelTraits>::useMoles();
 
-    static constexpr auto dim = GetPropType<TypeTag, Properties::ModelTraits>::dim();
-    static constexpr auto transportCompIdx = 1;
-
 public:
-    FreeFlowSubProblem(std::shared_ptr<const FVGridGeometry> fvGridGeometry, std::shared_ptr<CouplingManager> couplingManager)
-    : ParentType(fvGridGeometry, "Freeflow"), eps_(1e-6), couplingManager_(couplingManager)
+    FreeFlowSubProblem(std::shared_ptr<const FVGridGeometry> fvGridGeometry,
+                       std::shared_ptr<CouplingManager> couplingManager)
+    : ParentType(fvGridGeometry, "Freeflow"),
+    eps_(1e-6),
+    couplingManager_(couplingManager)
     {
         velocity_ = getParamFromGroup<Scalar>(this->paramGroup(), "Problem.Velocity");
         pressure_ = getParamFromGroup<Scalar>(this->paramGroup(), "Problem.Pressure");
@@ -171,29 +171,29 @@ public:
      *
      * \param globalPos The global position
      */
-     PrimaryVariables initialAtPos(const GlobalPosition& globalPos) const
-     {
-         PrimaryVariables values(0.0);
-
-         // This is only an approximation of the actual hydorostatic pressure gradient.
-         // Air is compressible (the density depends on pressure), thus the actual
-         // vertical pressure profile is non-linear.
-         // This discrepancy can lead to spurious flows at the outlet if the inlet
-         // velocity is chosen too small.
-         FluidState fluidState;
-         updateFluidStateForBC_(fluidState, pressure_);
-         const Scalar density = FluidSystem::density(fluidState, 0);
-         values[Indices::pressureIdx] = pressure_ - density*this->gravity()[1]*(this->gridGeometry().bBoxMax()[1] - globalPos[1]);
-
-
-         // gravity has negative sign
-         values[Indices::conti0EqIdx + 1] = moleFraction_;
-         values[Indices::velocityXIdx] = 4.0 * velocity_ * (globalPos[1] - this->gridGeometry().bBoxMin()[1])
-                                                         * (this->gridGeometry().bBoxMax()[1] - globalPos[1])
-                                                         / (height_() * height_());
-
-         return values;
-     }
+    PrimaryVariables initialAtPos(const GlobalPosition& globalPos) const
+    {
+        PrimaryVariables values(0.0);
+
+        // This is only an approximation of the actual hydorostatic pressure gradient.
+        // Air is compressible (the density depends on pressure), thus the actual
+        // vertical pressure profile is non-linear.
+        // This discrepancy can lead to spurious flows at the outlet if the inlet
+        // velocity is chosen too small.
+        FluidState fluidState;
+        updateFluidStateForBC_(fluidState, pressure_);
+        const Scalar density = FluidSystem::density(fluidState, 0);
+        values[Indices::pressureIdx] = pressure_ - density*this->gravity()[1]*(this->gridGeometry().bBoxMax()[1] - globalPos[1]);
+
+
+        // gravity has negative sign
+        values[Indices::conti0EqIdx + 1] = moleFraction_;
+        values[Indices::velocityXIdx] = 4.0 * velocity_ * (globalPos[1] - this->gridGeometry().bBoxMin()[1])
+                                                        * (this->gridGeometry().bBoxMax()[1] - globalPos[1])
+                                                        / (height_() * height_());
+
+        return values;
+    }
 
     /*!
      * \brief Return the sources within the domain.
@@ -203,9 +203,6 @@ public:
     NumEqVector sourceAtPos(const GlobalPosition &globalPos) const
     { return NumEqVector(0.0); }
 
-    void setTimeLoop(TimeLoopPtr timeLoop)
-    { timeLoop_ = timeLoop; }
-
     /*!
      * \brief Returns the intrinsic permeability of required as input parameter for the Beavers-Joseph-Saffman boundary condition
      */
@@ -230,6 +227,9 @@ public:
     const CouplingManager& couplingManager() const
     { return *couplingManager_; }
 
+    void setTimeLoop(TimeLoopPtr timeLoop)
+    { timeLoop_ = timeLoop; }
+
 private:
     bool onLeftBoundary_(const GlobalPosition &globalPos) const
     { return globalPos[0] < this->gridGeometry().bBoxMin()[0] + eps_; }
@@ -264,6 +264,7 @@ private:
         const Scalar enthalpy = FluidSystem::enthalpy(fluidState, paramCache, 0);
         fluidState.setEnthalpy(0, enthalpy);
     }
+
     // the height of the free-flow domain
     const Scalar height_() const
     { return this->gridGeometry().bBoxMax()[1] - this->gridGeometry().bBoxMin()[1]; }
diff --git a/exercises/solution/exercise-coupling-ff-pm/models/main.cc b/exercises/solution/exercise-coupling-ff-pm/models/main.cc
index 6416f2b7e9dc7432d1d0f4926a92f1d791892fec..759563c785d15b1605539dd3edd6875b874a41a4 100644
--- a/exercises/solution/exercise-coupling-ff-pm/models/main.cc
+++ b/exercises/solution/exercise-coupling-ff-pm/models/main.cc
@@ -19,7 +19,7 @@
 /*!
  * \file
  *
- * \brief A test problem for the coupled Stokes/Darcy problem (1p)
+ * \brief A test problem for the coupled Freeflow/PorousMedium problem (1p)
  */
 #include <config.h>
 
@@ -153,7 +153,8 @@ int main(int argc, char** argv)
     GetPropType<FreeflowTypeTag, Properties::IOFields>::initOutputModule(freeflowVtkWriter);
     freeflowVtkWriter.write(0.0);
 
-    VtkOutputModule<PorousMediumGridVariables, GetPropType<PorousMediumTypeTag, Properties::SolutionVector>> porousMediumVtkWriter(*porousMediumGridVariables, sol[porousMediumIdx], porousMediumName);
+    VtkOutputModule<PorousMediumGridVariables, GetPropType<PorousMediumTypeTag, Properties::SolutionVector>> porousMediumVtkWriter(*porousMediumGridVariables,
+                                                                                                                                   sol[porousMediumIdx], porousMediumName);
     using PorousMediumVelocityOutput = GetPropType<PorousMediumTypeTag, Properties::VelocityOutput>;
     porousMediumVtkWriter.addVelocityOutput(std::make_shared<PorousMediumVelocityOutput>(*porousMediumGridVariables));
     GetPropType<PorousMediumTypeTag, Properties::IOFields>::initOutputModule(porousMediumVtkWriter);
diff --git a/exercises/solution/exercise-coupling-ff-pm/models/params.input b/exercises/solution/exercise-coupling-ff-pm/models/params.input
index dc2a61319d61d8739eb10201b4446330dff26a86..54e2032126b497dc23725d058006cd512f362592 100644
--- a/exercises/solution/exercise-coupling-ff-pm/models/params.input
+++ b/exercises/solution/exercise-coupling-ff-pm/models/params.input
@@ -38,7 +38,6 @@ Swr = 0.005
 Snr = 0.01
 PcLowSw = Swr * 5.0
 PcHighSw = 1.0 - Snr * 5.0
-PorousMediumTemperature = 293.15
 Temperature = 293.15
 
 [Problem]
@@ -49,7 +48,11 @@ ExportFluxes = true
 
 [Newton]
 MaxSteps = 12
-MaxRelativeShift = 1e-5
+MaxRelativeShift = 1e-7
+TargetSteps = 5
 
 [Vtk]
 AddVelocity = 1
+
+[Assembly]
+NumericDifference.BaseEpsilon = 1e-8
diff --git a/exercises/solution/exercise-coupling-ff-pm/models/porousmediumsubproblem.hh b/exercises/solution/exercise-coupling-ff-pm/models/porousmediumsubproblem.hh
index 5f1b26ab2661ac2888af422eb5aa004e693af0a2..458c9db98a5511d2fa40d311ad64ca916be3b25e 100644
--- a/exercises/solution/exercise-coupling-ff-pm/models/porousmediumsubproblem.hh
+++ b/exercises/solution/exercise-coupling-ff-pm/models/porousmediumsubproblem.hh
@@ -34,6 +34,7 @@
 #include <dumux/porousmediumflow/problem.hh>
 
 namespace Dumux {
+
 /*!
  * \brief The porous medium flow sub problem
  */
@@ -56,7 +57,6 @@ class PorousMediumSubProblem : public PorousMediumFlowProblem<TypeTag>
 
     // copy some indices for convenience
     using Indices = typename GetPropType<TypeTag, Properties::ModelTraits>::Indices;
-
     // grid and world dimension
     static constexpr int dim = GridView::dimension;
     static constexpr int dimworld = GridView::dimensionworld;
diff --git a/exercises/solution/exercise-coupling-ff-pm/models/properties.hh b/exercises/solution/exercise-coupling-ff-pm/models/properties.hh
index ca80958014b76f08cce119f27791ab30af5d2838..5416e7641043d97788c7c10a57796f7e39a4aaea 100644
--- a/exercises/solution/exercise-coupling-ff-pm/models/properties.hh
+++ b/exercises/solution/exercise-coupling-ff-pm/models/properties.hh
@@ -19,19 +19,19 @@
 /*!
  * \file
  *
- * \brief The coupled exercise properties file or the models case.
+ * \brief Properties file for coupled models example.
  */
 #ifndef DUMUX_EXERCISE_COUPLED_MODELS_PROPERTIES_HH
 #define DUMUX_EXERCISE_COUPLED_MODELS_PROPERTIES_HH
 
-// Both domains
+// Both Domains
 #include <dune/grid/yaspgrid.hh>
+#include <dumux/multidomain/staggeredtraits.hh>
+#include <dumux/multidomain/boundary/stokesdarcy/couplingmanager.hh>
 
 #include <dumux/io/gnuplotinterface.hh>
 #include <dumux/material/fluidsystems/1padapter.hh>
 #include <dumux/material/fluidsystems/h2oair.hh>
-#include <dumux/multidomain/staggeredtraits.hh>
-#include <dumux/multidomain/boundary/stokesdarcy/couplingmanager.hh>
 
 // Porous medium flow domain
 #include <dumux/discretization/cctpfa.hh>
@@ -47,7 +47,7 @@
 
 #include "porousmediumsubproblem.hh"
 
-// Free flow domain
+// Free-flow
 #include <dumux/discretization/staggered/freeflow/properties.hh>
 #include <dumux/freeflow/compositional/navierstokesncmodel.hh>
 
@@ -81,24 +81,12 @@ struct CouplingManager<TypeTag, TTag::PorousMediumOnePNC>
 
 // Set the problem property
 template<class TypeTag>
-struct Problem<TypeTag, TTag::FreeflowNC> { using type = Dumux::FreeFlowSubProblem<TypeTag> ; };
-template<class TypeTag>
 struct Problem<TypeTag, TTag::PorousMediumOnePNC> { using type = Dumux::PorousMediumSubProblem<TypeTag>; };
-
-// Set the grid type
-template<class TypeTag>
-struct Grid<TypeTag, TTag::FreeflowNC> { using type = Dune::YaspGrid<2, Dune::EquidistantOffsetCoordinates<GetPropType<TypeTag, Properties::Scalar>, 2> >; };
 template<class TypeTag>
-struct Grid<TypeTag, TTag::PorousMediumOnePNC> { using type = Dune::YaspGrid<2>; };
+struct Problem<TypeTag, TTag::FreeflowNC> { using type = Dumux::FreeFlowSubProblem<TypeTag> ; };
 
 // The fluid system
 template<class TypeTag>
-struct FluidSystem<TypeTag, TTag::FreeflowNC>
-{
-    using H2OAir = FluidSystems::H2OAir<GetPropType<TypeTag, Properties::Scalar>>;
-    using type = FluidSystems::OnePAdapter<H2OAir, H2OAir::gasPhaseIdx>;
-};
-template<class TypeTag>
 struct FluidSystem<TypeTag, TTag::PorousMediumOnePNC>
 {
     using H2OAir = FluidSystems::H2OAir<GetPropType<TypeTag, Properties::Scalar>>;
@@ -108,18 +96,35 @@ struct FluidSystem<TypeTag, TTag::PorousMediumOnePNC>
     using type = H2OAir;
 #endif
 };
+template<class TypeTag>
+struct FluidSystem<TypeTag, TTag::FreeflowNC>
+{
+    using H2OAir = FluidSystems::H2OAir<GetPropType<TypeTag, Properties::Scalar>>;
+    using type = FluidSystems::OnePAdapter<H2OAir, H2OAir::gasPhaseIdx>;
+};
 
-// Do not replace one equation with a total mass balance
+// Use moles
 template<class TypeTag>
-struct ReplaceCompEqIdx<TypeTag, TTag::FreeflowNC> { static constexpr int value = 3; };
+struct UseMoles<TypeTag, TTag::PorousMediumOnePNC> { static constexpr bool value = true; };
+template<class TypeTag>
+struct UseMoles<TypeTag, TTag::FreeflowNC> { static constexpr bool value = true; };
+
+// Do not replace one equation with a total mass balance
 template<class TypeTag>
 struct ReplaceCompEqIdx<TypeTag, TTag::PorousMediumOnePNC> { static constexpr int value = 3; };
+template<class TypeTag>
+struct ReplaceCompEqIdx<TypeTag, TTag::FreeflowNC> { static constexpr int value = 3; };
 
-// Use formulation based on mass fractions
+// Set the grid type
 template<class TypeTag>
-struct UseMoles<TypeTag, TTag::FreeflowNC> { static constexpr bool value = true; };
+struct Grid<TypeTag, TTag::PorousMediumOnePNC> { using type = Dune::YaspGrid<2>; };
 template<class TypeTag>
-struct UseMoles<TypeTag, TTag::PorousMediumOnePNC> { static constexpr bool value = true; };
+struct Grid<TypeTag, TTag::FreeflowNC> { using type = Dune::YaspGrid<2, Dune::EquidistantOffsetCoordinates<GetPropType<TypeTag, Properties::Scalar>, 2> >; };
+
+//! Use a model with constant tortuosity for the effective diffusivity
+template<class TypeTag>
+struct EffectiveDiffusivityModel<TypeTag, TTag::PorousMediumOnePNC>
+{ using type = DiffusivityConstantTortuosity<GetPropType<TypeTag, Properties::Scalar>>; };
 
 #if EXNUMBER >= 1
 //! Set the default formulation to pw-Sn: This can be over written in the problem.
@@ -141,11 +146,6 @@ struct SpatialParams<TypeTag, TTag::PorousMediumOnePNC> {
 };
 #endif
 
-//! Use a model with constant tortuosity for the effective diffusivity
-template<class TypeTag>
-struct EffectiveDiffusivityModel<TypeTag, TTag::PorousMediumOnePNC>
-{ using type = DiffusivityConstantTortuosity<GetPropType<TypeTag, Properties::Scalar>>; };
-
 template<class TypeTag>
 struct EnableGridGeometryCache<TypeTag, TTag::FreeflowNC> { static constexpr bool value = true; };
 template<class TypeTag>
@@ -153,6 +153,6 @@ struct EnableGridFluxVariablesCache<TypeTag, TTag::FreeflowNC> { static constexp
 template<class TypeTag>
 struct EnableGridVolumeVariablesCache<TypeTag, TTag::FreeflowNC> { static constexpr bool value = true; };
 
-} //end namespace Dumux::Properties
+} // end namespace Dumux::Properties
 
 #endif
diff --git a/exercises/solution/exercise-coupling-ff-pm/turbulence/freeflowsubproblem.hh b/exercises/solution/exercise-coupling-ff-pm/turbulence/freeflowsubproblem.hh
index 66dbbc3b192757f538d136af30c090d0496fc083..11042b533bd30eee26438880cc87db821bbc9773 100644
--- a/exercises/solution/exercise-coupling-ff-pm/turbulence/freeflowsubproblem.hh
+++ b/exercises/solution/exercise-coupling-ff-pm/turbulence/freeflowsubproblem.hh
@@ -24,9 +24,9 @@
 #define DUMUX_FREEFLOW_TURBULENCE_SUBPROBLEM_SOL_HH
 
 
-#include <dumux/common/timeloop.hh>
 #include <dumux/common/properties.hh>
 #include <dumux/common/boundarytypes.hh>
+#include <dumux/common/timeloop.hh>
 #include <dumux/common/numeqvector.hh>
 #include <dumux/multidomain/boundary/stokesdarcy/couplingdata.hh>
 
@@ -41,7 +41,6 @@
 #endif
 
 namespace Dumux {
-
 /*!
  * \brief The free-flow sub problem
  */
@@ -101,7 +100,8 @@ public:
         refTemperature_ = getParamFromGroup<Scalar>(this->paramGroup(), "Problem.RefTemperature");
 
         diffCoeffAvgType_ = StokesDarcyCouplingOptions::stringToEnum(DiffusionCoefficientAveragingType{},
-                                                                     getParamFromGroup<std::string>(this->paramGroup(), "Problem.InterfaceDiffusionCoefficientAvg"));
+                                                                     getParamFromGroup<std::string>(this->paramGroup(),
+                                                                     "Problem.InterfaceDiffusionCoefficientAvg"));
 #if EXNUMBER >=1
         FluidSystem::init();
         Dumux::TurbulenceProperties<Scalar, dimWorld, true> turbulenceProperties;
@@ -109,7 +109,7 @@ public:
         const auto phaseIdx = 0;
         fluidState.setPressure(phaseIdx, refPressure_);
         fluidState.setTemperature(this->spatialParams().temperatureAtPos({}));
-        fluidState.setMassFraction(phaseIdx, phaseIdx, refMoleFrac_);
+        fluidState.setMassFraction(phaseIdx, phaseIdx, 1.0);
         Scalar density = FluidSystem::density(fluidState, phaseIdx);
         Scalar kinematicViscosity = FluidSystem::viscosity(fluidState, phaseIdx) / density;
         Scalar diameter = this->gridGeometry().bBoxMax()[1] - this->gridGeometry().bBoxMin()[1];
@@ -202,12 +202,11 @@ public:
             values.setCouplingNeumann(Indices::momentumYBalanceIdx);
             values.setBeaversJoseph(Indices::momentumXBalanceIdx);
         }
-
         return values;
     }
 
     /*!
-     * \brief Evaluate the boundary conditions for a dirichlet values at the boundary.
+     * \brief Evaluate the boundary conditions for dirichlet values at the boundary.
      *
      * \param element The finite element
      * \param scvf the sub control volume face
@@ -217,6 +216,20 @@ public:
     {
         const auto globalPos = scvf.ipGlobal();
         PrimaryVariables values(initialAtPos(globalPos));
+
+#if EXNUMBER == 1
+        if (onUpperBoundary_(globalPos) || onLowerBoundary_(globalPos))
+        {
+            values[Indices::turbulentKineticEnergyIdx] = 0.0;
+            values[Indices::dissipationIdx] = 0.0;
+        }
+#elif EXNUMBER  >= 2
+        if (onLowerBoundary_(globalPos))
+        {
+            values[Indices::turbulentKineticEnergyIdx] = 0.0;
+            values[Indices::dissipationIdx] = 0.0;
+        }
+#endif
         return values;
     }
 
@@ -328,20 +341,6 @@ public:
         values[Indices::dissipationIdx] = dissipation_;
 #endif
 
-#if EXNUMBER == 1
-        if (onLowerBoundary_(globalPos) || onUpperBoundary_(globalPos))
-        {
-            values[Indices::turbulentKineticEnergyIdx] = 0.0;
-            values[Indices::dissipationIdx] = 0.0;
-        }
-#elif EXNUMBER >= 2
-        if (onLowerBoundary_(globalPos) || onUpperBoundary_(globalPos))
-        {
-            values[Indices::turbulentKineticEnergyIdx] = 0.0;
-            values[Indices::dissipationIdx] = 0.0;
-        }
-#endif
-
 #if EXNUMBER >= 2
         if(onLowerBoundary_(globalPos))
             values[Indices::velocityXIdx] = 0.0;
diff --git a/exercises/solution/exercise-coupling-ff-pm/turbulence/main.cc b/exercises/solution/exercise-coupling-ff-pm/turbulence/main.cc
index ad0470ad0a2755d83bee78e4ef1fc25a43fecac4..15009b12f81bd032e9366b9a51baaf5caf9b06f0 100644
--- a/exercises/solution/exercise-coupling-ff-pm/turbulence/main.cc
+++ b/exercises/solution/exercise-coupling-ff-pm/turbulence/main.cc
@@ -63,7 +63,7 @@ int main(int argc, char** argv)
 
     // Define the sub problem type tags
     using FreeflowTypeTag = Properties::TTag::FreeflowModel;
-    using PorousMediumTypeTag = Properties::TTag::PorousMediumModel;
+    using PorousMediumTypeTag = Properties::TTag::PorousMediumFlowModel;
 
     // try to create a grid (from the given grid file or the input file)
     // for both sub-domains
@@ -156,7 +156,8 @@ int main(int argc, char** argv)
     GetPropType<FreeflowTypeTag, Properties::IOFields>::initOutputModule(freeflowVtkWriter);
     freeflowVtkWriter.write(0.0);
 
-    VtkOutputModule<PorousMediumGridVariables, GetPropType<PorousMediumTypeTag, Properties::SolutionVector>> porousMediumVtkWriter(*porousMediumGridVariables, sol[porousMediumIdx], porousMediumName);
+    using PorousMediumSolutionVector = GetPropType<PorousMediumTypeTag, Properties::SolutionVector>;
+    VtkOutputModule<PorousMediumGridVariables, PorousMediumSolutionVector> porousMediumVtkWriter(*porousMediumGridVariables, sol[porousMediumIdx], porousMediumName);
     using PorousMediumVelocityOutput = GetPropType<PorousMediumTypeTag, Properties::VelocityOutput>;
     porousMediumVtkWriter.addVelocityOutput(std::make_shared<PorousMediumVelocityOutput>(*porousMediumGridVariables));
     GetPropType<PorousMediumTypeTag, Properties::IOFields>::initOutputModule(porousMediumVtkWriter);
diff --git a/exercises/solution/exercise-coupling-ff-pm/turbulence/params.input b/exercises/solution/exercise-coupling-ff-pm/turbulence/params.input
index 0769cc10d4f3977f7a3b0bd6c25bcac56c122cc9..a4d7be6eb5bb45d2f97bf9a47c5679e5c11876ff 100644
--- a/exercises/solution/exercise-coupling-ff-pm/turbulence/params.input
+++ b/exercises/solution/exercise-coupling-ff-pm/turbulence/params.input
@@ -46,7 +46,7 @@ Swr = 0.005
 Snr = 0.01
 PcLowSw = Swr * 5.0
 PcHighSw = 1.0 - Snr * 5.0
-PorousMediumTemperature = 298.15 # [K]
+Temperature = 298.15 # [K]
 
 [Problem]
 Name = ex_coupling_turbulence_ff-pm
diff --git a/exercises/solution/exercise-coupling-ff-pm/turbulence/porousmediumsubproblem.hh b/exercises/solution/exercise-coupling-ff-pm/turbulence/porousmediumsubproblem.hh
index f7fd1fd835a6f3cebae701e5f09dc8ff578618c2..bc2bd75a63131f1ac416822e15b1ce2e5858ec7a 100644
--- a/exercises/solution/exercise-coupling-ff-pm/turbulence/porousmediumsubproblem.hh
+++ b/exercises/solution/exercise-coupling-ff-pm/turbulence/porousmediumsubproblem.hh
@@ -24,15 +24,15 @@
 #ifndef DUMUX_POROUSMEDIUMFLOW2P2C_SUBPROBLEM_HH
 #define DUMUX_POROUSMEDIUMFLOW2P2C_SUBPROBLEM_HH
 
-#include <dumux/common/properties.hh>
+#include <dumux/porousmediumflow/problem.hh>
 #include <dumux/common/boundarytypes.hh>
+#include <dumux/common/properties.hh>
 #include <dumux/common/timeloop.hh>
 #include <dumux/common/numeqvector.hh>
-
-#include <dumux/porousmediumflow/problem.hh>
 #include <dumux/multidomain/boundary/stokesdarcy/couplingdata.hh>
 
 namespace Dumux {
+
 /*!
  * \brief The porous medium sub problem
  */
@@ -91,7 +91,6 @@ public:
                                                                      getParamFromGroup<std::string>(this->paramGroup(), "Problem.InterfaceDiffusionCoefficientAvg"));
     }
 
-
     template<class SolutionVector, class GridVariables>
     void postTimeStep(const SolutionVector& curSol,
                       const GridVariables& gridVariables,
@@ -172,7 +171,7 @@ public:
                         const FVElementGeometry& fvGeometry,
                         const ElementVolumeVariables& elemVolVars,
                         const ElementFluxVariablesCache& elemFluxVarsCache,
-                        const SubControlVolumeFace& scvf)  const
+                        const SubControlVolumeFace& scvf) const
     {
         NumEqVector values(0.0);
 
@@ -267,6 +266,7 @@ private:
     std::shared_ptr<CouplingManager> couplingManager_;
     DiffusionCoefficientAveragingType diffCoeffAvgType_;
 };
-} //end namespace Dumux
+
+} //end namespace dUMUX
 
 #endif
diff --git a/exercises/solution/exercise-coupling-ff-pm/turbulence/properties.hh b/exercises/solution/exercise-coupling-ff-pm/turbulence/properties.hh
index 3b6b962a2ced962347b31f66e892739f3bafa869..2b9e6518b325b2d1d62c617bea6aed9da8ec4679 100644
--- a/exercises/solution/exercise-coupling-ff-pm/turbulence/properties.hh
+++ b/exercises/solution/exercise-coupling-ff-pm/turbulence/properties.hh
@@ -18,7 +18,6 @@
  *****************************************************************************/
 /*!
  * \file
- *
  * \brief The coupled exercise properties file or the turbulent case.
  */
 #ifndef DUMUX_EXERCISE_COUPLING_TURBULENCE_PROPERTIES_HH
@@ -26,18 +25,18 @@
 
 // Both domains
 #include <dune/grid/yaspgrid.hh>
-
-#include <dumux/material/fluidsystems/1padapter.hh>
-#include <dumux/material/fluidsystems/h2oair.hh>
 #include <dumux/multidomain/staggeredtraits.hh>
 #include <dumux/multidomain/boundary/stokesdarcy/couplingmanager.hh>
 
+#include <dumux/material/fluidsystems/h2oair.hh>
+#include <dumux/material/fluidsystems/1padapter.hh>
+
 // Porous medium flow domain
 #include <dumux/porousmediumflow/2p2c/model.hh>
 #include <dumux/discretization/cctpfa.hh>
 
-#include "../2pspatialparams.hh"
 #include"porousmediumsubproblem.hh"
+#include "../2pspatialparams.hh"
 
 // Free-flow domain
 #include <dumux/discretization/staggered/freeflow/properties.hh>
@@ -54,7 +53,7 @@ namespace Dumux::Properties {
 
 // Create new type tags
 namespace TTag {
-struct PorousMediumModel { using InheritsFrom = std::tuple<TwoPTwoCNI, CCTpfaModel>; };
+struct PorousMediumFlowModel { using InheritsFrom = std::tuple<TwoPTwoCNI, CCTpfaModel>; };
 #if EXNUMBER >= 1
 struct FreeflowModel { using InheritsFrom = std::tuple<SSTNCNI, StaggeredFreeFlowModel>; };
 #else
@@ -66,11 +65,11 @@ struct FreeflowModel { using InheritsFrom = std::tuple<NavierStokesNCNI, Stagger
 template<class TypeTag>
 struct CouplingManager<TypeTag, TTag::FreeflowModel>
 {
-    using Traits = StaggeredMultiDomainTraits<TypeTag, TypeTag, Properties::TTag::PorousMediumModel>;
+    using Traits = StaggeredMultiDomainTraits<TypeTag, TypeTag, Properties::TTag::PorousMediumFlowModel>;
     using type = Dumux::StokesDarcyCouplingManager<Traits>;
 };
 template<class TypeTag>
-struct CouplingManager<TypeTag, TTag::PorousMediumModel>
+struct CouplingManager<TypeTag, TTag::PorousMediumFlowModel>
 {
     using Traits = StaggeredMultiDomainTraits<Properties::TTag::FreeflowModel, Properties::TTag::FreeflowModel, TypeTag>;
     using type = Dumux::StokesDarcyCouplingManager<Traits>;
@@ -78,17 +77,19 @@ struct CouplingManager<TypeTag, TTag::PorousMediumModel>
 
 // Set the problem property
 template<class TypeTag>
-struct Problem<TypeTag, TTag::FreeflowModel> { using type = Dumux::FreeFlowSubProblem<TypeTag> ; };
+struct Problem<TypeTag, TTag::PorousMediumFlowModel> { using type = Dumux::PorousMediumSubProblem<TypeTag>; };
 template<class TypeTag>
-struct Problem<TypeTag, TTag::PorousMediumModel> { using type = Dumux::PorousMediumSubProblem<TypeTag>; };
+struct Problem<TypeTag, TTag::FreeflowModel> { using type = Dumux::FreeFlowSubProblem<TypeTag> ; };
 
 // Set the grid type
 template<class TypeTag>
-struct Grid<TypeTag, TTag::FreeflowModel> { using type = Dune::YaspGrid<2, Dune::TensorProductCoordinates<GetPropType<TypeTag, Properties::Scalar>, 2> >; };
+struct Grid<TypeTag, TTag::PorousMediumFlowModel> { using type = Dune::YaspGrid<2, Dune::TensorProductCoordinates<GetPropType<TypeTag, Properties::Scalar>, 2> >; };
 template<class TypeTag>
-struct Grid<TypeTag, TTag::PorousMediumModel> { using type = Dune::YaspGrid<2, Dune::TensorProductCoordinates<GetPropType<TypeTag, Properties::Scalar>, 2> >; };
+struct Grid<TypeTag, TTag::FreeflowModel> { using type = Dune::YaspGrid<2, Dune::TensorProductCoordinates<GetPropType<TypeTag, Properties::Scalar>, 2> >; };
 
-// The fluid system
+// the fluid system
+template<class TypeTag>
+struct FluidSystem<TypeTag, TTag::PorousMediumFlowModel> { using type = FluidSystems::H2OAir<GetPropType<TypeTag, Properties::Scalar>>; };
 template<class TypeTag>
 struct FluidSystem<TypeTag, TTag::FreeflowModel>
 {
@@ -96,28 +97,27 @@ struct FluidSystem<TypeTag, TTag::FreeflowModel>
   static constexpr auto phaseIdx = H2OAir::gasPhaseIdx; // simulate the air phase
   using type = FluidSystems::OnePAdapter<H2OAir, phaseIdx>;
 };
+
+// Use formulation based on mole fractions
+template<class TypeTag>
+struct UseMoles<TypeTag, TTag::PorousMediumFlowModel> { static constexpr bool value = true; };
 template<class TypeTag>
-struct FluidSystem<TypeTag, TTag::PorousMediumModel> { using type = FluidSystems::H2OAir<GetPropType<TypeTag, Properties::Scalar>>; };
+struct UseMoles<TypeTag, TTag::FreeflowModel> { static constexpr bool value = true; };
 
 // The gas component balance (air) is replaced by the total mass balance
 template<class TypeTag>
-struct ReplaceCompEqIdx<TypeTag, TTag::PorousMediumModel> { static constexpr int value = 3; };
+struct ReplaceCompEqIdx<TypeTag, TTag::PorousMediumFlowModel> { static constexpr int value = 3; };
 template<class TypeTag>
 struct ReplaceCompEqIdx<TypeTag, TTag::FreeflowModel> { static constexpr int value = 3; };
 
-// Use formulation based on mass fractions
-template<class TypeTag>
-struct UseMoles<TypeTag, TTag::FreeflowModel> { static constexpr bool value = true; };
-template<class TypeTag>
-struct UseMoles<TypeTag, TTag::PorousMediumModel> { static constexpr bool value = true; };
 
 //! Set the default formulation to pw-Sn: This can be over written in the problem.
 template<class TypeTag>
-struct Formulation<TypeTag, TTag::PorousMediumModel>
+struct Formulation<TypeTag, TTag::PorousMediumFlowModel>
 { static constexpr auto value = TwoPFormulation::p1s0; };
 
 template<class TypeTag>
-struct SpatialParams<TypeTag, TTag::PorousMediumModel>
+struct SpatialParams<TypeTag, TTag::PorousMediumFlowModel>
 { using type = TwoPSpatialParams<GetPropType<TypeTag, GridGeometry>, GetPropType<TypeTag, Scalar>>; };
 
 template<class TypeTag>
@@ -127,6 +127,6 @@ struct EnableGridFluxVariablesCache<TypeTag, TTag::FreeflowModel> { static const
 template<class TypeTag>
 struct EnableGridVolumeVariablesCache<TypeTag, TTag::FreeflowModel> { static constexpr bool value = true; };
 
-} //end namespace Dumux::Properties
+} // end namespace Dumux::Properties
 
 #endif