PMEAL / OpenPNM

A Python package for performing pore network modeling of porous media
http://openpnm.org
MIT License
444 stars 174 forks source link

effect of pore size on relative permeabilities #677

Closed JiyuanZhang221 closed 7 years ago

JiyuanZhang221 commented 7 years ago

I was assuming the size of pores/throats should have negligible effect on relative permeabilities. Unfortunately, I found the simulated relative permeabilities are identical when inputting different pore sizes following the example "Example: Regenerating Data from J.T. Gostick et al. / JPS 173 (2007) 277–290". Does anybody have a clue about the reasons?

-- coding: utf-8 --

import OpenPNM as op

workspace = op.Base.Workspace() workspace.clear()

N=1 nx=15 #No. of pores in x-direction ny=15 #No. of pores in y-direction nz=15 #No. of pores in z-direction I=3 # Scaling factor used for controlling pore size Lc = 0.005 # Lattice spacing

1 setting up network

sgl = op.Network.Cubic([nx, ny, nz], spacing=Lc, connectivity=6,name='relativeperm2') sgl.add_boundaries()

2 set up geometries

Ps = sgl.pores('boundary') Ts = sgl.find_neighbor_throats(pores=Ps,mode='not_intersection') boun = op.Geometry.Boundary(network=sgl,pores=Ps,throats=Ts,name='boun')

Ps = sgl.pores('boundary',mode='difference') Ts = sgl.find_neighbor_throats(pores=Ps,mode='intersection',flatten=True) geo = op.Geometry.Cube_and_Cuboid(network=sgl,pores=Ps,throats=Ts,name='geo')

geo.models.add(propname='pore.diameter', model=op.Geometry.models.pore_diameter.normal, shape=1.0, scale=0.00001, loc=0.00006*I,seeds='pore.seed') pore_dia=geo['pore.diameter']
geo.models.add(propname='throat.diameter', model=op.Geometry.models.throat_misc.neighbor, pore_prop='pore.diameter',mode='mean')

geo.models.add(propname='throat.length', model=op.Geometry.models.throat_length.straight)

geo.models.add(propname='throat.area', model=op.Geometry.models.throat_area.cuboid)

geo.models.add(propname='pore.area', model=op.Geometry.models.pore_area.cubic)

geo.models.add(propname='pore.volume', model=op.Geometry.models.pore_volume.cube)

geo.models.add(propname='throat.volume', model=op.Geometry.models.throat_volume.cuboid)

boun.regenerate()

air = op.Phases.Air(network = sgl, name = 'air') water = op.Phases.Water(network = sgl, name = 'water')

reset pore contact angle

water['pore.contact_angle'] = 120 air['throat.contact_angle']=110 air['throat.surface_tension']=.072 air['pore.surface_tension']=0.072

1 create physics objects associated with our phases

Ps = sgl.pores() Ts = sgl.throats() phys_water = op.Physics.Standard(network=sgl,phase=water,pores=Ps,throats=Ts,dynamic_data=True,name='standard_water_physics') phys_air = op.Physics.Standard(network=sgl,phase=air,pores=Ps,throats=Ts,dynamic_data=True,name='standard_air_physics')

2 calculating physics properties (capillary pressure, hydraulic conductance, etc)

phys_water.regenerate() phys_air.regenerate()

inlets = sgl.pores('top_boundary') used_inlets = [inlets[x] for x in range(0, len(inlets), 2)]

using every other pore in the bottom and boundary as an inlet

prevents extremely small diffusivity and permeability values in the z direction

used_inlets = [inlets[x] for x in range(0, len(inlets), 2)] OP_1 = op.Algorithms.OrdinaryPercolation(network=sgl,invading_phase=water,defending_phase=air)

OP_1.run(inlets = used_inlets,npts=50) OP_1.run() OP_1.evaluate_trapping(p_outlets=sgl.pores('bottom_boundary')) sw_eva=OP_1.evaluate_late_pore_filling(Pc=500)

OP_1.evaluate_late_throat_filling(Pc=1000)

OP_1.plot_drainage_curve() OP_1.plot_primary_drainage_curve()

sat = [] perm_air = {'0': [], '1': [], '2': []} diff_air = {'0': [], '1': [], '2': []} perm_water = {'0': [], '1': [], '2': []} diff_water = {'0': [], '1': [], '2': []}

max_inv_seq = max(OP_1['throat.inv_seq'])

num_seq = 10

pore_oc=(num_seq+1)[len(Ps)[0]] throat_oc=(num_seq+1)[len(Ts)[0]] sat1=(num_seq+1)*[0]

netname1='Stokes_alg_single_phase_air' netname2='Stokes_alg_single_phase_water' netname3='Stokes_alg_multi_phase_air' netname4='Stokes_alg_multi_phase_water'

for x in range(num_seq+1):

netname0=str(x)
netname1=netname1+netname0
netname2=netname2+netname0
netname3=netname3+netname0
netname4=netname4+netname0

OP_1.return_results(sat = x/num_seq)

#printing out so we know how far along we are
print('seq = '+str(round(max_inv_seq*(x/num_seq)))+' Seq out of '+str(round(max_inv_seq))+' total sequences')

final_pores = water['pore.occupancy']
pore_volumes = sgl['pore.volume']
final_throats = water['throat.occupancy']
throat_volumes = sgl['throat.volume']
pore_oc[x]=final_pores
throat_oc[x]=final_throats
saturation = (sum(final_pores*pore_volumes) + sum(final_throats*throat_volumes))/(sum(pore_volumes) + sum(throat_volumes))
sat.append(saturation)

#adding multiphase conductances
phys_air.add_model(model=op.Physics.models.multiphase.conduit_conductance,
           propname='throat.conduit_hydraulic_conductance',
           throat_conductance='throat.hydraulic_conductance')
phys_water.add_model(model=op.Physics.models.multiphase.conduit_conductance,
           propname='throat.conduit_hydraulic_conductance',
           throat_conductance='throat.hydraulic_conductance')

phys_air.add_model(model=op.Physics.models.multiphase.late_throat_filling,
           propname='throat.capillary_pressure',
           Pc=5000,throat_entry_pressure='throat.capillary_pressure')

phys_water.add_model(model=op.Physics.models.multiphase.late_throat_filling,
           propname='throat.capillary_pressure',
           Pc=5000,throat_entry_pressure='throat.capillary_pressure')

phys_air.add_model(model=op.Physics.models.multiphase.late_pore_filling,
           propname='throat.capillary_pressure',
           Pc=5000,throat_entry_pressure='throat.capillary_pressure')

phys_water.add_model(model=op.Physics.models.multiphase.late_pore_filling,
           propname='throat.capillary_pressure',
           Pc=5000,throat_entry_pressure='throat.capillary_pressure')
#bounds = [['front', 'back'], ['left', 'right'], ['top', 'bottom']]
bounds=[['top', 'bottom']]

for bound_increment in range(len(bounds)):
   #single phase
    Stokes_alg_single_phase_air = op.Algorithms.StokesFlow(name=netname1,network=sgl,phase=water)
    BC2_pores = sgl.pores(labels=bounds[bound_increment][1]+'_boundary')
    BC1_pores = sgl.pores(labels=bounds[bound_increment][0]+'_boundary')
    Stokes_alg_single_phase_air.set_boundary_conditions(bctype='Dirichlet',bcvalue=.6,pores=BC1_pores)
    Stokes_alg_single_phase_air.set_boundary_conditions(bctype='Dirichlet',bcvalue=.2,pores=BC2_pores)
    Stokes_alg_single_phase_air.run(conductance = 'hydraulic_conductance')
    effective_permeability_air_single = Stokes_alg_single_phase_air.calc_eff_permeability()
    Stokes_alg_multi_phase_air = op.Algorithms.StokesFlow(name=netname3,network=sgl,phase=air)
    Stokes_alg_multi_phase_water = op.Algorithms.StokesFlow(name=netname4,network=sgl,phase=water)
    Stokes_alg_multi_phase_air.set_boundary_conditions(bctype='Dirichlet',bcvalue=.6,pores=BC1_pores)
    Stokes_alg_multi_phase_water.set_boundary_conditions(bctype='Dirichlet',bcvalue=.6,pores=BC1_pores)
    Stokes_alg_multi_phase_air.set_boundary_conditions(bctype='Dirichlet',bcvalue=.2,pores=BC2_pores)
    Stokes_alg_multi_phase_water.set_boundary_conditions(bctype='Dirichlet',bcvalue=.2,pores=BC2_pores)
    #run algorithms with proper conduit conductance
    Stokes_alg_multi_phase_air.run(conductance = 'conduit_hydraulic_conductance')
    Stokes_alg_multi_phase_water.run(conductance = 'conduit_hydraulic_conductance')
    #calc effective properties
    effective_permeability_air_multi = Stokes_alg_multi_phase_air.calc_eff_permeability()
    effective_permeability_water_multi = Stokes_alg_multi_phase_water.calc_eff_permeability()
    relative_eff_perm_air = effective_permeability_air_multi/effective_permeability_air_single
    relative_eff_perm_water = effective_permeability_water_multi/effective_permeability_air_single
    perm_air[str(bound_increment)].append(relative_eff_perm_air)
    perm_water[str(bound_increment)].append(relative_eff_perm_water)

result=[sat,perm_air['0'], perm_water['0']]

jgostick commented 7 years ago

Hi @JiyuanZhang221, it's difficult to know why you're seeing this...it's quite possible a bug in the example, but may also be due to a bug in your script. For instance, how are you changing pores sizes, are you regenerating the geometry objects after you change the size distribution info.

Perhaps you could condense your script down to the bare minimum lines and share it here?

JiyuanZhang221 commented 7 years ago

Thanks @jgostick. Actually my primary target is to investigate the effect of pore/throat deformation at varying confining pressures on relative permeability. So I used a data file as an initial input to the pore sizes. Pore sizes were repeatedly updated according to certain rules and relative permeability were then re-calculated in a loop.

For consideration of simplicity, I attached a condensed version of the script for your reference (without calculation of deformation). In this condensed script, you can change the scaling factor (I) to see that the average pore sizes have very minor effect on relative permeability.

Thanks again.

jgostick commented 7 years ago

Hi @JiyuanZhang221 ... perhaps you could email me a script file that I could just run? I'm getting errors after cutting & pasting the above code. (jgostick at gmail)

xu-kai-xu commented 4 years ago

@jgostick , hi, I did an easy simulation by using two sets of pores to construct a pnm. The two pore diameters are 5e-5 m and 2.5 e-5 m, respectively. Then I use it to simulation single-phase flow. As intuitive feeling, I think the constant-pressure surface should not be vertical to the flow direction (z-axis) near the pore-size-changing areas. But the results show that the pressrue is the same at the plane cross the flow direction.

About the following two pictures:

image

visualization for the pore size.

image show the pressure and flow rate along y-axis.

It seems flow rate is ok. But I suspect the pressure result.

jgostick commented 4 years ago

I do not understand your question but the data you show above all looks good.

The pressure should be relatively constant at any z-position.

The flow should be 0 in every pore, except for boundary pores, because the mass balances to 0 (in - out). Your flow rate axis is scaled by 10^-22, so the red points are basically just showing numerical noise...the flows are indeed 0 as expected.

xu-kai-xu commented 4 years ago

Thank you. Still, I don't understand why the pore size will not change the pressure distribution. About the flow rate, I understand it.

ma-sadeghi commented 4 years ago

Thank you. Still, I don't understand why the pore size will not change the pressure distribution. About the flow rate, I understand it.

If you radically change the pore size distribution, you might see some noticeable change.

Try on a smaller network so the heterogeneity gets more pronounced.

jgostick commented 4 years ago

The only thing that pore size changes is the flow in the throats., not the pressure in the pores. There are several ways to think about this. Consider the pores at the same z-position...if there were a pressure difference between them, then flow would move between them, but there is no reason for flow to move laterally when the pressure is lower toward the outlet. Or, you can say that the flow rate through small throats is low, so no flow will go there, so therefore there won't be a pressure difference. Also, remember that the steady-state solution is at infinite time, so after a long time, the pressure will all pores will be equilibrated with it's neighbors. So if there is an 'instantaneous' low pressure somewhere, then fluid will rush in there to fill that location, and it's pressure will be close to it's neighbors.

xu-kai-xu commented 4 years ago

I guess I understand. I forget the steady-state condition. Thank you!

Then, I made some changes, the new pore size distribution is as follows:

p-cubic-half-5-f

From the figure, I can see pressure difference obviously in lateral direction. This may show the flow at the x-y plane.