In general OpenGeoSys is not using any intrinsic units, i.e. OGS assumes that a selfconsistent set of units (SI for example) is used for all quantities. However, there are some exceptions to the rule as for instance some empirical laws like some expressions for the vapor diffusion and latent heat that assume the temperature to be given in Kelvin. Therefore, we recommend using SI base units.
For more information on how to come up with a selfconsistent unit scheme see this PDF.
Most processes in OpenGeoSys are derived and implemented with a 3D setting in mind. I.e., the conceptual model of the processes exists in a 3D world and boundaries are twodimensional. The canonical (a) Neumann boundary conditions and (b) volumetric source terms for those processes are to (a) assign a flux (aka current density, i.e. quantity $Q$ per time per 2D area) to the boundary of the domain or to (b) add a quantity rate density ($Q$ per time per 3D volume) to the simulation domain, respectively. Here, $Q$ is usually a different quantity than the primary unknown of the process. E.g., for energy balances in OpenGeoSys the primary unknown usually is the temperature $T$, but $Q$ is energy.
Note: The connection between the primary unknown and these canonical BCs and STs is made via the implemented PDE. And if a PDE is implemented slightly differently (e.g. PDE (B) is PDE (A) multiplied by a factor of density), the corresponding quantities in the BCs and STs generally also change (by a factor of density in this example).
Lowerdimensional realizations of the mentioned processes essentially are slices through a 3D world. Hence, the parameterization does not change. Despite solving a problem whose solution varies only in one or two dimensions, you still have to apply Neumann BCs and volumetric STs as quantity per time per 2D area and quantity per time per 3D volume, respectively. Sometimes, these lower dimensional realizations are referred to as crosssectional models.
In contrast, for real twodimensional and onedimensional problems the currents of Neumann boundary conditions and volumetric source terms can be seen as normalized by their one and twodimensional boundaries, i.e. quantity per time per length^{dim} in the case of volumetric STs.
Attention: The units given below apply to 3D and said crosssectional models! That applies to all of OpenGeoSys’s processes; there are currently no processes in OpenGeoSys that are manifest 1D or 2D.
Units:
Quantity name  Dimension symbol  SI base units 

time  $\mathrm{T}$  $\mathrm{s}$ 
length  $\mathrm{L}$  $\mathrm{m}$ 
mass  $\mathrm{M}$  $\mathrm{kg}$ 
temperature  $\mathrm{\Theta}$  $\mathrm{K}$ 
Temperature:
Displacement:
Pressure:
Pressure (Liquid Flow Process):
The liquid flow process uses a volumebased PDE, not a massbased one. Therefore, the list above does not apply to the liquid flow process, but you have to divide all units of natural boundary conditions and source terms by density, leading to the following list:
OpenGeoSys uses the positive sign convention for tensile stresses, whereas compressional stresses carry a negative sign.
OpenGeoSys assumes plane strain conditions as default for 2D simulations.
OpenGeoSys uses internally a vector of primary vector components (“global component vector”). Although, which set of components is used, varies from process to process, the order of primary variable components for most THMbased processes is: $T$, $p$, $u_{x}$, $u_{y}$, $u_{z}$. However, there are some exceptions and as there are also different primary variables for each process, you will find a list containing the process variables as they appear in the global component vector. This order is used, e.g., to display the per component convergence of the nonlinear solver in the output.
To map the elasticity/stiffness tensor OpenGeoSys internally uses a Kelvin mapping with an adapted component ordering for computational reasons [1]. For 2D, the KelvinVector of the stress tensor looks like $\sigma=(\sigma_{xx},\sigma_{yy},\sigma_{zz},\sqrt{2}\sigma_{xy})$ whereas the 3D version reads as $\sigma=(\sigma_{xx},\sigma_{yy},\sigma_{zz},\sqrt{2}\sigma_{xy}, \sqrt{2}\sigma_{yz},\sqrt{2}\sigma_{xz})$.
For Kelvin mapping also see the conversion function documentation.
The input and output of symmetric tensors consists of the full (symmetric) tensor elements (without the factor $\sqrt{2}$), retaining the same order. I.e., Input and output components are be $\sigma=(\sigma_{xx},\sigma_{yy},\sigma_{zz},\sigma_{xy})$ and $\sigma=(\sigma_{xx},\sigma_{yy},\sigma_{zz},\sigma_{xy},\sigma_{yz},\sigma_{xz})$, respectively.
A staggered scheme solves coupled problems by alternating on separate physical domains (e.g. thermal and mechanical) in contrast to monolithic schemes which solve all domains simultaneously (e.g. thermomechanical). Thus staggered schemes add another level of iterations, however this may pay off since the subproblems are smaller and they enable a finer tuning of the specific solvers.
For hydromechanical processes the fixedstress split has been implemented, since it turned out advantageous [1].
For sake of brevity, we do not describe the scheme itself, but intend to provide guidance for its stabilization parameter.
On this parameter depends how many coupling iterations are needed and thus how long it takes to obtain a solution.
The optimal value of this parameter is not apriori known, only that it lies within a certain interval [2]
\begin{equation}
\frac{1}{2}\frac{\alpha^2}{K_\mathrm{1D}} \le \beta_\mathrm{FS} \le \frac{\alpha^2}{K_\mathrm{ph}},
\end{equation}
where $\alpha$ denotes Biot’s coefficient, and $K_\mathrm{1D}$ and $K_\mathrm{ph}$ are the bulk moduli described next.
By $K_\mathrm{ph}$ we mean the bulk modulus adjusted to the spatial dimension, so in three dimensions it coincides with the drained bulk modulus, whereas in lower dimensions it becomes a constrained bulk modulus (2D plane strain, 1D uniaxial strain).
Assuming isotropic, linear elasticity we have
\begin{eqnarray}
K_\mathrm{3D} &=& \lambda + \frac{2}{3}\mu, \
K_\mathrm{2D} &=& \lambda + \frac{2}{2}\mu, \
K_\mathrm{1D} &=& \lambda + \frac{2}{1}\mu. \
\end{eqnarray}
OGS sets the stabilization parameter, which corresponds to a coupling compressibility, via the coupling_scheme_parameter
\begin{equation}
\beta_\mathrm{FS} = p_\mathrm{FS} \frac{\alpha^2}{K_\mathrm{3D}},
\end{equation}
by default to $p_\mathrm{FS}=\frac{1}{2}$.
For isotropic, linear elasticity we provide the interval [2] and the recommended value [3] in dependence on Poisson’s ratio $\nu$ (note $\frac{\lambda}{\mu}=\frac{2\nu}{12\nu}$).
2D  3D  

$p_\mathrm{FS}^\mathrm{min}$  $\frac{1}{6}\frac{1+\nu}{1\nu}$  same as 2D 
$p_\mathrm{FS}^\mathrm{MW}$  $\frac{1+\nu}{3}$  $\frac{1}{2}$ 
$p_\mathrm{FS}^\mathrm{max}$  $\frac{2(1+\nu)}{3}$  $1$ 
[1] Kim, J. and Tchelepi, H.A. and Juanes, R. (2009): Stability, Accuracy and Efficiency of Sequential Methods for Coupled Flow and Geomechanics. SPE International, vol. 16, p. 249262, DOI:https://doi.org/10.2118/119084PA https://onepetro.org/SJ/articleabstract/16/02/249/204235
[2] Storvik, E. and Both, J.W. and Kumar, K. and Nordbotten, J.M. and Radu, F.A. (2019): On the optimization of the fixedstress splitting for Biots equations. International Journal for Numerical Methods in Engineering, vol. 120, p. 179194, DOI:10.1002/nme.6130 https://onlinelibrary.wiley.com/doi/full/10.1002/nme.6130
[3] Mikelic, A. and Wheeler, M.F. (2013): Convergence of iterative coupling for coupled flow and geomechanics. Computational Geosciences, vol. 17, p. 455461, DOI:https://doi.org/10.1007/s105960129318y https://link.springer.com/content/pdf/10.1007/s105960129318y.pdf
This article was written by Joerg Buchwald. If you are missing something or you find an error please let us know.
Generated with Hugo 0.101.0
in CI job 319943

Last revision: February 16, 2023
Commit: [PL] Extract another loop and collection of IVs de82977c1
 Edit this page on