Skip to content
Snippets Groups Projects
Commit 7ede763e authored by Martin Schneider's avatar Martin Schneider
Browse files

[doc][ex][ffpm] Update description in readme

parent 2e77fe11
No related branches found
No related tags found
1 merge request!213Feature/updates ff pm
......@@ -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
......@@ -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,
......@@ -386,7 +384,7 @@ 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,
In the problem constructor, uncomment the code calculating 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.
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment