<< Click to Display Table of Contents >> richards |
{ RICHARDS.PDE
A solution of Richards' equation in 1D.
Constant negative head at surface, unit gradient at bottom.
This problem runs slowly, because the very steep wave front
requires small cells and small timesteps to track accurately.
submitted by Neil Soicher of University of Hawaii.
}
title "1-D Richard's equation"
coordinates
cartesian1('y')
variables
h (1)
definitions
thr = 0.2
ths = 0.58
alpha = .08
n = 1.412
ks = 10
{Using Van Genuchten parameters for water content (wc),
water capacitance (C=d(wc)/dh), effecive saturation (se),
and hydraulic Conductivity (k) }
m = 1-1/n
wc = if h<0 then thr+(ths-thr)*(1+(abs(alpha*h))^n)^(-m) else ths
C = ((1-n)*abs(-alpha*h)^n*(1+abs(-alpha*h)^n)^((1/n)-2)*(ths-thr))/h
se = (wc-thr)/(ths-thr)
k = ks*sqrt(se)*(1-(1-se^(1/m))^m)^2
initial values
h = 199*exp(-(y-100)^2)-200
equations
h : dy(k*(dy(h)+1)) = C*dt(h)
boundaries
region 1
start(0)
line to (100) point value(h) = -1
front(h+150,1)
time 0 to 2
monitors
for cycle=10
elevation(h) from (0) to (100) as "pressure"
elevation(c) from (0) to (100) as "capacitance"
elevation(k) from (0) to (100) log as "conductivity"
grid(y)
plots
for t=0.001 by 0.001 to 0.01 0.1 by 0.1 to endtime
elevation(h) from (0) to (100) as "pressure"
elevation(c) from (0) to (100) as "capacitance"
elevation(k) from (0) to (100) log as "conductivity"
grid(y)
history(K) at (90) (95) (99) (100)
history(C) at (90) (95) (99) (100)
end "IE3vuxq/blOMIRLitV++FYgmXZuPz8D1+wvzXgpATJSsnTmsgWgSZOaLi+YOaMBsdxHjOXQUBxLPMVWceTZ+tzU0r6xbZ0Y9YaEBD8IG48nPCNezshKtEYPOYKh4ucdlKJWvPO8XzbScXAA9wKDowS86YuXXbtMkiY/S4U2KCpa"