mbakker7 / timml

An analytic element model for steady multi-layer flow
MIT License
38 stars 21 forks source link

changetrace in stripareasink needs to be fixed for water to leave aquifer during tracing #36

Open jentjr opened 4 years ago

jentjr commented 4 years ago
# set up the aquifer with a hydraulic conductivity of 30 ft/day, thickness of 50 ft, and a porosity of 0.25
ml = timml.ModelMaq(kaq=30, z=[50, 0], npor=0.25)
# Constant head boundary
rf = timml.Constant(ml, xr=0, yr=5000, hr=100, layer=0)
# set the hydraulic gradient to 0.005
uf = timml.Uflow(ml, slope=0.005, angle=0)
# set the recharge to 0.0027 ft/day
p = timml.StripAreaSink(ml, xleft=0, xright=5000, N=0.0027, layer=0)
# extraction wells
w1 = timml.Well(ml, xw=3500, yw=2500, Qw=4000, rw=0.25, layers=0, label='Well 1')
w2 = timml.Well(ml, xw=3400, yw=2700, Qw=4000, rw=0.25, layers=0, label='Well 2')
w3 = timml.Well(ml, xw=3400, yw=2300, Qw=4000, rw=0.25, layers=0, label='Well 3')

ml.solve()

# Plot a 5 year capture zone
ml.contour(win=[0, 5000, 0, 5000], ngr=50, layers=0, figsize=(10,10), labels=True)
w1.plotcapzone(nt=20, zstart=0, tmax=5*365.25, color='r')
w2.plotcapzone(nt=20, zstart=0, tmax=5*365.25)
w3.plotcapzone(nt=20, zstart=0, tmax=5*365.25)
---------------------------------------------------------------------------
AttributeError                            Traceback (most recent call last)
<ipython-input-55-7632c4335bcc> in <module>
      1 # Plot a 5 year capture zone
      2 ml.contour(win=[0, 5000, 0, 5000], ngr=50, layers=0, figsize=(10,10), labels=True)
----> 3 w1.plotcapzone(nt=20, zstart=0, tmax=5*365.25, color='r')
      4 w2.plotcapzone(nt=20, zstart=0, tmax=5*365.25)
      5 w3.plotcapzone(nt=20, zstart=0, tmax=5*365.25)

/usr/local/src/timml/timml/well.py in plotcapzone(self, nt, zstart, hstepmax, vstepfrac, tmax, nstepmax, silent, color, orientation, win, newfig, figsize, return_traces, metadata)
    244             nstepmax=nstepmax, silent=silent, color=color,
    245             orientation=orientation, win=win, newfig=newfig,
--> 246             figsize=figsize, return_traces=return_traces, metadata=metadata)
    247         if return_traces:
    248             return traces

/usr/local/src/timml/timml/util.py in tracelines(self, xstart, ystart, zstart, hstepmax, vstepfrac, tmax, nstepmax, silent, color, orientation, win, newfig, figsize, return_traces, metadata)
    219                 win=win,
    220                 returnlayers=True,
--> 221                 metadata=metadata,
    222             )
    223             if return_traces:

/usr/local/src/timml/timml/trace.py in timtraceline(ml, xstart, ystart, zstart, hstepmax, vstepfrac, tmax, nstepmax, win, silent, returnlayers, metadata)
     54         layerlist.append(modellayer)
     55         v0 = (
---> 56             ml.velocomp(x0, y0, z0, aq, [layer, ltype]) * direction
     57         )  # wordt nog gebruikt
     58         vx, vy, vz = v0

/usr/local/src/timml/timml/model.py in velocomp(self, x, y, z, aq, layer_ltype)
    274             qztop = qzlayer[layer]
    275             if layer == 0:
--> 276                 qztop += self.qztop(x, y)
    277             vz = (qzbot + (z - aq.zaqbot[layer]) / aq.Haq[layer] * \
    278                  (qztop - qzbot)) / aq.nporaq[layer]        

/usr/local/src/timml/timml/model.py in qztop(self, x, y, aq)
    112         if aq.ltype[0] == 'a':  # otherwise recharge cannot be added
    113             for e in aq.elementlist:
--> 114                 rv += e.qztop(x, y)
    115         return rv
    116 

/usr/local/src/timml/timml/stripareasink.py in qztop(self, x, y)
    161     def qztop(self, x, y):
    162         rv = 0.0
--> 163         if np.sqrt((x - self.xc) ** 2 + (y - self.yc) ** 2) <= self.R:
    164             rv = -self.parameters[0, 0]  # minus cause the parameter is the infiltration rate
    165         return rv

AttributeError: 'StripAreaSink' object has no attribute 'R'
mbakker7 commented 4 years ago

There are two problems here:

  1. The error in qztop, which should be fixed.
  2. StripAreaSink elements can not, in general, be combined with two-dimensional elements (wells) without the addition of line-elements along the boundary. The docs should be improved and an example should be added.
mbakker7 commented 4 years ago

qztop is fixed in https://github.com/mbakker7/timml/commit/15ead67fb8eb84f3b91ca5dbe657231fe7f87f6f

changetrace needs to be fixed, and an example needs to be added