Solving BL in python

This is another try to solve BL equation in FiPy. I gave up some time ago without actually giving it a proper trial. Now I'm going to try it out here with the new version of FiPy that I have installed in Conda.



In [1]:

    
phi = 0.4 # porosity
nw = 2.0
no = 2.0
krw0 = 0.4
kro0 = 1.0
swc = 0.15
sor = 0.2
mu_w = 0.001 # Pa.s
mu_o = 0.002 # Pa.s



In [ ]:

    
def kro(sw, kro0, sor, swc, no):
    res = ((swc<=sw) & (sw<=1-sor))*kro0*((1-sw-sor)/(1-sor-swc))**no \
    +((0.0<sw) & (sw<swc))*(1+(kro0-1)/swc*sw) \
    +(sw>1-sor)*0.0 \
    +(sw<=0.0)*1.0
    return res

def krw(sw, krw0, sor, swc, nw):
    res = ((swc<=sw) & (sw<=1-sor))*krw0*((sw-swc)/(1-sor-swc))**nw \
    +((1-sor<sw) & (sw<1.0))*(-(1-krw0)/sor*(1.0-sw)+1.0) \
    +(sw<=swc)*0.0 \
    +(sw>=1.0)*1.0
    return res

def dkrodsw(sw, kro0, sor, swc, no):
    res = ((swc<=sw) & (sw<=1-sor))*(-no*kro0/(1-sor-swc)*((1-sw-sor)/(1-sor-swc))**(no-1)) \
    +((0.0<sw) & (sw<swc))*(kro0-1)/swc \
    +((sw>1-sor) | (sw<=0.0))*0.0
    return res

def dkrwdsw(sw, krw0, sor, swc, nw):
    res = ((swc<=sw) & (sw<=1-sor))*nw*krw0/(1-sor-swc)*((sw-swc)/(1-sor-swc))**(nw-1) \
    +((1-sor<sw) & (sw<1.0))*(1-krw0)/sor \
    +((sw<swc) | (sw>=1.0))*0.0
    return res



In [2]:

    
import fipy



In [ ]:

    
nx = 50
L = 0.12 # domain length
u_inj = 1.0e-5 # m/s injection velocity
mesh = fipy.Grid1D(nx = nx, Lx = L)
sw_inj = 1.0 # left boundary condition
sw_init = swc
sw0 = swc

p_out = 1.0e5 # right boundary pressure
p_in  = 2.0e5 # left boundary pressure (making it simple)
p0 = p_out

sw = fipy.CellVariable(name="saturation", mesh=mesh, value=sw0, hasOld=True)
p  = fipy.CellVariable(name="pressure", mesh=mesh, value=p0, hasOld=True)



In [3]:

    
? fipy.Grid1D