Skip to content

Commit

Permalink
Merge pull request #1 from PMEAL/feature/UofT_work
Browse files Browse the repository at this point in the history 
Feature/uof t work
  • Loading branch information
@jhinebau
jhinebau committed Jul 10, 2014
2 parents 19bb653 + 1e19b1e commit c666e23
Show file tree
Hide file tree
Showing 22 changed files with 3,184 additions and 195 deletions.
64 changes: 33 additions & 31 deletions OpenPNM/Algorithms/__LinearSolver__.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -278,7 +278,7 @@ def rate(self,pores='',throats=''):
pores2[-sp.in1d(p1,pores)] = p1[-sp.in1d(p1,pores)]
X1 = self._result[pores1]
X2 = self._result[pores2]
g = self._conductance[throats]
g = self._fluid[self._conductance][throats]
R = sp.sum(sp.multiply(g,(X1-X2)))
return(R)

Expand All @@ -296,7 +296,7 @@ def _do_one_inner_iteration(self):
self._logger.info('Writing result to '+self.__class__.__name__+'[\''+self._conductance+'\']')
self[self._quantity] = self._result

def _calc_eff_prop(self):
def _calc_eff_prop(self,check_health=False):
try:
self[self._quantity]
except:
Expand All @@ -310,41 +310,43 @@ def _calc_eff_prop(self):
outlets = sp.where(self['pore.bcval_Dirichlet']==sp.amin(BCs))[0]

#Analyze input and output pores
#Check for coplanarity
# if self._net.iscoplanar(inlets) == False:
# raise Exception('The inlet pores do not define a plane. Effective property will be approximation')
# if self._net.iscoplanar(outlets) == False:
# raise Exception('The outlet pores do not define a plane. Effective property will be approximation')
#Ensure pores are on a face of domain (only 1 non-self neighbor each)
PnI = self._net.find_neighbor_pores(pores=inlets,mode='not_intersection',excl_self=True)
if sp.shape(PnI) != sp.shape(inlets):
self._logger.warning('The inlet pores have too many neighbors. Internal pores appear to be selected.')
PnO = self._net.find_neighbor_pores(pores=outlets,mode='not_intersection',excl_self=True)
if sp.shape(PnO) != sp.shape(outlets):
self._logger.warning('The outlet pores have too many neighbors. Internal pores appear to be selected.')
if check_health:
#Check for coplanarity
if self._net.iscoplanar(inlets) == False:
raise Exception('The inlet pores do not define a plane. Effective property will be approximation')
if self._net.iscoplanar(outlets) == False:
raise Exception('The outlet pores do not define a plane. Effective property will be approximation')
#Ensure pores are on a face of domain (only 1 non-self neighbor each)
PnI = self._net.find_neighbor_pores(pores=inlets,mode='not_intersection',excl_self=True)
if sp.shape(PnI) != sp.shape(inlets):
self._logger.warning('The inlet pores have too many neighbors. Internal pores appear to be selected.')
PnO = self._net.find_neighbor_pores(pores=outlets,mode='not_intersection',excl_self=True)
if sp.shape(PnO) != sp.shape(outlets):
self._logger.warning('The outlet pores have too many neighbors. Internal pores appear to be selected.')

#Fetch area and length of domain
#A = self._net.domain_area(face=inlets)
A = (40.5e-6*10)**2
L = self._net.domain_length(face_1=inlets,face_2=outlets)
x = self._result
#x = self._result
#Find flow through inlet face
Pin = []
Pn = []
for pore in inlets:
pore_value = x[pore]
neighbors = self._net.find_neighbor_pores(pore, excl_self = True)
for neighbor in neighbors:
neighbor_value = x[neighbor]
if(sp.absolute(neighbor_value - pore_value) > .000001):
Pin.append(pore)
Pn.append(neighbor)

Ts = self._net.find_connecting_throat(Pin,Pn)
g = self._fluid[self._conductance][Ts]
xin = x[Pin]
xout = x[Pn]
flow = g*(xin - xout)
#Pin = []
#Pn = []
#for pore in inlets:
# pore_value = x[pore]
# neighbors = self._net.find_neighbor_pores(pore, excl_self = True)
# for neighbor in neighbors:
# neighbor_value = x[neighbor]
# if(sp.absolute(neighbor_value - pore_value) > .000001):
# Pin.append(pore)
# Pn.append(neighbor)
#
#Ts = self._net.find_connecting_throat(Pin,Pn)
#g = self._fluid[self._conductance][Ts]
#xin = x[Pin]
#xout = x[Pn]
#flow = g*(xin - xout)
flow = self.rate(pores=inlets)
D = sp.sum(flow)*L/A/sp.absolute(BCs[0]-BCs[1])
return D

Expand Down
45 changes: 20 additions & 25 deletions OpenPNM/Algorithms/__OrdinaryPercolation__.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -22,21 +22,21 @@

class OrdinaryPercolation(GenericAlgorithm):
r"""
Simulates a capillary drainage experiment by looping through a list of capillary pressures
This function produces a pore_properties array called 'Pc_invaded' which contains the pressure at which a given pore was invaded. This list can be useful for reproducing the simulation for plotting or determining late pore filling.
Simulates a capillary drainage experiment by looping through a list of
capillary pressures.
Parameters
----------
loglevel : integer, optional
Level of the logger (10=Debug, 20=INFO, 30=Warning, 40=Error, 50=Critical)
loggername : string, optional
Set the name of the logger to be output on the console. Defaults to class name.
Note
----
Some info about OP?
network : OpenPNM Network Object
The network upon which the simulation will be run
name : string, optional
The name to assign to the Algorithm Object
Notes
-----
To run this algorithm, use 'setup()' to provide the necessary simulation
parameters, and then use 'run()' to execute it. Use 'update()' to send
the results of the simulation out of the algorithm object.
"""

def __init__(self, **kwargs):
Expand All @@ -46,14 +46,9 @@ def __init__(self, **kwargs):
super(OrdinaryPercolation,self).__init__(**kwargs)
self._logger.debug("Create Drainage Percolation Algorithm Object")

def setup(self,
invading_fluid = None,
defending_fluid = None,
inlets = [0],
npts = 25,
capillary_pressure = 'capillary_pressure',
AL=True,
**params):
def setup(self,invading_fluid = None,defending_fluid = None,inlets = [0],npts = 25,capillary_pressure = 'capillary_pressure',AL=True,**params):
r'''
'''
# Parse params
self._fluid_inv = invading_fluid
self._fluid_def = defending_fluid
Expand All @@ -64,6 +59,8 @@ def setup(self,
self._p_cap = capillary_pressure

def run(self):
r'''
'''
#See if setup has been run
try: capillary_pressure = self._p_cap
except:
Expand Down Expand Up @@ -94,8 +91,6 @@ def _do_outer_iteration_stage(self):
self._t_seq = sp.searchsorted(sp.unique(self._t_inv),self._t_inv)
self.set_pore_data(prop='inv_seq',data=self._p_seq)
self.set_throat_data(prop='inv_seq',data=self._t_seq)
#Remove temporary arrays and adjacency matrices
del self._net._adjacency_matrix['csr']['invaded']

def _do_one_inner_iteration(self,inv_val):
r"""
Expand All @@ -107,11 +102,10 @@ def _do_one_inner_iteration(self,inv_val):
#Generate a tlist containing boolean values for throat state
Tinvaded = self._t_cap<=inv_val
#Fill adjacency matrix with invasion state info
self._net.create_adjacency_matrix(data=Tinvaded,prop='invaded',sprsfmt='csr',dropzeros=True)
clusters = sprs.csgraph.connected_components(self._net._adjacency_matrix['csr']['invaded'])[1]
clusters = self._net.find_clusters(Tinvaded)
#Find all pores with at least 1 invaded throat (invaded)
Pinvaded = sp.zeros_like(clusters,dtype=bool)
nums = self._net.get_throat_indices('all')
nums = self._net.throats()
temp = self._net.find_connected_pores(nums)
temp = temp[Tinvaded]
temp = sp.hstack((temp[:,0],temp[:,1]))
Expand Down Expand Up @@ -207,6 +201,7 @@ def update(self, Pc=0, seq = None, occupancy='occupancy'):
temp = sp.array(t_inv,dtype=sp.float64,ndmin=1)
self._fluid_inv.set_throat_data(prop=occupancy,data=temp)
#Apply occupancy to defending fluid
if self._fluid_def != None:
temp = sp.array(~p_inv,dtype=sp.float64,ndmin=1)
self._fluid_def.set_pore_data(prop=occupancy,data=temp)
temp = sp.array(~t_inv,dtype=sp.float64,ndmin=1)
Expand Down
Loading

0 comments on commit c666e23

Please sign in to comment.