Dear all,
i set a problem with a cylindrical2D grid, it is composed, in the
vertical direction of two layers, in the lower one there is a velocity
from the center well to the outside, in the upper, the velocity is
negligible. A concentration is injected in a well located at the center.
as it is in a porou medium, there is dispersion, which is proportional
to velociaty. This dispersion has two different compounds in the radial
and in the vertical direction. For this, i use two cellVariable (D2 and
D2) and i set the coeff in the diffusion term as (D1.faceValue,D2.faceValue)
The BC at the outlet for the concentration is a null gradient
I tried two BC for the injection well : 1, where=mesh.faceLeft and
where=[0] it does not change
my problem :
when i run the solve command, iget this :
ValueError: The truth value of an array with more than one element is
ambiguous. Use a.any() or a.all()
the file is enclosed
thanks
olivier
# -*- coding: utf-8 -*-
"""
Created on Tue Sep 23 07:33:14 2014
@author: olive
"""
""" a model to calculate radial injeciton of heat and a non reactive tracer
in a 2 layer domain
data from hecht mendez"""
from fipy import *
from fipy.terms.cellTerm import CellTerm
rho_s = 2650
c_s = 880
rho_w = 1e3
c_w = 4.18e3
lb = 2 #W/m/K
n=0.3
R = (n*rho_w*c_w+(1-n)*rho_s*c_s)/(n*rho_w*c_w)
D_t = lb/(n*rho_w*c_w)
dff=1e-10*86400 # diffusion coeff m2/d
a_L = .1
a_Z = .01
# domain units m days
a=logspace(-2,0,50);a=r_[0,a];dr=a[1:]-a[:-1];r=a[1:]/2+a[:-1]/2;nr=len(r)
b1=5.;b2=2. # thickness of geol layers bottom(1) fast, top(2) slow
a=logspace(-2,-1,20);dy0=a/sum(a)
dy=r_[b1*dy0[-1::-1],b2*dy0]
y=cumsum(dy)-dy/2;ny=len(y)
mesh=CylindricalGrid2D(dx=dr,nx=nr,dy=dy,ny=ny)
R,Y=mesh.getCellCenters()
Q=120 # m3/j or 5 m3/h
v=Q/(2*pi*r*b1*n) #m/d
a,v1=meshgrid(v,ones(ny));v1=ravel(v1)
V = CellVariable(name="velocity",mesh=mesh,value=v1)
V.setValue(v1*1e-4,where=Y>b1)
C = CellVariable(name ="Tracer",mesh=mesh,value=0.,hasOld=True)
#C.faceGrad.constrain(V.faceValue-V.faceValue*C.faceValue,
where=mesh.facesRight)
C.constrain(1.,where=[0])
D1 = CellVariable(name="dispers",mesh=mesh,value=V.value*a_L+dff)
D2 = CellVariable(name="dispers",mesh=mesh,value=V.value*a_Z+dff)
eq = TransientTerm() == ImplicitDiffusionTerm(coeff =
(D1.faceValue,D2.faceValue)) - \
PowerLawConvectionTerm(coeff = (V,0.))
dt = .1
for i in range(25):
C.updateOld()
eq.solve(var=C,dt=dt)_______________________________________________
fipy mailing list
[email protected]
http://www.ctcms.nist.gov/fipy
[ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
