Added indexing of pores

gen/generate.py

@@ -15,7 +15,7 @@ def main():
     print("Min. number of throats per pore: %d" % min(pore_throat_counts))
     print("Avg. number of throats per pore: %f" % (sum(pore_throat_counts) / len(pore_throat_counts)))
     
-    target  = [2, 2, 2]
+    target  = [3, 3, 3]
     cutoff  = 0.5 * max([basenet.ub[i] - basenet.lb[i] \
                         for i in range(len(basenet.ub))])
     

gen/pool.py

@@ -3,13 +3,15 @@ import multiprocessing as mp
 class Dummy:
     def __init__(self, n):
         self.vals = [i for i in range(n)]
-    def square(self, i):
-        return self.vals[i]**2
+    def square(self, i, offset, values):
+        return self.vals[i]**2 + offset + max(values)
 
 def main():
     dummy = Dummy(10)
+    offset = 1
+    values = [-1, 0, 3, 1]
     with mp.Pool(processes=4) as p:
-        res = p.map(func=dummy.square, iterable=[i for i in range(6)])
+        res = p.starmap(func=dummy.square, iterable=[(i, offset, values) for i in range(6)])
         print(res)
 
 if __name__ == '__main__':

netflow/netgen.py

@@ -12,7 +12,7 @@ class Pore:
     """Class of a pore connected to other pores via throats."""
     id = 0
     def __init__(self, pos: List[float], r: float, label: str = LABELS[0], \
-        throats: Set = None, id: int = -1):
+        throats: Set = None, index: int = -1, id: int = -1):
         if id == -1: # no id given
             self.id = Pore.id; Pore.id = Pore.id+1
         else: # id provided
@@ -22,6 +22,8 @@ class Pore:
         self.throats = throats # throats set
         self.r = r # radius
         self.label = label # label string like '', 'in', or 'out'
+        # Identifier for pores in list used in generate_dendrogram
+        self.index = index
     def __repr__(self): return str(self.__class__) + ': ' + \
         str({k:self.__dict__[k] for k in self.__dict__ if k != 'throats'}) + \
         ' {0:d} throats'.format(len(self.throats)) # dont flush throat objects
@@ -105,7 +107,8 @@ class CellList:
     
     find_candidates() assumes pore is located in interior of domain (non-periodic)
     """
-    def __init__(self, pores, domainSize: List[float], cellSize: float, Lmax: float):
+    def __init__(self, pores: List[Pore], domainSize: List[float], \
+                         cellSize: float, Lmax: float):
         from math import ceil
         self.dim = len(domainSize)
         self.Lmax = Lmax
@@ -164,6 +167,19 @@ class CellList:
                             candidates[nbor] = l
         return candidates
     
+    def find_best_candidate_in_cell(self, cellIdx: int, porePos: List[float],
+                        poreIdx: int, throatLen: float) -> Tuple[int, float]:
+        best_candidate_index = -1
+        best_length = self.Lmax
+        for nbor in self.poresSorted[cellIdx]:
+            if (nbor.index == poreIdx): continue
+            l = distance(nbor.pos, porePos)
+            diff = abs(l - throatLen)
+            if (l < self.Lmax and diff < abs(best_length - throatLen)):
+                best_candidate_index = nbor.index
+                best_length = l
+        return (best_candidate_index, best_length)
+
     def cell_stats(self):
         cellCounts = [len(self.poresSorted[i]) for i in range(self.totalCells)]
         n    = sum(cellCounts)
@@ -561,7 +577,7 @@ def generate_dendrogram(basenet: Network, targetsize: List[int], \
     Periodic throats are marked with a periodicity label, e.g., '1 0 -1' in 3d.
     If cutoff == float('nan'), pores are uniformly distributed.
     """
-    from random import seed, random, randint
+    from random import seed, random, randint, shuffle
     from itertools import product
     seed(sd)
     d = len(targetsize) # number of spatial dimensions
@@ -585,11 +601,13 @@ def generate_dendrogram(basenet: Network, targetsize: List[int], \
     if (cutoff != cutoff): # uniform
         print("\b"*21 + "uniform pore distribution")
         # loop over network sections & add uniformly distributed pores drawn from basenet
+        runningIndex = 0
         for s in product(*[range(k) for k in targetsize]):
             for k, pore in enumerate(basepores):
                 pos = [random()*l for l in L]
                 pores.append(Pore(pos=pos, r=pore.r, label=LABELS[0], 
-                    throats=pore.throats.copy()))
+                    throats=pore.throats.copy(), index=runningIndex))
+                runningIndex += 1
     else: # dendrogram-based
         # extract cluster hierarchy from basenet
         centroids = [pore.pos for pore in basepores]
@@ -626,12 +644,16 @@ def generate_dendrogram(basenet: Network, targetsize: List[int], \
             for k, pore in enumerate(basepores):
                 pos = [c-lb + si*(ub-lb) \
                     for c,si,lb,ub in zip(centroids[k],s,basenet.lb,basenet.ub)]
-                pores.insert(randint(0,len(pores)), Pore(pos=pos,
-                    r=pore.r, label=LABELS[0], throats=pore.throats.copy()))
+                pores.append(Pore(pos=pos, r=pore.r, label=LABELS[0],
+                    throats=pore.throats.copy()))
+        # Finally, shuffle pores (random order)
+        shuffle(pores)
+                
         print("\b"*21 + "left {0:d} of {1:d} clusters (incl. {2:d} pores) untouched".\
             format(sum([int(not j) for j in touched]), len(touched), len(basepores)))
-        # flip pores outside back into domain
-        for pore in pores:
+        # flip pores outside back into domain and set respective index of all pores
+        for i, pore in enumerate(pores):
+            pore.index = i
             for k in range(d):
                 if pore.pos[k] < 0:
                     pore.pos[k] = L[k] + pore.pos[k]
@@ -652,14 +674,8 @@ def generate_dendrogram(basenet: Network, targetsize: List[int], \
     throats = []
     k = 0 # count of unrealised throats
     
-    """
     # DEBUG
-    import sys
-    max_candidates = 0
-    min_candidates = sys.maxsize
-    avg_candidates = 0
-    """
-    # Spawn parallel 
+    #max_throat_diff = 0.
     for i, pore in enumerate(pores[:n]):
         # (re)run triangulation if maximal number of neighbors exceeds
         # threshold given by 2x init. maximum (large nnbrs -> expensive search)
@@ -672,16 +688,6 @@ def generate_dendrogram(basenet: Network, targetsize: List[int], \
         # TODO: Search neighboring cells in parallel
         candidates = cellList.find_candidates(pore)
         
-        """
-        # DEBUG
-        nCandidates = len(candidates)
-        
-        avg_candidates += nCandidates
-        if (nCandidates > max_candidates):
-            max_candidates = nCandidates
-        if (nCandidates < min_candidates):
-            min_candidates = nCandidates
-        """
         # establish throat connections between pore and best candidate
         while (len(pore.throats) > 0):
             """
@@ -696,6 +702,11 @@ def generate_dendrogram(basenet: Network, targetsize: List[int], \
             if (len(candidates) == 0):
                 k = k+1; continue # no candidates left, give up
             nbor, lt = min(candidates.items(), key=lambda pl: abs(pl[1] - lt))
+            """
+            diff = abs(lt - ltc)
+            if (max_throat_diff < diff):
+                max_throat_diff = diff
+            """
             # Remove nbor from candidates
             del candidates[nbor]
             nborc = nbor # memorize potential periodic copy
@@ -727,14 +738,8 @@ def generate_dendrogram(basenet: Network, targetsize: List[int], \
         
     print("\b"*24 + "{0:d} throats in total, {1:d} unrealised".\
         format(len(throats), k))
-    """
-    # DEBUG
-    avg_candidates /= n
-    print("Candidates statistics")
-    print("Max. candidates per pore: %d" % max_candidates)
-    print("Min. candidates per pore: %d" % min_candidates)
-    print("Avg. candidates per pore: %f" % avg_candidates)
-    """
+    #print("Max. throat length difference: %e" % max_throat_diff)
+    
     # assemble and return network
     network = Network(lb=[0.0 for k in range(d)],
         ub=L, Lmax=basenet.Lmax, label='from_' + basenet.label)
@@ -757,27 +762,30 @@ def __add_buffer_layers(pores: List[Pore], targetsize: List[float],
     for k, pore in enumerate(pores): pore.label = str(k) + d * ' 0'
     # add pore buffer layers
     copies = {pore:set() for pore in pores} # {pore:{set of periodic copies}}
+    runningIndex = n
     for k in range(d):
         # left bound
         pores_layer = [pore for pore in pores \
             if targetsize[k]-Lbuffer <= pore.pos[k] < targetsize[k]]
         for pore in pores_layer:
-            pcopy = Pore(pore.pos, pore.r, throats=pore.throats, id=pore.id)
+            pcopy = Pore(pore.pos, pore.r, throats=pore.throats, index=runningIndex, id=pore.id)
             pcopy.pos[k] = pcopy.pos[k]-targetsize[k]
             lbl = pore.label.split(); lbl[k + 1] = '-1'
             pcopy.label = ' '.join(lbl) # periodicity label
             pores.append(pcopy)
             copies[pores[int(lbl[0])]].add(pcopy)
+            runningIndex += 1
         # right bound
         pores_layer = [pore for pore in pores \
             if 0 <= pore.pos[k] < Lbuffer]
         for pore in pores_layer:
-            pcopy = Pore(pore.pos, pore.r, throats=pore.throats, id=pore.id)
+            pcopy = Pore(pore.pos, pore.r, throats=pore.throats, index=runningIndex, id=pore.id)
             pcopy.pos[k] = pcopy.pos[k]+targetsize[k]
             lbl = pore.label.split(); lbl[k + 1] = '1'
             pcopy.label = ' '.join(lbl) # periodicity label
             pores.append(pcopy)
             copies[pores[int(lbl[0])]].add(pcopy)
+            runningIndex += 1
     # reset labels of pores inside domain (without buffer layers)
     for pore in pores[:n]: pore.label = LABELS[0]
     return copies