#------------------------------------------------------------------------------- # Name: wellfun-Theis.py # Purpose: Calculate groundwater lowering by abstraction from confined aquifer # Source: Huisman, 1972. Groundwater Recorvery. Macmillan. # Created: 21-01-2012 # Copyright: (c) Grondwaterformules.nl #------------------------------------------------------------------------------- #!/usr/bin/env python from scipy import * from scipy import special from scipy.special import expn as ExpInt from pylab import * def Theis(r,t,Q0,kD,S): u2 = (S*r**2)/(4*kD*t) s = Q0/(4*pi*kD)*ExpInt(1,u2) return s def plot_s(r,t,s): # Define colors orange = '#FF9900' blue = '#2F64B2' # Create plots graph = subplot(111) # Set ticklines xticklines = graph.get_xticklines() yticklines = graph.get_yticklines() for line in xticklines+yticklines: line.set_linewidth(3) for line in xticklines: line.set_visible(False) # Set gridlines gridlines = graph.get_xgridlines() gridlines.extend( graph.get_ygridlines() ) for line in gridlines: line.set_linestyle('-') grid(True) # Set ticklabels xticklabels = graph.get_xticklabels() yticklabels = graph.get_yticklabels() for label in xticklabels+yticklabels: label.set_color('k') label.set_fontsize('medium') label.set_visible(True) # Create annotation ylabel('verlaging (m)') xlabel('afstand tot de put (m)') # Plot lowring in time for ti in range(len(t)): si = s[ti,:] graph.plot(r, -si, color=blue, lw=2.0) # Set axes limits # If you see a plot without lines, delete or comment next line: graph.set_ylim(-5, 0) # Filename filename = "Theis.png" # Save figure to file fig = figure(1) mydpi = 100.0 fig.set_dpi(mydpi) fig.set_size_inches(6.0,4.0) fig.savefig(filename,dpi=mydpi, facecolor='w', edgecolor='w') # Show the plot on screen show() if __name__ == '__main__': r = arange(0.1,1000.0,0.1) # r is een numpy array #r = array([0.1,10.0,100.0,500.0,1000.0]) t = array([0,1,2,3,4,5]) t = 10**t Q0 = 1000000.0/365.0 kD = 1500.0 S = 0.0001 # Calculate lowering in time (Hantush) s = zeros((len(t),len(r))) for ti in range(len(t)): for ri in range(len(r)): s[ti,ri] = Theis(r[ri],t[ti],Q0,kD,S) # Plot lowering of groundwater plot_s(r,t,s)