MODFLOW

The governing partial differential equation for a confined aquifer used in MODFLOW is: :\frac{\partial}{\partial x} \left[ K_{xx} \frac{\partial h}{\partial x} \right] + \frac{\partial}{\partial y} \left[ K_{yy} \frac{\partial h}{\partial y} \right] + \frac{\partial}{\partial z} \left[ K_{zz} \frac{\partial h}{\partial z} \right] + W = S_{S} \frac{\partial h}{\partial t} where • K_{xx}, K_{yy} and K_{zz} are the values of hydraulic conductivity along the x, y, and z coordinate axes (L/T) • h is the potentiometric head (L) • W is a volumetric flux per unit volume representing sources and/or sinks of water, where negative values are extractions, and positive values are injections (T−1) • S_{S} is the specific storage of the porous material (L−1); and • t\, is time (T) Finite difference The finite difference form of the partial differential in a discretized aquifer domain (represented using rows, columns and layers) is: :\begin{align} & \mathit{CR}_{i,j-\tfrac{1}{2},k}\left(h^m_{i,j-1,k}-h^m_{i,j,k}\right) + \mathit{CR}_{i,j+\tfrac{1}{2},k}\left(h^m_{i,j+1,k}-h^m_{i,j,k}\right) + \\ & \mathit{CC}_{i-\tfrac{1}{2},j,k}\left(h^m_{i-1,j,k}-h^m_{i,j,k}\right) + \mathit{CC}_{i+\tfrac{1}{2},j,k}\left(h^m_{i+1,j,k}-h^m_{i,j,k}\right) + \\ & \mathit{CV}_{i,j,k-\tfrac{1}{2}}\left(h^m_{i,j,k-1}-h^m_{i,j,k}\right) + \mathit{CV}_{i,j,k+\tfrac{1}{2}}\left(h^m_{i,j,k+1}-h^m_{i,j,k}\right) + \\ & P_{i,j,k}\,h^m_{i,j,k} + Q_{i,j,k} = \mathit{SS}_{i,j,k}\left(\Delta r_j \Delta c_i \Delta v_k\right) \frac{h^m_{i,j,k}-h^{m-1}_{i,j,k}}{t^m-t^{m-1}} \end{align} where :h^m_{i,j,k}\, is the hydraulic head at cell i,j,k at time step m :CV, CR and CC are the hydraulic conductances, or branch conductances between node i,j,k and a neighboring node :P_{i,j,k}\, is the sum of coefficients of head from source and sink terms :Q_{i,j,k}\, is the sum of constants from source and sink terms, where Q_{i,j,k} is flow out of the groundwater system (such as pumping) and Q_{i,j,k}>0.0\, is flow in (such as injection) :\mathit{SS}_{i,j,k}\, is the specific storage :\Delta r_j\Delta c_i\Delta v_k\, are the dimensions of cell i,j,k, which, when multiplied, represent the volume of the cell; and :t^m\, is the time at time step m This equation is formulated into a system of equations to be solved as: :\begin{align} &\mathit{CV}_{i,j,k-\tfrac{1}{2}} h^m_{i,j,k-1} + \mathit{CC}_{i-\tfrac{1}{2},j,k} h^m_{i-1,j,k} + \mathit{CR}_{i,j-\tfrac{1}{2},k} h^m_{i,j-1,k} \\ &+ \left( - \mathit{CV}_{i,j,k-\tfrac{1}{2}} - \mathit{CC}_{i-\tfrac{1}{2},j,k} - \mathit{CR}_{i,j-\tfrac{1}{2},k} - \mathit{CR}_{i,j+\tfrac{1}{2},k} - \mathit{CC}_{i+\tfrac{1}{2},j,k} - \mathit{CV}_{i,j,k+\tfrac{1}{2}} + \mathit{HCOF}_{i,j,k}\right) h^m_{i,j,k} \\ &+ \mathit{CR}_{i,j+\tfrac{1}{2},k} h^m_{i,j+1,k} + \mathit{CC}_{i+\tfrac{1}{2},j,k} h^m_{i+1,j,k} + \mathit{CV}_{i,j,k+\tfrac{1}{2}} h^m_{i,j,k+1} = \mathit{RHS}_{i,j,k} \end{align} where :\begin{align} \mathit{HCOF}_{i,j,k} &= P_{i,j,k} - \frac{\mathit{SS}_{i,j,k}\Delta r_j \Delta c_i \Delta_k}{t^m-t^{m-1}} \\ \mathit{RHS}_{i,j,k} &= -Q_{i,j,k} - \mathit{SS}_{i,j,k}\Delta r_j \Delta c_i \Delta v_k \frac{h^{m-1}_{i,j,k}}{t^m-t^{m-1}} \end{align} or in matrix form as: :A\mathbf{h}=\mathbf{q} where :A is a matrix of the coefficients of head for all active nodes in the grid :\mathbf{h} is a vector of head values at the end of time step m for all nodes in the grid; and :\mathbf{q} is a vector of the constant terms, RHS, for all nodes of the grid. Limitations • The water must have a constant density, dynamic viscosity (and consequently temperature) throughout the modelling domain (SEAWAT is a modified version of MODFLOW which is designed for density-dependent groundwater flow and transport) • The principal components of anisotropy of the hydraulic conductivity used in MODFLOW is displayed on the right. This tensor does not allow non-orthogonal anisotropies, as could be expected from flow in fractures. Horizontal anisotropy for an entire layer can be represented by the coefficient "TRPY" (Data Item 3 Page 153). == Versions ==