The function grd.S()
computes the deficit S for all sub surface layers (not the overland flow layer)
grd.S(ddistx, X, ddist)
The output is a scalar
It follows the expression dS/dt = -dD/dt
from Skaugen,Peerebom and Nilsson(2015)
## Not run: ------------------------------------ # grd.S() ## ---------------------------------------------