OilfieldToolsNet/opm-simulators
4
0
Fork 0
mirror of https://github.com/OPM/opm-simulators.git synced 2025-02-25 18:55:30 -06:00
Code Issues Packages Projects Releases Wiki Activity

update the tutorial for the coupled models

Browse Source
This commit is contained in:
Andreas Lauser
2012-02-13 18:44:10 +00:00
committed by Andreas Lauser
parent 7e14c1964f
commit 25aba6a3f2
5 changed files with 337 additions and 362 deletions
Download Patch File Download Diff File Expand all files Collapse all files

479
doc/handbook/tutorial-coupled.tex
View File

@@ -1,4 +1,4 @@
\section[Fully-coupled model]{Solving a problem using a Fully-Coupled Model}\label{tutorial-coupled} \section[Fully-Implicit Model]{Solving a Problem Using a Fully-Coupled Model}\label{tutorial-coupled}
The process of setting up a problem using \Dumux can be roughly divided into four parts: The process of setting up a problem using \Dumux can be roughly divided into four parts:
\begin{enumerate} \begin{enumerate}
@@ -49,7 +49,7 @@ The solved equations are the mass balances of water and oil:
0 0
\end{align} \end{align}
\subsection{The main file} \subsection{The Main File}
Listing \ref{tutorial-coupled:mainfile} shows the main application file Listing \ref{tutorial-coupled:mainfile} shows the main application file
\texttt{tutorial/tutorial\_coupled.cc} for the coupled two-phase \texttt{tutorial/tutorial\_coupled.cc} for the coupled two-phase
@@ -62,101 +62,91 @@ above.
\end{lst} \end{lst}
From line \ref{tutorial-coupled:include-begin} to line From line \ref{tutorial-coupled:include-begin} to line
\ref{tutorial-coupled:include-end} the \Dune and \Dumux files, which \ref{tutorial-coupled:include-end} the headers required are included.
contain some required functions and classes, and the problem file are included.
At line \ref{tutorial-coupled:set-type-tag} the type tag of the At line \ref{tutorial-coupled:set-type-tag} the type tag of the
problem, which is going to be simulated, is set. All other data types problem, which is going to be simulated, is specified. All other data
can be retrieved via the \Dumux property system and only depend on this types can be retrieved via the \Dumux property system and only depend
single type tag. Retrieving them is done between line on this single type tag. For a more thourough introduction to the
\ref{tutorial-coupled:retrieve-types-begin} and \Dumux property system, see chapter~\ref{sec:propertysystem}.
\ref{tutorial-coupled:retrieve-types-end}. For an introduction to the
property system, see section \ref{sec:propertysystem}.
The first thing which should be done at run time is to initialize the After this \Dumux' default startup routine is called on line
message passing interface using \Dune's \texttt{MPIHelper} class. Line \ref{tutorial-coupled:call-start}. To provide a error message, the
\ref{tutorial-coupled:init-mpi} is essential if the simulation is usage message which is displayed to the user if the simulation is
intended to be run on more than one processor at the same time using MPI. Next, called incorrectly, is printed via the custom function which is defined
the command-line arguments are parsed starting at line on line~\ref{tutorial-coupled:usage-function}.
\ref{tutorial-coupled:parse-args-begin} down to line
\ref{tutorial-coupled:parse-args-end}. In this example, it is checked if and
from which time on a previous run of the simulation should be restarted. Furthermore, we
parse the time when the simulation ends and the initial time-step size.
After this, a grid is created in line \subsection{The Problem Class}
\ref{tutorial-coupled:create-grid} and the problem is instantiated for
its leaf grid view in line \ref{tutorial-coupled:instantiate-problem}.
Finally, on line \ref{tutorial-coupled:initTimeManager} the time
manager is created with the parsed starting parameters. If requested by
the user, a state written to disk by a previous simulation run can be
restored if via the restart flag.
The simulation procedure is started in line
\ref{tutorial-coupled:execute}.
\subsection{The problem class}
When solving a problem using \Dumux, the most important file is the When solving a problem using \Dumux, the most important file is the
so-called \textit{problem file} as shown in listing so-called \textit{problem file} as shown in
\ref{tutorial-coupled:problemfile} of listing~\ref{tutorial-coupled:problemfile}.
\texttt{tutorialproblem\_coupled.hh}.
\begin{lst}[File tutorial/tutorialproblem\_coupled.hh]\label{tutorial-coupled:problemfile} \mbox{} \begin{lst}[File tutorial/tutorialproblem\_coupled.hh]\label{tutorial-coupled:problemfile} \mbox{}
\lstinputlisting[basicstyle=\ttfamily\scriptsize,numbers=left, \lstinputlisting[basicstyle=\ttfamily\scriptsize,numbers=left,
numberstyle=\tiny, numbersep=5pt, firstline=27]{../../tutorial/tutorialproblem_coupled.hh} numberstyle=\tiny, numbersep=5pt, firstline=28]{../../tutorial/tutorialproblem_coupled.hh}
\end{lst} \end{lst}
First, a new type tag is created for the problem in line First, a new type tag is created for the problem in line
\ref{tutorial-coupled:create-type-tag}. In this case, the new type \ref{tutorial-coupled:create-type-tag}. In this case, the new type
tag inherits all properties from the \texttt{BoxTwoP} type tag, tag inherits all properties from the \texttt{BoxTwoP} type tag, which
which means that for this problem the two-phase box model is chosen as means that for this problem the two-phase box model is chosen as
discretization scheme. Further, it inherits from the spatial parameters type tag, discretization scheme. Further, it inherits from the spatial
which is defined in the problem-dependent spatial parameters file (line \ref{tutorial-coupled:define-spatialparameters-typetag}). parameters type tag, which is defined in the problem-dependent spatial
On line \ref{tutorial-coupled:set-problem}, a parameters file (line
problem class is attached to the new type tag, while the grid which \ref{tutorial-coupled:define-spatialparameters-typetag}). On line
is going to be used is defined in line \ref{tutorial-coupled:set-grid} -- \ref{tutorial-coupled:set-problem}, a problem class is attached to the
in this case that is \texttt{SGrid}. Since there's no uniform new type tag, while the grid which is going to be used is defined in
mechanism to allocate grids in \Dune, the \texttt{Grid} property also contains line \ref{tutorial-coupled:set-grid} -- in this case that is
a static \texttt{create()} method which provides just that. The \texttt{SGrid} uses three variables of \texttt{Dune::YaspGrid}. Since there's no uniform mechanism to
Type \texttt{Dune::FieldVector} to define the lower left corner of the domain allocate grids in \Dune, \Dumux features the concept of grid creators.
(\texttt{L}), the upper right corner of the domain (\texttt{H}) and the number In this case the generic \texttt{CubeGridCreator} which creates a
of cells in $x$ and $y$ direction (\texttt{N}). structztured hexahedron grid of a specified size and resolution. For
this grid creator physical domain of the grid is specified via the
run-time parameters \texttt{Grid.upperRightX},
\texttt{Grid.upperRightY}, \texttt{Grid.numberOfCellsX} and
\texttt{Grid.numberOfCellsY}. These parameters can be specified via
the command-line or in a parameter file,
see~\ref{tutorial-coupled:runtime-parameters}.
Next, the appropriate fluid system, which specifies both, information about Next, the appropriate fluid system, which specifies the thermodynamic
the fluid mixture as well as about the pure substances, has to be chosen. relations of the fluid phases, has to be chosen. By default, the
The two-phase model uses by default the \texttt{FluidSystem2P}, which assumes two-phase model uses the \texttt{TwoPImmiscibleFluidSystem}, which
immiscibility of the phases, but requires the wetting and non-wetting phases assumes immiscibility of the phases, but requires that the components
to be explicitly set. In this case, liquid water which uses the relations from used for the wetting and non-wetting phases to be explicitly set. In
this case, liquid water which uses the relations from
IAPWS'97~\cite{IAPWS1997} is chosen as the wetting phase on line IAPWS'97~\cite{IAPWS1997} is chosen as the wetting phase on line
\ref{tutorial-coupled:wettingPhase} and liquid oil is chosen as the \ref{tutorial-coupled:wettingPhase} and liquid oil is chosen as the
non-wetting phase on line \ref{tutorial-coupled:nonwettingPhase}. non-wetting phase on line \ref{tutorial-coupled:nonwettingPhase}. The
The last property, which is set in line last property, which is set in line \ref{tutorial-coupled:gravity},
\ref{tutorial-coupled:gravity}, is optional and tells the model not to tells the model not to use gravity.
use gravity.
Parameters which are specific to a set-up, such as boundary and initial Parameters which are specific to a physical set-up to be simulated,
conditions, source terms or temperature within the domain, and which are such as boundary and initial conditions, source terms or temperature
required to solve the differential equations of the models are within the domain, and which are required to solve the differential
specified via a \textit{problem} class. If the two-phase box model is equations of the models are specified via a \textit{problem} class. If
used, this class should be derived from \texttt{TwoPBoxProblem} as done the two-phase box model is used, this class should be derived from
in line \ref{tutorial-coupled:def-problem}. \texttt{TwoPBoxProblem} as done in line
\ref{tutorial-coupled:def-problem}.
The problem class always has at least five methods: The problem class always has at least five methods:
\begin{itemize} \begin{itemize}
\item A method \texttt{boundaryTypes()} specifying the type of \item A method \texttt{boundaryTypes()} specifying the type of
boundary conditions at each vertex. boundary conditions at each vertex.
\item A method \texttt{dirichlet()} specifying the actual values for \item A method \texttt{dirichlet()} specifying the actual values for
the Dirichlet conditions at each Dirichlet vertex. the \textsc{Dirichlet} conditions at each \textsc{Dirichlet} vertex.
\item A method \texttt{neumann()} specifying the actual values for \item A method \texttt{neumann()} specifying the actual values for
the Neumann conditions, which are usually evaluated at the the \textsc{Neumann} conditions, which are usually evaluated at the
integration points of the Neumann boundary faces. integration points of the \textsc{Neumann} boundary faces.
\item A method for source or sink terms called \texttt{source()}, usually evaluated at \item A method for source or sink terms called \texttt{source()}, usually evaluated at
the center of a sub-control volume. the center of a control volume.
\item A method called \texttt{initial()} for specifying the initial \item A method called \texttt{initial()} for specifying the initial
conditions at each vertex. conditions at each vertex.
\end{itemize} \end{itemize}
For the definition of the the boundary condition types and of the values of the Dirichlet boundaries, For the definition of the the boundary condition types and of the
two parameters are required: values of the \textsc{Dirichlet} boundaries, two parameters are
available:
\begin{description} \begin{description}
\item [values:] A vector which stores the result of the method. What \item [values:] A vector which stores the result of the method. What
the values in this vector mean is dependent on the method: For the values in this vector mean is dependent on the method: For
@@ -165,19 +155,23 @@ two parameters are required:
condition types. It has as many entries as the model has primary variables / equations. condition types. It has as many entries as the model has primary variables / equations.
For the typical case, in which all equations have the same boundary For the typical case, in which all equations have the same boundary
condition type at a certain position, there are two methods that set the appropriate conditions condition type at a certain position, there are two methods that set the appropriate conditions
for all primary variables / equations: Either \texttt{setAllDirichlet()} or \texttt{setAllNeumann()}. for all primary variables / equations: \texttt{setAllDirichlet()} and \texttt{setAllNeumann()}.
\item [vertex:] The boundary condition and the Dirichlet values are specified at the vertex, which represents a \item [vertex:] The boundary condition and the Dirichlet values are
sub-control volume. This avoids the specification of two different boundary condition types for one equation at one sub-control volume. specified for a vertex, which represents a control volume in the box
Be aware that the second parameter is a Dune grid entity with the codimension dim. discretization. This avoids the specification of two different
boundary condition types for one equation at different control
volume. Be aware that the second parameter is a Dune grid entity
with the codimension \texttt{dim}.
\end{description} \end{description}
To ensure that no boundaries are undefined, a small safeguard value \texttt{eps\_} To ensure that no boundaries are undefined, a small safeguard value
is usually added when comparing spatial coordinates. The left boundary is \texttt{eps\_} is usually added when comparing spatial
hence not assigned where the first coordinate is equal to zero, but where it is coordinates. The left boundary is hence detected by comparing the
first coordinate is equal with zero, but by testing whether it is
smaller than a very small value \texttt{eps\_}. smaller than a very small value \texttt{eps\_}.
Methods which make statements about boundary segments of the grid (i.e. Methods which make statements about boundary segments of the grid
\texttt{neumann()}) get six parameters: (i.e. \texttt{neumann()}) get called with six arguments:
\begin{description} \begin{description}
\item[values:] A vector \texttt{neumann()}, in which the mass fluxes per area unit \item[values:] A vector \texttt{neumann()}, in which the mass fluxes per area unit
over the boundary segment are specified. over the boundary segment are specified.
@@ -193,115 +187,150 @@ Methods which make statements about boundary segments of the grid (i.e.
\end{description} \end{description}
Similarly, the \texttt{initial()} and \texttt{source()} methods Similarly, the \texttt{initial()} and \texttt{source()} methods
specify properties of sub-control volumes and thus only get specify properties of control volumes and thus only get
\texttt{values}, \texttt{element}, \texttt{fvElemGeom} and \texttt{values}, \texttt{element}, \texttt{fvElemGeom} and
\texttt{scvIdx} as parameters. \texttt{scvIdx} as arguments.
In addition to these five methods, there might be some model-specific In addition to these five methods, there might be some model-specific
methods. If the isothermal two-phase model is used, this includes methods. If the isothermal two-phase model is used, this includes for
for example a \texttt{temperature()} method which returns the temperature in Kelvin example a \texttt{temperature()} method which returns the temperature
of the fluids and the rock matrix in the domain. This temperature is in \textsc{Kelvin} of the fluids and the rock matrix in the
then used by the model to calculate fluid properties which possibly domain. This temperature is then used by the model to calculate fluid
depend on it, e.g. density. The \texttt{bboxMax()} method that is used here is a vector properties which possibly depend on it, e.g. density. The
providing information about the spatial extends of the simulation domain. This method \texttt{bboxMax()} (``\textbf{max}imum coordinated of the grid's
and the respective \texttt{bboxMin()} method \textbf{b}ounding \textbf{b}ox'') method that is used here to
can be found in the base class \texttt{Dumux::BoxProblem<TypeTag>}. determine the extend of the phsical domain, returns a vector with the
maximum values of each coordinate of the grid's physical. This method
and the analogous \texttt{bboxMin()} method are provided by the base
class \texttt{Dumux::BoxProblem<TypeTag>}.
\subsection{Defining Fluid Properties}\label{tutorial-coupled:description-fluid-class}
\subsection{Defining fluid properties}\label{tutorial-coupled:description-fluid-class} The \Dumux distribution includes some common substances which can be
used out of the box. The properties of the pure substances (such as
the components nitrogen, water, or the pseudo-component air) are
provided by header files located in the folder
\verb+dumux/material/components+.
The \Dumux distribution includes some common substances which can be used Most often, when two or more components are considered, fluid
out of the box. The properties of the pure substances (such as the components interactions such as solubility effects come into play and properties
nitrogen, water, or the pseudo-component air) are stored in header files in of mixtures such as density or enthalpy are of interest. These
the folder \verb+dumux/material/components+. Each of these files interactions are defined by {\em fluid systems}, which are located in
defines a class with the same name as the component but starting with a capital \verb+dumux/material/fluidsystems+. A more thorough overview of the
letter, e.g. \texttt{Water}, and are derived from \texttt{Component}. \Dumux fluid framework can be found in chapter~\ref{sec:fluidframework}.
Most often, when two or more components are considered, fluid interactions
such as solubility effects come into play and properties of mixtures such as
density or enthalpy are of interest. These interactions are defined in
a specific \verb+fluidsystem+ in the folder \verb+dumux/material/fluidsystems+.
It features methods returning fluid properties like density, enthalpy, viscosity,
etc. by accessing the pure components as well as binary coefficients such as
Henry or diffusion coefficients, which are stored in
\verb+dumux/material/binarycoefficients+. New fluids which are not yet
available in the \Dumux distribution can be defined analogously.
% In this example, a class for the definition of a two-phase system is used. This allows for the choice % In this example, a class for the definition of a two-phase system is used. This allows for the choice
% of the two components oil and water and for access of the parameters that are relevant for the two-phase model. % of the two components oil and water and for access of the parameters that are relevant for the two-phase model.
\subsection{The definition of the parameters that are dependent on space}\label{tutorial-coupled:description-spatialParameters} \subsection{Definiting Spatially Dependent Parameters}\label{tutorial-coupled:description-spatialParameters}
In \Dumux, the properties of the porous medium such as the \textit{intrinsic In \Dumux, many properties of the porous medium can depend on the
permeability}, the \textit{porosity}, the \textit{heat capacity} as spatial location. Such properties are the \textit{intrinsic
well as the \textit{heat conductivity} can be defined in space using a permeability}, the parameters of the \textit{capillary pressure} and
so-called \texttt{spatial parameters} class. However, because the soil the \textit{relative permeability}, the \textit{porosity}, the
also has an effect on the material laws of the fluids (e.g. \textit{capillarity}), \textit{heat capacity} as well as the \textit{heat conductivity}. Such
their selection and the definition of their attributes (e.g. \textit{residual parameters are define using a so-called \textit{spatial parameters}
saturations} or \textit{van Genuchten parameters}) are also accomplished in the spatial parameters. class.
The base class \texttt{Dumux::BoxSpatialParameters<TypeTag>} contains a general The the box discretization is to be used, the spatial paramters class
averaging procedure for vertex-centered box-methods. should be derived from the base class
\texttt{Dumux::BoxSpatialParameters<TypeTag>}. Listing
Listing \ref{tutorial-coupled:spatialparametersfile} shows the file \ref{tutorial-coupled:spatialparametersfile} shows the file
\verb+tutorialspatialparameters_coupled.hh+: \verb+tutorialspatialparameters_coupled.hh+:
\begin{lst}[File tutorial/tutorialspatialparameters\_coupled.hh]\label{tutorial-coupled:spatialparametersfile} \mbox{} \begin{lst}[File tutorial/tutorialspatialparameters\_coupled.hh]\label{tutorial-coupled:spatialparametersfile} \mbox{}
\lstinputlisting[basicstyle=\ttfamily\scriptsize,numbers=left, \lstinputlisting[basicstyle=\ttfamily\scriptsize,numbers=left,
numberstyle=\tiny, numbersep=5pt, firstline=27]{../../tutorial/tutorialspatialparameters_coupled.hh} numberstyle=\tiny, numbersep=5pt, firstline=28]{../../tutorial/tutorialspatialparameters_coupled.hh}
\end{lst} \end{lst}
First, the spatial parameters type tag has to be created (line \ref{tutorial-coupled:define-spatialparameters-typetag}), First, the spatial parameters type tag is created on line
from which the problem type tag inherits its properties and which can be classified as spatial parameter properties \ref{tutorial-coupled:define-spatialparameters-typetag}. The type tag
(problem file line \ref{tutorial-coupled:create-type-tag}). These are e.g. the spatial parameters class itself (line \ref{tutorial-coupled:set-spatialparameters}) for the problem then derives from it. The \Dumux properties defined on
or a certain material law (line \ref{tutorial-coupled:rawlaw} \label{tutorial-coupled:materialLaw}). the type tag for the spatial parameters are for example, the spatial
\Dumux provides several material laws in the folder parameters class itself (line
\verb+dumux/material/fluidmatrixinteractions+. \ref{tutorial-coupled:set-spatialparameters}) or the capillary
The selected one -- here it is a relation according to a regularized version of Brooks \& Corey -- is included pressure/relative permability relations\footnote{Taken together, the
in line \ref{tutorial-coupled:rawLawInclude}. After the selection, capillary pressure and the relative permeability relations are
an adapter in line \ref{tutorial-coupled:eff2abs} translates between effective saturations, which are employed in the Brooks \& Corey parameterization called \textit{material law}.} which ought to be used by the
and which deduce the residual saturations, and the simulated saturations. simulation (line
As the applied raw law knows best which kind of parameters are necessary, \ref{tutorial-coupled:rawlaw} \label{tutorial-coupled:materialLaw}).
it provides a parameter class \texttt{RegularizedBrooksCoreyParams} that is \Dumux provides several material laws in the folder
accessible via the member \texttt{Params} and defined in line \verb+dumux/material/fluidmatrixinteractions+. The selected one --
\ref{tutorial-coupled:matLawObjectType}. The material law object here it is a relation according to a regularized version of
is now instantiated correctly as a private object \textsc{Brooks} \& \textsc{Corey} -- is included in line
in line \ref{tutorial-coupled:matParamsObject}. \ref{tutorial-coupled:rawLawInclude}. After the selection, an adapter
class is specified ob line \ref{tutorial-coupled:eff2abs} to
translates between effective and absolute saturations. This way,
residual saturations can be specified in a generic way. As only used
material law knows which the names of the parameters which it
requires, it provides a parameter class
\texttt{RegularizedBrooksCoreyParams} which is exported by the type
\texttt{Params} and defined in line
\ref{tutorial-coupled:matLawObjectType}. In this case, the spatial
parameters only require a single set of parameters which means that it
only requires a single material parameter object as can be seen on
line~\ref{tutorial-coupled:matParamsObject}.
In line \ref{tutorial-coupled:permeability} a method returning the In line \ref{tutorial-coupled:permeability}, a method returning the
intrinsic permeability can be found. As can be seen, the method has intrinsic permeability is specified. As can be seen, the method has
to be called with three arguments: to be called with three arguments:
(\texttt{Element}) is again the considered element, which also holds information \begin{description}
about its geometry and position, the second argument \item[\texttt{element}:] Just like for the problem itself, this
(\texttt{fvElemGeom}) holds information about the finite-volume geometry induced parameter describes the considered element by means of a \Dune
by the box-method, and the third defines the local index of the current sub-control entity. Elements provide information about its geometry and
volume. The intrinsic permeability is a tensor and is thus returned in form of position and can be mapped to a global index.
a $\texttt{dim} \times \texttt{dim}$-matrix where \texttt{dim} is the dimension \item[\texttt{fvElemGeom}:] Holds information about the finite-volume
of the problem. geometry induced by the box-method on the element.
\item[\texttt{scvIdx}:] Is the index of the sub-control volume of the
element which is considered. This is equivalent to the local index
of the vertex which corrosponts to the considered control volume in
the element.
\end{description}
The intrinsic permeability is a tensor and is thus the method returnes
a $\texttt{dim} \times \texttt{dim}$-matrix where \texttt{dim} is the
dimension of the grid.
The method \texttt{porosity()} defined in line The method \texttt{porosity()} defined in line
\ref{tutorial-coupled:porosity} is called with the same arguments as \ref{tutorial-coupled:porosity} is called with the same arguments as
the permeability function described before and returns the porosity \texttt{intrinsicPermeability()} and returns a scalar value for
dependent on the position in the domain. porosity dependent on the position in the domain.
Next, the method \texttt{materialLawParams()} defines in line Next, the method \texttt{materialLawParams()}, defined on line
\ref{tutorial-coupled:matLawParams} which \verb+materialLawParams+ object \ref{tutorial-coupled:matLawParams}, specifies
should be applied at this specific position. Although in this case only one object is returned, \verb+materialLawParams+ object applies at the specified
in general the problem may be heterogeneous, returning different objects at different positions in space. position. Although in this case only one object is returned, in
While the selection of the type of this object was already explained (line \ref{tutorial-coupled:rawLawInclude}), general, the problem may be heterogeneous, which necessitates
some specific parameter values of the applied material law, such as the Brooks \& Corey parameters, are still needed. This is returning different objects at different positions in space. While
done in the constructor body (line \ref{tutorial-coupled:setLawParams}). the selection of the type of this object was already explained (line
Depending on the type of the \texttt{materialLaw} object, the adequate \texttt{set}-methods \ref{tutorial-coupled:rawLawInclude}), some specific parameter values
are provided by the object to access all necessary parameters of the used material law, such as the \textsc{Brooks} \&
for the applied material law. The name of the access / set functions as well as the rest of the implementation \textsc{Corey} parameters, are still needed. This is done in the
of the material description can be found in constructor at line \ref{tutorial-coupled:setLawParams}. Depending on
the type of the \texttt{materialLaw} object, the \texttt{set}-methods
might be different than those given in this example. The name of the
access / set functions as well as the rest of the implementation of
the material description can be found in
\verb+dumux/material/fluidmatrixinteractions/2p+. \verb+dumux/material/fluidmatrixinteractions/2p+.
\subsection{The definition of run-time parameters}\label{tutorial-coupled:runtime-parameters}
Some parameters need to be specified at runtime. These can either be
specified directly via command line optioms or the can be collected in
a parameter file. The parameter file which \Dumux uses by default has
the name of the program with ``.input'' appended. The parameter file
for the coupled tutorial is thus named \verb+tutorial_coupled.input+
and defaults to the content shown in listing~\ref{tutorial-coupled:parameter-file}.
\begin{lst}[File tutorial/tutorial\_coupled.input]\label{tutorial-coupled:parameter-file} \mbox{}
\lstinputlisting[basicstyle=\ttfamily\scriptsize,numbers=left,numberstyle=\tiny]{../../tutorial/tutorial_coupled.input}
\end{lst}
\subsection{Exercises} \subsection{Exercises}
\label{tutorial-coupled:exercises} \label{tutorial-coupled:exercises}
The following exercises will give you the opportunity to learn how you The following exercises will give you the opportunity to learn how you
can change soil parameters, boundary conditions and fluid properties can change soil parameters, boundary conditions, run-time parameters
in \Dumux. and fluid properties in \Dumux.
\subsubsection{Exercise 1} \subsubsection{Exercise 1}
\renewcommand{\labelenumi}{\alph{enumi})} For Exercise 1 you have \renewcommand{\labelenumi}{\alph{enumi})} For Exercise 1 you have
@@ -309,8 +338,8 @@ to make only some small changes in the tutorial files.
\begin{enumerate} \begin{enumerate}
\item \textbf{Run the Model} \\ \item \textbf{Runing the Program} \\
To get an impression what the results should look like you can first run the original version of the coupled tutorial model by typing \texttt{./tutorial\_coupled 5e5 10}. The first number behind the simulation name defines the time span of the simulation run in seconds, the second number defines the initial time-step size. Note that the time-step size is automatically optimized during the simulation. For the visualization with paraview please refer to \ref{quick-start-guide}.\\ To get an impression what the results should look like you can first run the original version of the coupled tutorial model by typing \texttt{./tutorial\_coupled}. Note, that the time-step size is automatically adapted during the simulation. For the visualization of the results using paraview please refer to \ref{quick-start-guide}.\\
\item \textbf{Changing the Model Domain and the Boundary Conditions} \\ \item \textbf{Changing the Model Domain and the Boundary Conditions} \\
Change the size of the model domain so that you get a rectangle with Change the size of the model domain so that you get a rectangle with
@@ -326,20 +355,19 @@ To get an impression what the results should look like you can first run the ori
Compile the main file by typing \texttt{make tutorial\_coupled} and Compile the main file by typing \texttt{make tutorial\_coupled} and
run the model as explained above. run the model as explained above.
\item \textbf{Changing Fluids} \\ \item \textbf{Changing Fluids} \\
Now you can change the fluids. Use DNAPL instead of Oil and Brine instead of Water. To do that you have to select different components via the property system in the problem file: Now you can change the fluids. Use DNAPL instead of Oil and Brine instead of Water. To do that, you have to select different components via the property system in the problem file:
\begin{enumerate} \begin{enumerate}
\item Brine: The class \texttt{Dumux::Brine} acts as an adapter to the fluid system that alters a pure water class by adding some salt. Hence, the class \texttt{Dumux::Brine} uses a pure water class, such as \texttt{Dumux::H2O}, as a second template argument after the data type \texttt{<Scalar>} as a template argument (be sure to use the complete water class with its own template parameter). \item Brine: Brine is thermodynamically very similar to pure water but also considers a fixed amount of salt in the liquid phase. Hence, the class \texttt{Dumux::Brine} uses a pure water class, such as \texttt{Dumux::H2O}, as a second template argument after the data type \texttt{<Scalar>} as a template argument.
\item DNAPL: A standard set of chemical substances, such as Oil and Brine, is already included (via a list of \texttt{\#include ..} commandos) and hence easily accessible by default. This is not the case for the class \texttt{Dumux::SimpleDNAPL}, however, which is located in the folder \texttt{dumux/material/components/}. Try to include the file and select the component via the property system. \item DNAPL: A standard set of chemical substances, such as Oil and Brine, is already included in the problem (via a list of \texttt{\#include ..} statements) and hence easily accessible by default. However, this is not the case for the class \texttt{Dumux::SimpleDNAPL}, which describes a simple \textbf{d}ense \textbf{n}on-\textbf{a}queous \textbf{p}hase \textbf{l}iquid and is located in the folder \texttt{dumux/material/components/}. Try to include the file and select the component as the non-wetting phase via the property system.
\end{enumerate} \end{enumerate}
If you want to take a closer look on how the fluid classes are defined and which substances are already available please browse through the files in the directory If you want to take a closer look on how the fluid classes are defined and which substances are already available please browse through the files in the directory
\texttt{/dumux/material/components}. \texttt{/dumux/material/components} and read chapter~\ref{sec:fluidframework}.
\item \textbf{Use the \Dumux fluid system} \\ \item \textbf{Use a full-fledged fluid system} \\
\Dumux usually organizes fluid mixtures via a \texttt{fluidsystem}, see also chapter \ref{sec:fluidframework}. In order to include a fluidsystem you first have to comment the lines \ref{tutorial-coupled:2p-system-start} to \ref{tutorial-coupled:2p-system-end} in the problem file. If you use eclipse, this can easily be done by pressing \textit{str + shift + 7} -- the same as to cancel the comment later on.\\ \Dumux usually describes fluid mixtures via \textit{fluid systems}, see also chapter \ref{sec:fluidframework}. In order to include a fluid system, you first have to comment out lines \ref{tutorial-coupled:2p-system-start} to \ref{tutorial-coupled:2p-system-end} in the problem file. If you use eclipse, this can easily be done by pressing \textit{Ctrl + Shift + 7} -- the same as to cancel the comment later on.\\
Now include the file \texttt{fluidsystems/h2oairsystem.hh} in the material folder, and set a property \texttt{FluidSystem} with the appropriate type, \texttt{Dumux::H2OAirFluidSystem<TypeTag>}. However, this rather complicated fluidsystem uses tabularized fluid data, which need to be initialized (i.e. the tables need to be filled with values) in the constructor body of the current problem by adding \texttt{GET\_PROP\_TYPE(TypeTag, FluidSystem)::init();}. As water flow replacing a gas is much faster, test your simulation only until 2e3 seconds and start with a time step of 1 second.\\ Now include the file \texttt{fluidsystems/h2oairsystem.hh} in the material folder, and set a property \texttt{FluidSystem} with the appropriate type, i.e. \texttt{Dumux::H2OAirFluidSystem<TypeTag>}. However, this is a rather complicated fluid system considers mixtures of components and also uses tabulated components that need to be initialized -- i.e. the tables need to be filled with values. Initializating the fluid system is normally done in the constructor of the problem by calling \texttt{GET\_PROP\_TYPE(TypeTag, FluidSystem)::init();}. As water flow replacing a gas is much faster, test your simulation only until $2000$ seconds and start with a time step of $1$ second.\\
Please reverse the changes of this example, as we still use bulk phases and hence do not need such an extensive fluid system. Please reverse the changes made in this part of the exercise, as we will continue to use immiscible phases from here on and hence do not need a complex fluid system.
\item \textbf{Changing Constitutive Relations} \\ \item \textbf{Changing Constitutive Relations} \\
Use an unregularized linear law with an entry pressure of $p_e = 0.0$ and maximal capillary pressure of e.g. $p_{c_{max}} = 2000.0$ instead of using a Use an unregularized linear law with an entry pressure of $p_e = 0.0$ and maximal capillary pressure of e.g. $p_{c_{max}} = 2000.0$ instead of using a
@@ -368,37 +396,43 @@ of the linear law and the respective \texttt{set}-functions can be found
\caption{Exercise 1f: Set-up of a model domain with a heterogeneity. $\Delta \text{x} = 20$ m $\Delta \text{y} = 20$ m.}\label{tutorial-coupled:exercise1_d} \caption{Exercise 1f: Set-up of a model domain with a heterogeneity. $\Delta \text{x} = 20$ m $\Delta \text{y} = 20$ m.}\label{tutorial-coupled:exercise1_d}
\end{figure} \end{figure}
domain. You can use the fluids of exercise 1c).\\ domain. You can use the fluids of exercise 1c).\\
Hint: The current position of the element can be obtained via \texttt{element.geometry().center();}.\\ \textbf{Hint:} The current position of the control volume can be obtained via \texttt{element.geometry().corner(scvIdx);}.\\
When does the front cross the material border? In paraview, the option \textit{View} $\rightarrow$ \textit{Animation View} is nice to get a rough feeling of the time-step sizes. When does the front cross the material border? In paraview, the
animation inspector (\textit{View} $\rightarrow$ \textit{Animation
View}) is a convenient way to get a rough feeling of the time-step
sizes.
\end{enumerate} \end{enumerate}
\subsubsection{Exercise 2} \subsubsection{Exercise 2}
For this exercise you should create a new problem file analogous to For this exercise you should create a new problem file analogous to
the file \texttt{tutorialproblem\_coupled.hh} (e.g. with the name the file \texttt{tutorialproblem\_coupled.hh} (e.g. with the name
\texttt{ex2\_tutorialproblem\_coupled.hh} and new spatial parameters \texttt{ex2\_tutorialproblem\_coupled.hh} and new spatial parameters
just like \texttt{tutorialspatialparameters\_coupled.hh}. The new problem file needs to just like \texttt{tutorialspatialparameters\_coupled.hh}. The new
problem file needs to
be included in the file \texttt{tutorial\_coupled.cc}.\\ be included in the file \texttt{tutorial\_coupled.cc}.\\
The new files should contain the definition of new classes with The new files should contain the definition of new classes with names
names that relate to the file name, such as that relate to the file name, such as
\texttt{Ex2TutorialProblemCoupled}. Make sure that you also adjust the guardian \texttt{Ex2TutorialProblemCoupled}. Make sure that you also adjust the
macros in lines \ref{tutorial-coupled:guardian1} and \ref{tutorial-coupled:guardian1} guardian macros in lines \ref{tutorial-coupled:guardian1} and
in the header files (e.g. change \\ \ref{tutorial-coupled:guardian1}
in the header files (e.g. change \\
\texttt{DUMUX\_TUTORIALPROBLEM\_COUPLED\_HH} to \texttt{DUMUX\_TUTORIALPROBLEM\_COUPLED\_HH} to
\texttt{DUMUX\_EX2\_TUTORIALPROBLEM\_COUPLED\_HH}). Besides also adjusting the guardian macros, \texttt{DUMUX\_EX2\_TUTORIALPROBLEM\_COUPLED\_HH}). Besides also
the new problem file should define and use a new type tag for the problem as well as a new problem class adjusting the guardian macros, the new problem file should define and
e.g. \texttt{Ex2TutorialProblemCoupled}. Make sure to assign your newly defined spatial use a new type tag for the problem as well as a new problem class
parameter class to the \texttt{SpatialParameters} property for the new e.g. \texttt{Ex2TutorialProblemCoupled}. Make sure to assign your
newly defined spatial parameter class to the
\texttt{SpatialParameters} property for the new
type tag. \\ type tag. \\
After this, change the \texttt{create()} method of the \texttt{Grid} After this, change the run-time parameters so that they match the
property so that it matches the domain described domain described by figure \ref{tutorial-coupled:ex2_Domain}. Adapt
by figure \ref{tutorial-coupled:ex2_Domain}. Adapt the problem class the problem class so that the boundary conditions are consistent with
so that the boundary conditions are consistent with figure figure \ref{tutorial-coupled:ex2_BC}. Initially, the domain is fully
\ref{tutorial-coupled:ex2_BC}. Initially the domain is fully saturated saturated with water and the pressure is $p_w = 5 \times
with water and the pressure is $p_w = 5 \times 10^5 \text{Pa}$. Oil 10^5\;\text{Pa}$. Oil infiltrates from the left side. Create a grid
infiltrates from the left side. Create a grid with $20$ cells in with $20$ cells in $x$-direction and $10$ cells in $y$-direction. The
$x$-direction and $10$ cells in $y$-direction. The simulation time simulation time should be set to $10^6\;\text{s}$ with an
should be set to $1\times 10^6 \text{ s}$ with an initial time step size of initial time step size of $100\;\text{s}$.
$100 \text{ s}$.
Now include your new problem file in the main file and replace the Now include your new problem file in the main file and replace the
\texttt{TutorialProblemCoupled} type tag by the one you've created and \texttt{TutorialProblemCoupled} type tag by the one you've created and
@@ -406,30 +440,30 @@ compile the program.
\begin{figure}[ht] \begin{figure}[ht]
\psfrag{K1}{K $= 10^{-7}\text{ m}^2$} \psfrag{K1}{K $= 10^{-7}\;\text{m}^2$}
\psfrag{phi1}{$\phi = 0.2$} \psfrag{phi1}{$\phi = 0.2$}
\psfrag{Lin}{Brooks-Corey Law} \psfrag{Lin}{Brooks-Corey Law}
\psfrag{Lin2}{$\lambda = 1.8$, $p_e = 1000$} \psfrag{Lin2}{$\lambda = 1.8$, $p_e = 1000\;\text{Pa}$}
\psfrag{K2}{K $= 10^{-9}\text{ m}^2$} \psfrag{K2}{K $= 10^{-9}\;\text{m}^2$}
\psfrag{phi2}{$\phi = 0.15$} \psfrag{phi2}{$\phi = 0.15$}
\psfrag{BC1}{Brooks-Corey Law} \psfrag{BC1}{\textsc{Brooks}-\textsc{Corey} Law}
\psfrag{BC2}{$\lambda = 2$, $p_e = 1500$} \psfrag{BC2}{$\lambda = 2$, $p_e = 1500\;\text{Pa}$}
\psfrag{H1y}{50 m} \psfrag{H1y}{$50\;\text{m}$}
\psfrag{H2y}{15 m} \psfrag{H2y}{$15\;\text{m}$}
\psfrag{H3y}{20 m} \psfrag{H3y}{$20\;\text{m}$}
\psfrag{L1x}{100 m} \psfrag{L1x}{$100\;\text{m}$}
\psfrag{L2x}{50 m} \psfrag{L2x}{$50\;\text{m}$}
\psfrag{L3x}{25 m} \psfrag{L3x}{$25\;\text{m}$}
\centering \centering
\includegraphics[width=0.8\linewidth,keepaspectratio]{EPS/Ex2_Domain.eps} \includegraphics[width=0.8\linewidth,keepaspectratio]{EPS/Ex2_Domain.eps}
\caption{Set-up of the model domain and the soil parameters}\label{tutorial-coupled:ex2_Domain} \caption{Set-up of the model domain and the soil parameters}\label{tutorial-coupled:ex2_Domain}
\end{figure} \end{figure}
\begin{figure}[ht] \begin{figure}[ht]
\psfrag{pw}{$p_w = 5 \times 10^5$ [\text{Pa}]} \psfrag{pw}{$p_w = 5 \times 10^5\;\text{Pa}$}
\psfrag{S}{$S_n = 1.0$} \psfrag{S}{$S_n = 1.0$}
\psfrag{qw}{$q_w = 2 \times 10^{-4}$ [kg/$\text{m}^2$s]} \psfrag{qw}{$q_w = 2 \times 10^{-4} \;\text{kg}/\text{m}^2\text{s}$}
\psfrag{qo}{$q_n = 0.0$ [kg/$\text{m}^2$s]} \psfrag{qo}{$q_n = 0.0 \;\text{kg}/\text{m}^2\text{s}$}
\psfrag{no flow}{no flow} \psfrag{no flow}{no flow}
\centering \centering
\includegraphics[width=0.8\linewidth,keepaspectratio]{EPS/Ex2_Boundary.eps} \includegraphics[width=0.8\linewidth,keepaspectratio]{EPS/Ex2_Boundary.eps}
@@ -437,24 +471,33 @@ compile the program.
\end{figure} \end{figure}
\begin{itemize} \begin{itemize}
\item Increase the simulation time to e.g. $4\times 10^7 \text{ s}$. Investigate the saturation: Is the value range reasonable? \item Increase the simulation time to e.g. $4\times 10^7 \;\text{s}$. Investigate the saturation: Is the value range reasonable?
\item What happens if you increase the resolution of the grid? \item What happens if you increase the resolution of the grid?
\end{itemize} \end{itemize}
\subsubsection{Exercise 3: Parameter file input.} \subsubsection{Exercise 3: Parameter File Input.}
As you have experienced, compilation takes quite some time. Therefore, \Dumux 2.1 provides a simple method to read in parameters (such as simulation end time or modelling parameters) via \texttt{Paramter Input Files}. The tests in the Test-folder \texttt{/test/} already use this system.\\
If you look at the Application in \texttt{/test/boxmodels/2p/}, you see that the main file looks rather empty: The parameter file \texttt{test\_2p.input} is read by a standard start procedure, which is called in the main function. This should be adapted for your problem at hand. The program run has to be called with the parameter file as argument.
In the code, parameters can be read via the macro \texttt{GET\_RUNTIME\_PARAM(TypeTag, Scalar, MyWonderfulGroup.MyWonderfulParameter);}. In \texttt{test\_2p}, \texttt{MyWonderfulGroup} is the group \texttt{SpatialParameters} - any type of groups is applicable, if the group definition in the parameter file is enclosed in square brackets. The parameters are then listed thereafter. Try and use as much parameters as possible via the input file, such as lens dimension, grid resolution, soil properties etc.
As you have experienced, compilation takes quite some time. Therefore,
\Dumux provides a simple method to read in parameters at run-time
via \textit{paramter input files}.\\
\subsubsection{Exercise 4: Create a new Fluid System.} In the code, parameters can be read via the macro
\texttt{GET\_RUNTIME\_PARAM(TypeTag, Scalar,
MyWonderfulGroup.MyWonderfulParameter);}. At the end of the
simulation a list of parameters is printed if the command line option
\texttt{-printParams 1} is passed to the simulation. Add some (for
example \texttt{Newton.MaxSteps} and \texttt{EnableGravity}) to the
parameter file \texttt{tutorial\_coupled.input} and observe what
happens if they are modified.
Create a new file for benzene called \texttt{benzene.hh} and implement \subsubsection{Exercise 4: Create a New Component.}
a new fluid system. (You may get a hint by looking at existing fluid
systems in the directory \verb+/dumux/material/fluidsystems+). \\ Create a new file for the benzene component called \texttt{benzene.hh}
and implement a new component. (You may get a hint by looking at
existing fluid systems in the directory \verb+/dumux/material/components+). \\
Use benzene as a new fluid and run the model of Exercise 2 with water Use benzene as a new fluid and run the model of Exercise 2 with water
and benzene. Benzene has a density of $889.51 \, \text{kg} / \text{m}^3$ and benzene. Benzene has a density of $889.51 \, \text{kg} /
and a viscosity of $0.00112 \, \text{Pa} \; \text{s}$. \text{m}^3$ and a viscosity of $0.00112 \, \text{Pa} \; \text{s}$.
\clearpage \newpage \clearpage \newpage
%%% Local Variables: %%% Local Variables:

83
examples/tutorial_coupled.cc
View File

@@ -3,7 +3,7 @@
/***************************************************************************** /*****************************************************************************
* Copyright (C) 2007-2008 by Klaus Mosthaf * * Copyright (C) 2007-2008 by Klaus Mosthaf *
* Copyright (C) 2007-2008 by Bernd Flemisch * * Copyright (C) 2007-2008 by Bernd Flemisch *
* Copyright (C) 2008-2009 by Andreas Lauser * * Copyright (C) 2008-2012 by Andreas Lauser *
* Institute for Modelling Hydraulic and Environmental Systems * * Institute for Modelling Hydraulic and Environmental Systems *
* University of Stuttgart, Germany * * University of Stuttgart, Germany *
* email: <givenname>.<name>@iws.uni-stuttgart.de * * email: <givenname>.<name>@iws.uni-stuttgart.de *
@@ -28,72 +28,29 @@
*/ */
#include "config.h" /*@\label{tutorial-coupled:include-begin}@*/ #include "config.h" /*@\label{tutorial-coupled:include-begin}@*/
#include "tutorialproblem_coupled.hh" /*@\label{tutorial-coupled:include-problem-header}@*/ #include "tutorialproblem_coupled.hh" /*@\label{tutorial-coupled:include-problem-header}@*/
#include <dumux/common/start.hh> /*@\label{tutorial-coupled:include-end}@*/
#include <dune/common/mpihelper.hh> //! Prints a usage/help message if something goes wrong or the user asks for help
#include <iostream> /*@\label{tutorial-coupled:include-end}@*/ void usage(const char *progName, const std::string &errorMsg) /*@\label{tutorial-coupled:usage-function}@*/
void usage(const char *progname)
{ {
std::cout << "usage: " << progname << " [--restart restartTime] tEnd dt\n"; std::cout
exit(1); << "\nUsage: " << progName << " [options]\n";
if (errorMsg.size() > 0)
std::cout << errorMsg << "\n";
std::cout
<< "\n"
<< "The List of Mandatory arguments for this program is:\n"
<< "\t-tEnd The end of the simulation [s]\n"
<< "\t-dtInitial The initial timestep size [s]\n"
<< "\t-Grid.upperRightX The x-coordinate of the grid's upper-right corner [m]\n"
<< "\t-Grid.upperRightY The y-coordinate of the grid's upper-right corner [m]\n"
<< "\t-Grid.numberOfCellsX The grid's x-resolution\n"
<< "\t-Grid.numberOfCellsY The grid's y-resolution\n"
<< "\n";
} }
int main(int argc, char** argv) int main(int argc, char** argv)
{ {
try { typedef TTAG(TutorialProblemCoupled) TypeTag; /*@\label{tutorial-coupled:set-type-tag}@*/
typedef TTAG(TutorialProblemCoupled) TypeTag; /*@\label{tutorial-coupled:set-type-tag}@*/ return Dumux::start<TypeTag>(argc, argv, usage); /*@\label{tutorial-coupled:call-start}@*/
typedef GET_PROP_TYPE(TypeTag, Grid) Grid; /*@\label{tutorial-coupled:retrieve-types-begin}@*/
typedef GET_PROP_TYPE(TypeTag, TimeManager) TimeManager;
typedef GET_PROP_TYPE(TypeTag, Problem) Problem; /*@\label{tutorial-coupled:retrieve-types-end}@*/
// Initialize MPI
Dune::MPIHelper::instance(argc, argv); /*@\label{tutorial-coupled:init-mpi}@*/
// parse the command line arguments
if (argc < 3) /*@\label{tutorial-coupled:parse-args-begin}@*/
usage(argv[0]);
// parse restart time if restart is requested
int argPos = 1;
bool restart = false;
double startTime = 0;
if (std::string("--restart") == argv[argPos]) {
restart = true;
++argPos;
// use restart time as start time
std::istringstream(argv[argPos++]) >> startTime;
}
// read the initial time step and the end time
if (argc - argPos != 2)
usage(argv[0]);
double tEnd, dt; /*@\label{tutorial-coupled:parse-tEn-and-dt}@*/
std::istringstream(argv[argPos++]) >> tEnd;
std::istringstream(argv[argPos++]) >> dt; /*@\label{tutorial-coupled:parse-args-end}@*/
// create the grid
Grid *gridPtr = GET_PROP(TypeTag, Grid)::create(); /*@\label{tutorial-coupled:create-grid}@*/
// create time manager responsible for global simulation control
TimeManager timeManager;
// instantiate the problem on the leaf grid
Problem problem(timeManager, gridPtr->leafView()); /*@\label{tutorial-coupled:instantiate-problem}@*/
timeManager.init(problem, startTime, dt, tEnd, restart); /*@\label{tutorial-coupled:initTimeManager}@*/
// run the simulation
timeManager.run(); /*@\label{tutorial-coupled:execute}@*/
return 0;
}
catch (Dune::Exception &e) { /*@\label{tutorial-coupled:catch-dune-exceptions}@*/
// Catch exceptions thrown somewhere in DUNE
std::cerr << "Dune reported error: " << e << std::endl;
}
catch (...) { /*@\label{tutorial-coupled:catch-other-exceptions}@*/
// Catch exceptions thrown elsewhere
std::cerr << "Unknown exception thrown!\n";
throw;
}
return 3;
} }

7
examples/tutorial_coupled.input Normal file
View File

@@ -0,0 +1,7 @@
tEnd=500000 # duration of the simulation [s]
dtInitial=10 # initial time step size [s]
[Grid]
upperRightX=300 # x-coordinate of the upper-right corner of the grid [m]
upperRightY=60 # y-coordinate of the upper-right corner of the grid [m]
numberOfCellsX=100 # x-resolution of the grid
numberOfCellsY=1 # y-resolution of the grid

107
examples/tutorialproblem_coupled.hh
View File

@@ -2,7 +2,7 @@
// vi: set et ts=4 sw=4 sts=4: // vi: set et ts=4 sw=4 sts=4:
/***************************************************************************** /*****************************************************************************
* Copyright (C) 2008-2009 by Melanie Darcis, Klaus Mosthaf * * Copyright (C) 2008-2009 by Melanie Darcis, Klaus Mosthaf *
* Copyright (C) 2009 by Andreas Lauser * * Copyright (C) 2009-2012 by Andreas Lauser *
* Institute for Modelling Hydraulic and Environmental Systems * * Institute for Modelling Hydraulic and Environmental Systems *
* University of Stuttgart, Germany * * University of Stuttgart, Germany *
* email: <givenname>.<name>@iws.uni-stuttgart.de * * email: <givenname>.<name>@iws.uni-stuttgart.de *
@@ -28,54 +28,36 @@
#ifndef DUMUX_TUTORIAL_PROBLEM_COUPLED_HH // guardian macro /*@\label{tutorial-coupled:guardian1}@*/ #ifndef DUMUX_TUTORIAL_PROBLEM_COUPLED_HH // guardian macro /*@\label{tutorial-coupled:guardian1}@*/
#define DUMUX_TUTORIAL_PROBLEM_COUPLED_HH // guardian macro /*@\label{tutorial-coupled:guardian2}@*/ #define DUMUX_TUTORIAL_PROBLEM_COUPLED_HH // guardian macro /*@\label{tutorial-coupled:guardian2}@*/
// the numerical model // The numerical model
#include <dumux/boxmodels/2p/2pmodel.hh> #include <dumux/boxmodels/2p/2pmodel.hh>
// the DUNE grid used // The DUNE grid used
#include <dune/grid/sgrid.hh> #include <dune/grid/yaspgrid.hh>
// spatialy dependent parameters // Spatialy dependent parameters
#include "tutorialspatialparameters_coupled.hh" #include "tutorialspatialparameters_coupled.hh"
// the components that are used // The components that are used
#include <dumux/material/components/h2o.hh> #include <dumux/material/components/h2o.hh>
#include <dumux/material/components/oil.hh> #include <dumux/material/components/oil.hh>
#include <dumux/common/structuredcubegridcreator.hh>
namespace Dumux namespace Dumux{
{ // Forward declaration of the problem class
// forward declaration of the problem class
template <class TypeTag> template <class TypeTag>
class TutorialProblemCoupled; class TutorialProblemCoupled;
namespace Properties namespace Properties {
{ // Create a new type tag for the problem
// create a new type tag for the problem
NEW_TYPE_TAG(TutorialProblemCoupled, INHERITS_FROM(BoxTwoP, TutorialSpatialParametersCoupled)); /*@\label{tutorial-coupled:create-type-tag}@*/ NEW_TYPE_TAG(TutorialProblemCoupled, INHERITS_FROM(BoxTwoP, TutorialSpatialParametersCoupled)); /*@\label{tutorial-coupled:create-type-tag}@*/
// Set the "Problem" property // Set the "Problem" property
SET_PROP(TutorialProblemCoupled, Problem) /*@\label{tutorial-coupled:set-problem}@*/ SET_PROP(TutorialProblemCoupled, Problem) /*@\label{tutorial-coupled:set-problem}@*/
{ typedef Dumux::TutorialProblemCoupled<TypeTag> type;}; { typedef Dumux::TutorialProblemCoupled<TypeTag> type;};
// Set the grid // Set grid and the grid creator to be used
SET_PROP(TutorialProblemCoupled, Grid) /*@\label{tutorial-coupled:set-grid}@*/ SET_TYPE_PROP(TutorialProblemCoupled, Grid, Dune::YaspGrid</*dim=*/2>); /*@\label{tutorial-coupled:set-grid}@*/
{ SET_TYPE_PROP(TutorialProblemCoupled, GridCreator, Dumux::CubeGridCreator<TypeTag>); /*@\label{tutorial-coupled:set-grid}@*/
typedef Dune::SGrid<2,2> type;
static type *create() /*@\label{tutorial-coupled:create-grid-method}@*/
{
typedef typename type::ctype ctype;
Dune::FieldVector<int, 2> cellRes; // vector holding resolution of the grid
Dune::FieldVector<ctype, 2> lowerLeft(0.0); // Coordinate of lower left corner of the grid
Dune::FieldVector<ctype, 2> upperRight; // Coordinate of upper right corner of the grid
cellRes[0] = 100;
cellRes[1] = 1;
upperRight[0] = 300;
upperRight[1] = 60;
return new Dune::SGrid<2,2>(cellRes,
lowerLeft,
upperRight);
}
};
// Set the wetting phase // Set the wetting phase
SET_PROP(TutorialProblemCoupled, WettingPhase) /*@\label{tutorial-coupled:2p-system-start}@*/ SET_PROP(TutorialProblemCoupled, WettingPhase) /*@\label{tutorial-coupled:2p-system-start}@*/
@@ -97,12 +79,10 @@ SET_BOOL_PROP(TutorialProblemCoupled, EnableGravity, false); /*@\label{tutorial-
} }
/*! /*!
* \ingroup TwoPBoxModel * \ingroup TwoPBoxModel
* *
* \brief Tutorial problem for a fully coupled twophase box model. * \brief Tutorial problem for a fully coupled twophase box model.
*/ */
// Definition of the actual problem
template <class TypeTag> template <class TypeTag>
class TutorialProblemCoupled : public TwoPProblem<TypeTag> /*@\label{tutorial-coupled:def-problem}@*/ class TutorialProblemCoupled : public TwoPProblem<TypeTag> /*@\label{tutorial-coupled:def-problem}@*/
{ {
@@ -134,23 +114,17 @@ public:
{ {
} }
// Specifies the problem name. This is used as a prefix for files //! Specifies the problem name. This is used as a prefix for files
// generated by the simulation. //! generated by the simulation.
const char *name() const const char *name() const
{ return "tutorial_coupled"; } { return "tutorial_coupled"; }
//! Returns true if a restart file should be written. //! Returns true if a restart file should be written.
/* The default behaviour is to write no restart file.
*/
bool shouldWriteRestartFile() const /*@\label{tutorial-coupled:restart}@*/ bool shouldWriteRestartFile() const /*@\label{tutorial-coupled:restart}@*/
{ { return false; }
return false;
}
//! Returns true if the current solution should be written to disk (i.e. as a VTK file) //! Returns true if the current solution should be written to disk
/*! The default behaviour is to write out the solution for //! as a VTK file
* every time step. Else, the user has to change the divisor in this function.
*/
bool shouldWriteOutput() const /*@\label{tutorial-coupled:output}@*/ bool shouldWriteOutput() const /*@\label{tutorial-coupled:output}@*/
{ {
return return
@@ -158,13 +132,13 @@ public:
(this->timeManager().timeStepIndex() % 1 == 0); (this->timeManager().timeStepIndex() % 1 == 0);
} }
// Return the temperature within a finite volume. We use constant //! Return the temperature within a finite volume. We use constant
// 10 degrees Celsius. //! 10 degrees Celsius.
Scalar temperature() const Scalar temperature() const
{ return 283.15; }; { return 283.15; };
// Specifies which kind of boundary condition should be used for //! Specifies which kind of boundary condition should be used for
// which equation for a finite volume on the boundary. //! which equation for a finite volume on the boundary.
void boundaryTypes(BoundaryTypes &BCtypes, const Vertex &vertex) const void boundaryTypes(BoundaryTypes &BCtypes, const Vertex &vertex) const
{ {
const GlobalPosition &pos = vertex.geometry().center(); const GlobalPosition &pos = vertex.geometry().center();
@@ -175,18 +149,19 @@ public:
} }
// Evaluate the Dirichlet boundary conditions for a finite volume //! Evaluate the Dirichlet boundary conditions for a finite volume
// on the grid boundary. Here, the 'values' parameter stores //! on the grid boundary. Here, the 'values' parameter stores
// primary variables. //! primary variables.
void dirichlet(PrimaryVariables &values, const Vertex &vertex) const void dirichlet(PrimaryVariables &values, const Vertex &vertex) const
{ {
values[Indices::pwIdx] = 200.0e3; // 200 kPa = 2 bar values[Indices::pwIdx] = 200.0e3; // 200 kPa = 2 bar
values[Indices::SnIdx] = 0.0; // 0 % oil saturation on left boundary values[Indices::SnIdx] = 0.0; // 0 % oil saturation on left boundary
} }
// Evaluate the boundary conditions for a Neumann boundary //! Evaluate the boundary conditions for a Neumann boundary
// segment. Here, the 'values' parameter stores the mass flux in //! segment. Here, the 'values' parameter stores the mass flux in
// normal direction of each phase. Negative values mean influx. //! [kg/(m^2 * s)] in normal direction of each phase. Negative
//! values mean influx.
void neumann(PrimaryVariables &values, void neumann(PrimaryVariables &values,
const Element &element, const Element &element,
const FVElementGeometry &fvElemGeom, const FVElementGeometry &fvElemGeom,
@@ -209,8 +184,8 @@ public:
} }
} }
// Evaluate the initial value for a control volume. For this //! Evaluate the initial value for a control volume. For this
// method, the 'values' parameter stores primary variables. //! method, the 'values' parameter stores primary variables.
void initial(PrimaryVariables &values, void initial(PrimaryVariables &values,
const Element &element, const Element &element,
const FVElementGeometry &fvElemGeom, const FVElementGeometry &fvElemGeom,
@@ -220,10 +195,10 @@ public:
values[Indices::SnIdx] = 1.0; values[Indices::SnIdx] = 1.0;
} }
// Evaluate the source term for all phases within a given //! Evaluate the source term for all phases within a given
// sub-control-volume. In this case, the 'values' parameter stores //! sub-control-volume. In this case, the 'values' parameter
// the rate mass generated or annihilate per volume unit. Positive //! stores the rate mass generated or annihilate per volume unit
// values mean that mass is created. //! in [kg / (m^3 * s)]. Positive values mean that mass is created.
void source(PrimaryVariables &values, void source(PrimaryVariables &values,
const Element &element, const Element &element,
const FVElementGeometry &fvElemGeom, const FVElementGeometry &fvElemGeom,

23
examples/tutorialspatialparameters_coupled.hh
View File

@@ -35,9 +35,7 @@
#include <dumux/material/fluidmatrixinteractions/2p/regularizedbrookscorey.hh> /*@\label{tutorial-coupled:rawLawInclude}@*/ #include <dumux/material/fluidmatrixinteractions/2p/regularizedbrookscorey.hh> /*@\label{tutorial-coupled:rawLawInclude}@*/
#include <dumux/material/fluidmatrixinteractions/2p/efftoabslaw.hh> #include <dumux/material/fluidmatrixinteractions/2p/efftoabslaw.hh>
namespace Dumux namespace Dumux {
{
//forward declaration //forward declaration
template<class TypeTag> template<class TypeTag>
class TutorialSpatialParametersCoupled; class TutorialSpatialParametersCoupled;
@@ -93,9 +91,8 @@ public:
// determine appropriate parameters depening on selected materialLaw // determine appropriate parameters depening on selected materialLaw
typedef typename MaterialLaw::Params MaterialLawParams; /*@\label{tutorial-coupled:matLawObjectType}@*/ typedef typename MaterialLaw::Params MaterialLawParams; /*@\label{tutorial-coupled:matLawObjectType}@*/
/*! Intrinsic permeability tensor K \f$[m^2]\f$ depending
//! Intrinsic permeability tensor K \f$[m^2]\f$ depending * on the position in the domain
/*! on the position in the domain
* *
* \param element The finite volume element * \param element The finite volume element
* \param fvElemGeom The finite-volume geometry in the box scheme * \param fvElemGeom The finite-volume geometry in the box scheme
@@ -107,12 +104,10 @@ public:
const Dune::FieldMatrix<Scalar, dim, dim> &intrinsicPermeability(const Element &element, /*@\label{tutorial-coupled:permeability}@*/ const Dune::FieldMatrix<Scalar, dim, dim> &intrinsicPermeability(const Element &element, /*@\label{tutorial-coupled:permeability}@*/
const FVElementGeometry &fvElemGeom, const FVElementGeometry &fvElemGeom,
int scvIdx) const int scvIdx) const
{ { return K_; }
return K_;
}
//! Define the porosity \f$[-]\f$ of the porous medium depending /*! Define the porosity \f$[-]\f$ of the porous medium depending
/*! on the position in the domain * on the position in the domain
* *
* \param element The finite volume element * \param element The finite volume element
* \param fvElemGeom The finite-volume geometry in the box scheme * \param fvElemGeom The finite-volume geometry in the box scheme
@@ -121,12 +116,10 @@ public:
* Alternatively, the function porosityAtPos(const GlobalPosition& globalPos) could be defined, where globalPos * Alternatively, the function porosityAtPos(const GlobalPosition& globalPos) could be defined, where globalPos
* is the vector including the global coordinates of the finite volume. * is the vector including the global coordinates of the finite volume.
*/ */
double porosity(const Element &element, /*@\label{tutorial-coupled:porosity}@*/ Scalar porosity(const Element &element, /*@\label{tutorial-coupled:porosity}@*/
const FVElementGeometry &fvElemGeom, const FVElementGeometry &fvElemGeom,
int scvIdx) const int scvIdx) const
{ { return 0.2; }
return 0.2;
}
/*! Return the parameter object for the material law (i.e. Brooks-Corey) /*! Return the parameter object for the material law (i.e. Brooks-Corey)
* depending on the position in the domain * depending on the position in the domain