Calculates infiltration into a layered soil
View definitions of variables used in this subroutine:
You may press here to view a diagram of the layer structure used by subroutine SOILW.
Every soil layer is associated with a node, which is defined as the top of that layer. Hence, the distance of the node JZ from the soil surface (i.e. its depth) is the distance of the top of the layer JZ from the soil surface. This distance has the variable name Z(JZ) [m]. The variable DZ(JZ) denotes the thickness of layer JZ. V(JZ) indicates a "volume" to be associated with the node JZ, and is typically taken to be the average thickness of the 2 layers surrounding the node.
The layers discussed above are mathematical constructs. A real soil profile, however, may have real structural layers within it, separated by discontinuities in soil hydraulic properties. This presents a challenge for Newton-Raphson transport solvers: the driving potentials, concentrations and conductivities are discontinuous about these boundaries.
To deal with this problem, every soil layer in Cupid is divided into an upper and lower half. The potential PN(IHR,JZ) is assigned to the the node between layers JZ-1 and JZ. Water contents, conductivities, and matric flux potentials in these halves are computed from the hydraulic properties of the layer and the matric potential of the nearest node. E.g., the water content in the lower half of layer JZ is computed as:
WNL(IHR,JZ)=WS(JZ)*(PE(JZ)/PN(IHR,JZ+1))**(1/BX(JZ))
A node must always be located where a discontinuity in soil properties occurs.
The soil water infiltration mode used in Cupid Version 7 is adapted from code supplied by G.S. Campbell. It uses a matric flux potential formalism, but is written explicitly in terms of the matric potential (see Campbell, Soil Physics with BASIC (1983) Chapter 8 for a more detailed discussion).
The goal is to obtain mass balance at each node in the soil:
F(JZ) = (LIQUID+VAPOR) MATRIC FLUX + GRAVITATIONAL FLUX + CHANGE IN STORGE - SOURCE = 0