refactored the term 'face' to 'cell' to match what I've found in literature

This commit is contained in:
Stephen Mardson McQuay 2010-11-08 15:07:48 -07:00
parent 075dd586a3
commit 93fec67cc7
3 changed files with 85 additions and 85 deletions

View File

@ -8,7 +8,7 @@ import scipy.spatial
from interp.baker import run_baker from interp.baker import run_baker
from interp.tools import log from interp.tools import log
from simplex import face, contains from simplex import cell, contains
from smcqdelaunay import * from smcqdelaunay import *
@ -32,12 +32,12 @@ class grid(object):
if q: if q:
self.q = np.array(q) self.q = np.array(q)
self.faces = {} self.cells = {}
self.faces_for_vert = defaultdict(list) self.cells_for_vert = defaultdict(list)
def get_containing_simplex(self, X): def get_containing_simplex(self, X):
if not self.faces: if not self.cells:
raise Exception("face connectivity is not set up") raise Exception("cell connectivity is not set up")
# get closest point # get closest point
(dist, indicies) = self.tree.query(X, 2) (dist, indicies) = self.tree.query(X, 2)
@ -49,26 +49,26 @@ class grid(object):
log.debug('distance = %0.4f' % dist[0]) log.debug('distance = %0.4f' % dist[0])
simplex = None simplex = None
checked_faces = [] checked_cells = []
faces_to_check = list(self.faces_for_vert[closest_point]) cells_to_check = list(self.cells_for_vert[closest_point])
attempts = 0 attempts = 0
while not simplex and faces_to_check: while not simplex and cells_to_check:
attempts += 1 attempts += 1
# if attempts > 20: # if attempts > 20:
# raise Exception("probably recursing to many times") # raise Exception("probably recursing to many times")
cur_face = faces_to_check.pop(0) cur_cell = cells_to_check.pop(0)
checked_faces.append(cur_face) checked_cells.append(cur_cell)
if cur_face.contains(X, self): if cur_cell.contains(X, self):
simplex = cur_face simplex = cur_cell
continue continue
for neighbor in cur_face.neighbors: for neighbor in cur_cell.neighbors:
if (neighbor not in checked_faces) and (neighbor not in faces_to_check): if (neighbor not in checked_cells) and (neighbor not in cells_to_check):
faces_to_check.append(neighbor) cells_to_check.append(neighbor)
if not simplex: if not simplex:
raise Exception('no containing simplex found') raise Exception('no containing simplex found')
@ -137,21 +137,21 @@ class grid(object):
assert( len(self.verts) == len(self.q) ) assert( len(self.verts) == len(self.q) )
for c, i in enumerate(zip(self.verts, self.q)): for c, i in enumerate(zip(self.verts, self.q)):
r += "%d vert(%s): q(%0.4f)" % (c,i[0], i[1]) r += "%d vert(%s): q(%0.4f)" % (c,i[0], i[1])
face_str = ", ".join([str(f.name) for f in self.faces_for_vert[c]]) cell_str = ", ".join([str(f.name) for f in self.cells_for_vert[c]])
r += " faces: [%s]" % face_str r += " cells: [%s]" % cell_str
r += "\n" r += "\n"
if self.faces: if self.cells:
for v in self.faces.itervalues(): for v in self.cells.itervalues():
r += "%s\n" % v r += "%s\n" % v
return r return r
class delaunay_grid(grid): class delaunay_grid(grid):
face_re = re.compile(r''' cell_re = re.compile(r'''
-\s+(?P<face>f\d+).*? -\s+(?P<cell>f\d+).*?
vertices:\s(?P<verts>.*?)\n.*? vertices:\s(?P<verts>.*?)\n.*?
neighboring\s faces:\s+(?P<neigh>[\sf\d]*) neighboring\s cells:\s+(?P<neigh>[\sf\d]*)
''', re.S|re.X) ''', re.S|re.X)
point_re = re.compile(r''' point_re = re.compile(r'''
@ -168,7 +168,7 @@ class delaunay_grid(grid):
def get_containing_simplex(self, X): def get_containing_simplex(self, X):
if not self.faces: if not self.cells:
self.construct_connectivity() self.construct_connectivity()
# get closest point # get closest point
@ -181,24 +181,24 @@ class delaunay_grid(grid):
log.debug('distance = %0.4f' % dist[0]) log.debug('distance = %0.4f' % dist[0])
simplex = None simplex = None
checked_faces = [] checked_cells = []
faces_to_check = self.faces_for_vert[closest_point] cells_to_check = self.cells_for_vert[closest_point]
attempts = 0 attempts = 0
while not simplex and faces_to_check: while not simplex and cells_to_check:
attempts += 1 attempts += 1
# if attempts > 20: # if attempts > 20:
# raise Exception("probably recursing to many times") # raise Exception("probably recursing to many times")
cur_face = faces_to_check.pop(0) cur_cell = cells_to_check.pop(0)
checked_faces.append(cur_face) checked_cells.append(cur_cell)
if cur_face.contains(X, self): if cur_cell.contains(X, self):
simplex = cur_face simplex = cur_cell
continue continue
for neighbor in cur_face.neighbors: for neighbor in cur_cell.neighbors:
if (neighbor not in checked_faces) and (neighbor not in faces_to_check): if (neighbor not in checked_cells) and (neighbor not in cells_to_check):
faces_to_check.append(neighbor) cells_to_check.append(neighbor)
if not simplex: if not simplex:
raise Exception('no containing simplex found') raise Exception('no containing simplex found')
@ -262,16 +262,16 @@ class delaunay_grid(grid):
R is a grid object that is the (a) containing simplex around point X R is a grid object that is the (a) containing simplex around point X
S is a connectivity-based nearest-neighbor lookup, limited to 3 extra verts S is a connectivity-based nearest-neighbor lookup, limited to 3 extra verts
""" """
if not self.faces: if not self.cells:
self.construct_connectivity() self.construct_connectivity()
# get closest point # get closest point
(dist, indicies) = self.tree.query(X, 2) (dist, indicies) = self.tree.query(X, 2)
simplex = None simplex = None
for face in self.faces_for_vert[indicies[0]]: for cell in self.cells_for_vert[indicies[0]]:
if face.contains(X, self): if cell.contains(X, self):
simplex = face simplex = cell
break break
if not simplex: if not simplex:
@ -312,27 +312,27 @@ class delaunay_grid(grid):
""" """
log.debug('start') log.debug('start')
qdelaunay_string = get_qdelaunay_dump_str(self) qdelaunay_string = get_qdelaunay_dump_str(self)
face_to_faces = [] cell_to_cells = []
for matcher in grid.face_re.finditer(qdelaunay_string): for matcher in grid.cell_re.finditer(qdelaunay_string):
d = matcher.groupdict() d = matcher.groupdict()
face_name = d['face'] cell_name = d['cell']
verticies = d['verts'] verticies = d['verts']
neighboring_faces = d['neigh'] neighboring_cells = d['neigh']
cur_face = face(face_name) cur_cell = cell(cell_name)
self.faces[face_name] = cur_face self.cells[cell_name] = cur_cell
for v in grid.vert_re.findall(verticies): for v in grid.vert_re.findall(verticies):
vertex_index = int(v[1:]) vertex_index = int(v[1:])
cur_face.add_vert(vertex_index) cur_cell.add_vert(vertex_index)
self.faces_for_vert[vertex_index].append(cur_face) self.cells_for_vert[vertex_index].append(cur_cell)
nghbrs = [(face_name, i) for i in neighboring_faces.split()] nghbrs = [(cell_name, i) for i in neighboring_cells.split()]
face_to_faces.extend(nghbrs) cell_to_cells.extend(nghbrs)
for rel in face_to_faces: for rel in cell_to_cells:
if rel[1] in self.faces: if rel[1] in self.cells:
self.faces[rel[0]].add_neighbor(self.faces[rel[1]]) self.cells[rel[0]].add_neighbor(self.cells[rel[1]])
log.debug('end') log.debug('end')

View File

@ -7,7 +7,7 @@ import numpy as np
from scipy.spatial import KDTree from scipy.spatial import KDTree
from interp.grid import grid from interp.grid import grid
from interp.grid.simplex import face from interp.grid.simplex import cell
from interp.tools import exact_func, exact_func_3D from interp.tools import exact_func, exact_func_3D
@ -58,45 +58,45 @@ class gmsh_grid(grid):
gmsh_file.readline() # $EndNodes gmsh_file.readline() # $EndNodes
gmsh_file.readline() # $Elements gmsh_file.readline() # $Elements
# temporary dict used to compute face connectivity # temporary dict used to compute cell connectivity
neighbors = {} neighbors = {}
element_count = int(gmsh_file.readline()) element_count = int(gmsh_file.readline())
for i in xrange(element_count): for i in xrange(element_count):
cur_line = gmsh_file.readline() cur_line = gmsh_file.readline()
cur_line = cur_line.split() cur_line = cur_line.split()
cur_face_index, node_type, rest = (int(cur_line[0]), cur_cell_index, node_type, rest = (int(cur_line[0]),
int(cur_line[1]), int(cur_line[1]),
[int(j) for j in cur_line[2:]]) [int(j) for j in cur_line[2:]])
if(node_type == THREE_NODE_TRIANGLE): if(node_type == THREE_NODE_TRIANGLE):
points_for_cur_face = [i-1 for i in rest[rest[0]+1:]] points_for_cur_cell = [i-1 for i in rest[rest[0]+1:]]
cur_face = face(cur_face_index) cur_cell = cell(cur_cell_index)
for cur_point in points_for_cur_face: for cur_point in points_for_cur_cell:
self.faces_for_vert[cur_point].append(cur_face) self.cells_for_vert[cur_point].append(cur_cell)
cur_face.verts = points_for_cur_face cur_cell.verts = points_for_cur_cell
self.faces[cur_face_index] = cur_face self.cells[cur_cell_index] = cur_cell
edges = [tuple(sorted(i)) for i in combinations(points_for_cur_face, EDGES_FOR_FACE_CONNECTIVITY)] edges = [tuple(sorted(i)) for i in combinations(points_for_cur_cell, EDGES_FOR_FACE_CONNECTIVITY)]
# edge is two verts # edge is two verts
for edge in edges: for edge in edges:
if edge in neighbors: if edge in neighbors:
neighbors[edge].append(cur_face_index) neighbors[edge].append(cur_cell_index)
else: else:
neighbors[edge] = [cur_face_index] neighbors[edge] = [cur_cell_index]
for k,v in neighbors.iteritems(): for k,v in neighbors.iteritems():
if len(v) > 1: if len(v) > 1:
self.faces[v[0]].add_neighbor(self.faces[v[1]]) self.cells[v[0]].add_neighbor(self.cells[v[1]])
self.faces[v[1]].add_neighbor(self.faces[v[0]]) self.cells[v[1]].add_neighbor(self.cells[v[0]])
def dump_to_blender_files(self, pfile = '/tmp/points.p', ffile = '/tmp/faces.p'): def dump_to_blender_files(self, pfile = '/tmp/points.p', ffile = '/tmp/cells.p'):
pickle.dump([(p[0], p[1], p[2]) for p in self.verts], open(pfile, 'w')) pickle.dump([(p[0], p[1], p[2]) for p in self.verts], open(pfile, 'w'))
pickle.dump([f.verts for f in self.faces.itervalues()], open(ffile, 'w')) pickle.dump([f.verts for f in self.cells.itervalues()], open(ffile, 'w'))
class gmsh_grid3D(grid): class gmsh_grid3D(grid):
@ -136,41 +136,41 @@ class gmsh_grid3D(grid):
gmsh_file.readline() # $EndNodes gmsh_file.readline() # $EndNodes
gmsh_file.readline() # $Elements gmsh_file.readline() # $Elements
# temporary dict used to compute face connectivity # temporary dict used to compute cell connectivity
neighbors = {} neighbors = {}
element_count = int(gmsh_file.readline()) element_count = int(gmsh_file.readline())
for i in xrange(element_count): for i in xrange(element_count):
cur_line = gmsh_file.readline() cur_line = gmsh_file.readline()
cur_line = cur_line.split() cur_line = cur_line.split()
cur_face_index, node_type, rest = (int(cur_line[0]), cur_cell_index, node_type, rest = (int(cur_line[0]),
int(cur_line[1]), int(cur_line[1]),
[int(j) for j in cur_line[2:]]) [int(j) for j in cur_line[2:]])
if(node_type == FOUR_NODE_TET): if(node_type == FOUR_NODE_TET):
points_for_cur_face = [i-1 for i in rest[rest[0]+1:]] points_for_cur_cell = [i-1 for i in rest[rest[0]+1:]]
cur_face = face(cur_face_index) cur_cell = cell(cur_cell_index)
for cur_point in points_for_cur_face: for cur_point in points_for_cur_cell:
self.faces_for_vert[cur_point].append(cur_face) self.cells_for_vert[cur_point].append(cur_cell)
cur_face.verts = points_for_cur_face cur_cell.verts = points_for_cur_cell
self.faces[cur_face_index] = cur_face self.cells[cur_cell_index] = cur_cell
edges = [tuple(sorted(i)) for i in combinations(points_for_cur_face, EDGES_FOR_VOLUME_CONNECTIVITY)] edges = [tuple(sorted(i)) for i in combinations(points_for_cur_cell, EDGES_FOR_VOLUME_CONNECTIVITY)]
for edge in edges: for edge in edges:
if edge in neighbors: if edge in neighbors:
neighbors[edge].append(cur_face_index) neighbors[edge].append(cur_cell_index)
else: else:
neighbors[edge] = [cur_face_index] neighbors[edge] = [cur_cell_index]
for k,v in neighbors.iteritems(): for k,v in neighbors.iteritems():
if len(v) > 1: if len(v) > 1:
self.faces[v[0]].add_neighbor(self.faces[v[1]]) self.cells[v[0]].add_neighbor(self.cells[v[1]])
self.faces[v[1]].add_neighbor(self.faces[v[0]]) self.cells[v[1]].add_neighbor(self.cells[v[0]])
def dump_to_blender_files(self, pfile = '/tmp/points.p', ffile = '/tmp/faces.p'): def dump_to_blender_files(self, pfile = '/tmp/points.p', ffile = '/tmp/cells.p'):
pickle.dump([(p[0], p[1], p[2]) for p in self.verts], open(pfile, 'w')) pickle.dump([(p[0], p[1], p[2]) for p in self.verts], open(pfile, 'w'))
pickle.dump([f.verts for f in self.faces.itervalues()], open(ffile, 'w')) pickle.dump([f.verts for f in self.cells.itervalues()], open(ffile, 'w'))

View File

@ -21,7 +21,7 @@ def contains(X, R):
return r return r
class face(object): class cell(object):
def __init__(self, name): def __init__(self, name):
self.name = name self.name = name
self.verts = [] self.verts = []
@ -35,7 +35,7 @@ class face(object):
def add_neighbor(self, n): def add_neighbor(self, n):
""" """
reference to another face object reference to another cell object
""" """
self.neighbors.append(n) self.neighbors.append(n)
@ -44,7 +44,7 @@ class face(object):
X = point of interest X = point of interest
G = corrensponding grid object (G.verts) G = corrensponding grid object (G.verts)
because of the way i'm storing things, because of the way i'm storing things,
a face simply stores indicies, and so one a cell simply stores indicies, and so one
must pass in a reference to the grid object must pass in a reference to the grid object
containing real verts. containing real verts.
@ -54,7 +54,7 @@ class face(object):
def __str__(self): def __str__(self):
neighbors = [str(i.name) for i in self.neighbors] neighbors = [str(i.name) for i in self.neighbors]
return '<face %s: verts: %s neighbors: [%s]>' %\ return '<cell %s: verts: %s neighbors: [%s]>' %\
( (
self.name, self.name,
self.verts, self.verts,