- Generated on Mon Feb 3 2020 19:58:04 for CLASSIC by
1.8.13
|
CLASSIC
Canadian Land Surface Scheme including Biogeochemical Cycles
|
Evaluates infiltration of water into soil under saturated conditions. More...
Functions/Subroutines | |
| subroutine | waterinfiltratesat (WMOVE, TMOVE, LZF, NINF, TRMDR, TPOND, ZPOND, R, TR, EVAP, PSIF, GRKINF, THLINF, THLIQX, TBARWX, DELZX, ZBOTX, FMAX, ZF, DZF, DTFLOW, THLNLZ, THLQLZ, DZDISP, WDISP, WABS, ITER, NEND, ISIMP, IGRN, IG, IGP1, IGP2, ILG, IL1, IL2, JL, N) |
Evaluates infiltration of water into soil under saturated conditions.
| subroutine waterinfiltratesat | ( | real, dimension (ilg,igp2), intent(inout) | WMOVE, |
| real, dimension (ilg,igp2), intent(inout) | TMOVE, | ||
| integer, dimension (ilg), intent(inout) | LZF, | ||
| integer, dimension (ilg), intent(inout) | NINF, | ||
| real, dimension (ilg), intent(inout) | TRMDR, | ||
| real, dimension (ilg), intent(inout) | TPOND, | ||
| real, dimension (ilg), intent(inout) | ZPOND, | ||
| real, dimension (ilg), intent(in) | R, | ||
| real, dimension (ilg), intent(in) | TR, | ||
| real, dimension (ilg), intent(in) | EVAP, | ||
| real, dimension (ilg,igp1), intent(in) | PSIF, | ||
| real, dimension(ilg,igp1), intent(in) | GRKINF, | ||
| real, dimension(ilg,igp1), intent(in) | THLINF, | ||
| real, dimension(ilg,igp1), intent(in) | THLIQX, | ||
| real, dimension(ilg,igp1), intent(in) | TBARWX, | ||
| real, dimension (ilg,igp1), intent(in) | DELZX, | ||
| real, dimension (ilg,igp1), intent(in) | ZBOTX, | ||
| real, dimension (ilg), intent(inout) | FMAX, | ||
| real, dimension (ilg), intent(inout) | ZF, | ||
| real, dimension(ilg), intent(inout) | DZF, | ||
| real, dimension(ilg), intent(inout) | DTFLOW, | ||
| real, dimension(ilg), intent(inout) | THLNLZ, | ||
| real, dimension(ilg), intent(inout) | THLQLZ, | ||
| real, dimension(ilg), intent(inout) | DZDISP, | ||
| real, dimension(ilg), intent(inout) | WDISP, | ||
| real, dimension(ilg), intent(inout) | WABS, | ||
| integer, dimension(ilg), intent(inout) | ITER, | ||
| integer, dimension(ilg), intent(inout) | NEND, | ||
| integer, dimension(ilg), intent(inout) | ISIMP, | ||
| integer, dimension (ilg), intent(in) | IGRN, | ||
| integer, intent(in) | IG, | ||
| integer, intent(in) | IGP1, | ||
| integer, intent(in) | IGP2, | ||
| integer, intent(in) | ILG, | ||
| integer, intent(in) | IL1, | ||
| integer, intent(in) | IL2, | ||
| integer, intent(in) | JL, | ||
| integer, intent(in) | N | ||
| ) |
| [in,out] | wmove | Water movement matrix \([m^3 m^{-2}]\) |
| [in,out] | tmove | Temperature matrix associated with ground water movement [C] |
| [in,out] | lzf | Index of soil layer in which wetting front is located |
| [in,out] | ninf | Number of levels involved in water movement |
| [in,out] | trmdr | Time remaining in current time step [s] |
| [in,out] | tpond | Temperature of ponded water [C] |
| [in,out] | zpond | Depth of ponded water \([m] (z_p)\) |
| [in] | r | Rainfall rate at ground surface \([m s^{-1}]\) |
| [in] | tr | Temperature of rainfall [C] |
| [in] | evap | Surface evaporation rate \([m s^{-1}]\) |
| [in] | psif | Soil water suction across the wetting front \([m] (\Psi_f)\) |
| [in] | grkinf | Hydraulic conductivity of soil behind the wetting front \([m s^{-1}] (K)\) |
| [in] | thlinf | Volumetric liquid water content behind the wetting front \([m^3 m^{-3}]\) |
| [in] | thliqx | Volumetric liquid water content of soil layer \([m^3 m^{-3}]\) |
| [in] | tbarwx | Temperature of water in soil layer [C] |
| [in] | delzx | Permeable depth of soil layer \([m] (\Delta z_{z,w})\) |
| [in] | zbotx | Depth of bottom of soil layer [m] |
| [in,out] | fmax | Maximum infiltration rate, defined as minimum value of GRKINF |
| [in,out] | zf | Depth of the wetting front \([m] (z_f)\) |
| [in] | igrn | Flag to indicate whether infiltration is occurring |
General calculations performed: The infiltration rate \(F_{inf}\) under saturated conditions is calculated as
\(F_{inf} = K_{inf} [(\Psi_f + z_f + z_p)/z_f ]\)
where \(K_{inf}\) is the hydraulic conductivity of the soil behind the wetting front, \(\Psi_f\) is the soil moisture suction across the wetting front, \(z_f\) is the depth of the wetting front, and \(z_p\) is the depth of water ponded on the surface. Since \(\Psi_f\) will vary by soil layer, and \(z_f\) and \(z_p\) will vary with time, the period of infiltration is divided into two-minute segments. A maximum iteration flag NEND is defined as the number of iteration segments plus 100, as a safeguard against runaway iterations. Since the surface infiltration rate will be limited by the lowest infiltration rate encountered, a maximum infiltration rate FMAX is defined as the minimum value of GRKINF, the hydraulic conductivity behind the wetting front, in all of the soil layers above and including that containing the wetting front.
In loop 200, a check is performed to determine whether the liquid water content of the current soil layer, THLIQX, equals or exceeds that behind the wetting front, THLINF. If so, the depth of the wetting front is reset to the depth of the soil layer. The water content of the soil layer is assigned to the level of the water movement matrix WMOVE corresponding to the soil layer index plus 1, and the temperature of the water in the layer is assigned to the same level of the matrix TMOVE. FMAX is recalculated, and the flags LZF and NINF are each incremented by 1. The flag ISIMP is set to 1, indicating that the following calculations are to be bypassed for this iteration.
If ISIMP is not 1, the time period DTFLOW applying to the current iteration loop is set to two minutes or to the remainder of the time step, whichever is less. If GRKINF for the current soil layer is not vanishingly small, the flag ISIMP is set to -2. Otherwise, it is deemed that infiltration is suppressed, and the temperature and depth of water ponded on the surface are simply updated. The temperature of the ponded water is calculated as the weighted average of the current pond temperature and the rainfall added to it. The ponded water remaining after rainfall and evaporation have taken place is calculated as ZPTEST. If ZPTEST is less than zero, it is deduced that evaporation must have removed the ponded water before the end of the current iteration loop. If this is the case, the time period of the current iteration loop is recalculated as the amount of time required for the difference between the evaporation and the rainfall rates to consume the ponded water, and the pond depth ZPOND is set to zero. Otherwise, ZPOND is set to ZPTEST. Finally, ISIMP is set to -1.
The 400 loop addresses the infiltration process under saturated conditions. Such infiltration is modelled as “piston flow”. First the current infiltration rate FINF is calculated using the equation given above. If the wetting front has passed the bottom of the lowest soil layer, PSIF Cis neglected. If FINF is greater than the rainfall rate R, FINF is set equal to R; if FINF is greater than CFMAX, FINF is set equal to FMAX. If the wetting front has not passed the bottom of the lowest soil layer, the change in depth of the wetting front DZF over the current time interval is calculated from the amount of infiltrating water WINF and the volumetric water content behind the wetting front, THLINF. The amount of soil water WDISP displaced by this movement of the wetting front is calculated from DZF and THLIQX, and the depth to which this water penetrates, DZDISP, is calculated as WDISP/(THLINF – THLIQX). The amount of soil water WABS entrained by the movement of WDISP itself is calculated as the product of DZDISP and THLIQX. The change in depth of the wetting front, behind which the liquid water content of the soil layer is THLINF, is now the sum of the depth represented by infiltration of water at the surface, DZF, and the depth represented by displacement of the pre-existing soil water, DZDISP. If this change in depth causes the wetting front to move beyond the bottom of the current soil layer, DTFLOW is recalculated as the amount of time required for the composite wetting front to reach the bottom of the soil layer, and DZF, WDISP, DZDISP and WABS are likewise recalculated. As in the case for ISIMP = -1, the temperature of the ponded water on the surface is calculated as the weighted average of the current pond temperature and the rainfall added to it. The ponded water remaining after rainfall, infiltration and evaporation have taken place is calculated as ZPTEST. If ZPTEST is less than zero, it is deduced that infiltration and evaporation must have removed the ponded water before the end of the current iteration loop. If this is the case, the time period of the current iteration loop is recalculated as the amount of time required for the infiltration and evaporation minus the rainfall to consume the ponded water. The pond depth ZPOND is set to zero; if the wetting front has not passed the bottom of the lowest soil layer, DZF, WDISP, DZDISP and WABS are recalculated. Otherwise, ZPOND is set to ZPTEST. Finally, the first layer of the water movement matrix WMOVE is incremented with the amount of the infiltrated water, and the first layer of the matrix TMOVE with the temperature of the infiltrated water. The last layer of WMOVE is incremented by WDISP+WABS, and the depth of the wetting front by DZF+DZDISP.
In the 500 loop, if ISIMP < 0 (i.e. water movement has occurred), the time remaining in the current time step is recalculated. If the wetting front is at a soil layer boundary, if the time remaining is non-zero, if the wetting front is still within the modelled soil column, and if the calculated water content behind the wetting front is greater than zero, FMAX is updated, and LZF and NINF are incremented by 1. The NINF level of the TMOVE matrix is set to the water temperature of the new soil layer.
At the end of the iteration pass, checks are done in the 600 loop to ascertain whether the number of iteration passes is still less than NEND, whether either ponded water still exists or rain is still falling, and whether there is still time remaining in the current time step. If these conditions are all fulfilled, the counter NPNTS representing the number of points in the current vector of mosaic tiles for which infiltration is still occurring is incremented by 1. Otherwise, the iteration flag for the current tile is changed from 1 to 0, signaling the end of the saturated infiltration calculations for that tile.
1.8.13