Personal tools

Difference between revisions of "Module PorousMedia"

From MohidWiki

Jump to: navigation, search
(Water content)
(Water content)
Line 73: Line 73:
 
   
 
   
 
<math>\left\{\begin{matrix}
 
<math>\left\{\begin{matrix}
h=-DWZ\cdot 0.5\,\,\,\,\,\,for\,\, the\,\, cells\,\, immediately\,\, above\,\, the\,\, water\,\, table  
+
\,\,\,\,h=-DWZ\cdot 0.5\,\,\,\,\,\,for\,\, the\,\, cells\,\, immediately\,\, above\,\, the\,\, water\,\, table  
 
\\  
 
\\  
 
\\\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\! \!\!\!\!\!\!h=-(-DZZ-h)\,\,\,\,\,\,for\,\, the\,\, other\,\, cells
 
\\\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\! \!\!\!\!\!\!h=-(-DZZ-h)\,\,\,\,\,\,for\,\, the\,\, other\,\, cells
Line 81: Line 81:
  
 
[[Image:Figure05.jpg|thumb|center|590px|Figure 1: Suction Head evaluation]]
 
[[Image:Figure05.jpg|thumb|center|590px|Figure 1: Suction Head evaluation]]
 +
 +
 +
As shown in the picture the suction head is calculated in order to maintain the same total head (H= in the cells in agreement with the field capacity definition.
  
 
===Water retention===
 
===Water retention===

Revision as of 13:09, 3 March 2011

Overview

Module Porous Media is based on the Finite Volume Method in order to take in account the non-linear behavior of the soil. This module allows the calculation of the water fluxes in a porous media. In particular the infiltration flow is calculated from the potential water column available for this process. The potential water column is calculated by the Module Basin as following:


PIC=ThroughFall-IWC


where:

PIC is the Potential Infiltration column (m)
ThroughFall is the percentage of the rain that reaches the soil (m)
IWC is the Initial Water column (m)


The infiltration flows is calculated by the Buckingham-Darcy equation (Jury et al,1991)

Main Processes

Water flow

In general the flow evaluation is based on the conservation equations of mass and momentum. In the case of soil it is assumed that the forces of inertia are almost zero; therefore the balance is reduced to the forces of pressure, gravity and viscous. The equation that describes the motion of the water in the soil is the Buckingham Darcy equation (Jury et al,1991).


Failed to parse (unknown error): \left\{\begin{matrix} J=-K\left ( \theta \right )\left ( \frac{\partial h}{\partial \theta } + 1 \right )\,\,\,\,[1]\\ \\ \\h=T+P+\psi\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,[2]\\ \end{matrix}\right .


where:

J is the water velocity (m/s)
H is the pressure head (m)
θ is the water content
K is the hydraulic conductivity
T is the topography (m)
P is the hydrostatic pressure (m)
Failed to parse (unknown error): \psi is the suction head (m)


The equation of Buckingham Darcy is applied at the cells interface in order to obtain the water velocity between two cell and subsequently the flow by multiplying the obtained velocity for the cross section.

Water content

The Water content is the amount of water stored in the soil, usually it is indicated with the lecter \theta. After the determination of aquifer level (water table) the water content is associated at each cells by the following criteria:

  • If the cell is located above the water table \theta=\theta_{s}
  • If the cell is located over the water table \theta=\theta_{ns}

where:

\theta=\theta_{s} is the water content in the saturated soil(m3/m3)
\theta=\theta_{ns} is the water content in the non saturated soil (m3/m3)

The water content in the non saturated soil can be associated with observed data, if these are not available it is imposed by default equal to the Field Capacity. In order to calculate this parameter by the equation [] the evaluation of the suction head is needed as following:


Failed to parse (unknown error): \left\{\begin{matrix} \,\,\,\,h=-DWZ\cdot 0.5\,\,\,\,\,\,for\,\, the\,\, cells\,\, immediately\,\, above\,\, the\,\, water\,\, table \\ \\\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\! \!\!\!\!\!\!h=-(-DZZ-h)\,\,\,\,\,\,for\,\, the\,\, other\,\, cells \end{matrix}\right


Figure 1: Suction Head evaluation


As shown in the picture the suction head is calculated in order to maintain the same total head (H= in the cells in agreement with the field capacity definition.

Water retention

The shape of water retention curves can be characterized by several models, one of them known as the van Genuchten model:

\theta(h) = \theta_r + \frac{\theta_s - \theta_r}{\left[ 1+(\alpha |h|)^n \right]^{1-1/n}}

where

\theta(h) is the water retention curve [L3L−3];
|h| is suction pressure ([L−1] or cm of water);
\theta_s saturated water content [L3L−3];
\theta_r residual water content [L3L−3];
\alpha is related to the inverse of the air entry suction, \alpha >0 ([L−1], or cm−1); and,
n is a measure of the pore-size distribution, n>1 (dimensionless).

Evapotranspiration

Some water may disappear from the soil because of the evaporation and transpiration processes, which become a sink in soil water profile. These two processes, currently named Evapotranspiration may be modeled using the Penmann Monteith equation.

 \overset{\text{Energy flux rate}}{\lambda_v E=\frac{\Delta R_n   +   \rho_a c_p  \left(  \delta q  \right) g_a }{\Delta  + \gamma \left (    1 + g_a / g_s    \right)}}~ \iff ~  \overset{\text{Volume flux rate}}{ET_o=\frac{\Delta R_n   +   \rho_a c_p  \left(  \delta q  \right) g_a } { \left(   \Delta  + \gamma \left (    1 + g_a / g_s    \right)    \right) \lambda_v }}

λv = Latent heat of vaporization. Energy required per unit mass of water vaporized. (J/g)
Lv = Volumetric latent heat of vaporization. Energy required per water volume vaporized. (Lv = 2453 MJ m-3)
E = Mass water evapotranspiration rate (g s-1 m-2)
ETo = Water volume evapotranspired (m3 s-1 m-2)
Δ = Rate of change of saturation specific humidity with air temperature. (Pa K-1)
Rn = Net irradiance (W m-2), the external source of energy flux
cp = Specific heat capacity of air (J kg-1 K-1)
ρa = dry air density (kg m-3)
δe = vapor pressure deficit, or specific humidity (Pa)
ga = Hydraulic conductivity of air, atmospheric conductance (m s-1)
gs = Conductivity of stoma, surface conductance (m s-1)
γ = Psychrometric constant (γ ≈ 66 Pa K-1)

All of these calculations about potential evapotranspiration are made in module Basin in MohiLand model. However, not all of the potential water that can be evaporated or transpirated will be in fact removed from the soil. The water that will really leave the soil through these processes is calculated in the Porous Media module. In Figure below it can be seen that the actual transpiration and evaporation are calculated in Porous Media module of MohidLand. The actual evaporation, which happens only at the soil surface, is calculated based on a pressure head limit imposed by the user. It allows the model to not evaporate any surface water, even if it is available, when the soil is very far way from the saturation point.

Evapotranspiration fluxogram in Mohid Land model

Other Features

How to Generate Info needed in Porous Media

SoilMap

Model needs to know soil ID in each cell and layer to pick hydraulic properties from that type of soil.

Options:

  • Constant value
  • Soil Grid. One possible option is to associate with soil shape file. In this case can use MOHID GIS going to menu [Tools]->[Shape to Grid Data] and provide: i) the grid (model grid), ii) the soil shape file and iii) the corespondence between soil codes and soil ID defined in data file.

Soil ID must be defined in Module FillMatrix standards for each soil horizon defined (grid example):

<beginhorizon>
KLB                       : 1
KUB                       : 10

<beginproperty>
NAME                      : SoilID
DEFAULTVALUE              : 1
INITIALIZATION_METHOD     : ASCII_FILE
REMAIN_CONSTANT           : 1
FILENAME                  : ..\..\GeneralData\PorousMedia\SoilID200m.dat
<endproperty>

<beginproperty>
<endproperty>
..
<endhorizon>

Remarks

All the soil ID's appearing in the soil grid(s) must be defined in the PorousMedia data file in terms of hydraulic properties:

<beginsoiltype>
ID                        : 1
THETA_S                   : 0.3859
THETA_R                   : 0.0476
N_FIT                     : 1.39
SAT_K                     : 3.5556e-6
ALPHA                     : 2.75
L_FIT                     : 0.50
THETA_CV_MIN              : 0.2844
THETA_CV_MAX              : 0.3791
<endsoiltype>

<beginsoiltype>
ID                        : 2
...
<endsoiltype>
...

Soil Bottom

The soil depth must be known by the model. This is computed by the model from terrain altitude (topography) and soil bottom altitude. As so, a soil bottom grid is needed.

Options:

  • Grid File.

Soil depth (and soil bottom altitude, the effective grid needed) can be defined with a constant depth or estimated from slope HOW TO SoilBottom LINK.

Define the grid just generated, in the porous media data file with:

BOTTOM_FILE  : ..\..\GeneralData\PorousMedia\BottomLevel.dat

Water Level

Options:

  • Grid File.

The water table altitude represents the initial altitude of the water table. It is recommended to do a spin-up run to estabilize water level and then do a continuous simulation starting with the final water table achieved. Use the following blocks with Module FillMatrix standards:

<beginwaterlevel>
NAME                      : waterlevel
INITIALIZATION_METHOD     : ASCII_FILE
DEFAULTVALUE              : 0
REMAIN_CONSTANT           : 0
FILENAME                  : ..\..\GeneralData\PorousMedia\WaterLevel0.50.dat
<endwaterlevel>

Impermeability

Impermeability values (0 - completely permeable, 1 - impermeable) must be provided.

Options:

  • Constant Value.
  • Grid File. One possible option is to associate with land use shape file. In this case can use MOHID GIS going to menu [Tools]->[Shape to Grid Data] and provide: i) the grid (model grid), ii) the land use shape file and iii) the corespondence between land use codes and Impermeability values.

Use the following blocks with Module FillMatrix standards:

<beginimpermeablefraction>
NAME                      : impermeablefraction
INITIALIZATION_METHOD     : ASCII_FILE
DEFAULTVALUE              : 0
REMAIN_CONSTANT           : 1
FILENAME                  : ..\..\GeneralData\PorousMedia\AreaImpermeavel.dat
<endimpermeablefraction>

Outputs

References

  • Jury,W.A.,Gardner,W.R.,Gardner,W.H., 1991,Soil Physics
  • Van Genuchten, M.T., A closed form equation for predicting the hydraulic conductivity of unsaturated soils
  • Wu,J.,Zhang, R., Gui,S.,1999, Modelling soil water movement with water uptake by roots, Plant and soil 215: 7-17

Data File

Keywords

Keywords read in the Data File

Keyword                   : Data Type         Default     !Comment

BOTTOM_FILE               : char              -           !Path to Bottom Topography File
START_WITH_FIELD          : logical           1           !Sets Theta initial Field Capacity
CONTINUOUS                : logical           0           !Continues from previous run
STOP_ON_WRONG_DATE        : logical           1           !Stops if previous run end is different from actual
                                                          !Start
OUTPUT_TIME               : sec. sec. sec.    -           !Output Time
TIME_SERIE_LOCATION       : char              -           !Path to File which defines Time Series
CONTINUOUS_OUTPUT_FILE    : logical           1           !Writes "famous" iter.log
CONDUTIVITYFACE           : integer           1           !Way to interpolate conducivity face
                                                          !1 - Average, 2 - Maximum, 3 - Minimum, 4 - Weigthed
HORIZONTAL_K_FACTOR       : real              1.0         !Factor for Horizontal Conductivity = Kh / Kv
CUT_OFF_THETA_LOW         : real              1e-6        !Disables calculation when Theta is near ThetaR
CUT_OFF_THETA_HIGH        : real              1e-15       !Set Theta = ThetaS when Theta > ThetaS - CUT_OFF_THETA_HIGH
MIN_ITER                  : integer           2           !Number of iterations below which the DT is increased
MAX_ITER                  : integer           3           !Number of iterations above which the DT is decreased
LIMIT_ITER                : integer           50          !Number of iterations of a time step (for restart)
THETA_TOLERANCE           : real              0.001       !Converge Parameter
INCREASE_DT               : real              1.25        !Increase of DT when iter < MIN_ITER
DECREASE_DT               : real              0.70        !Decrease of DT when iter > MAX_ITER

<beginproperty>
NAME                      : Theta / waterlevel 

see Module FillMatrix for more options

<endproperty>

Sample

Some keywords of the PorousMedia input file:

<beginsoiltype>
*Saturated water content
THETA_S                   : 0.55

*residual water content
THETA_R                   : 0.1

*N of Van Genuchten: 1.5(clay) to 4.5 (sandy)
N_FIT                     : 1.15

*Saturated hydraulic conductivity
SAT_K                     : 1.39e-8

*alpha of Van Genuchten: 0.005 (clay) - 0.035 (sandy)
ALPHA                     : 0.02

*L of Mualem - Van Genuchten (mostly 0.5)
L_FIT                     : 0.50
<endsoiltype>