Changeset 7713
- Timestamp:
- May 11, 2010, 10:59:06 AM (15 years ago)
- Location:
- anuga_core/source/anuga
- Files:
-
- 6 edited
-
fit_interpolate/general_fit_interpolate.py (modified) (1 diff)
-
fit_interpolate/search_functions.py (modified) (8 diffs)
-
fit_interpolate/test_meshquad.py (modified) (2 diffs)
-
fit_interpolate/test_search_functions.py (modified) (3 diffs)
-
geometry/quad.py (modified) (6 diffs)
-
geometry/test_geometry.py (modified) (3 diffs)
Legend:
- Unmodified
- Added
- Removed
-
anuga_core/source/anuga/fit_interpolate/general_fit_interpolate.py
r7711 r7713 36 36 import numpy as num 37 37 38 39 # tests fail if 2 is used40 MAX_VERTICES_PER_CELL = 13 # A value of 8 or lower can cause problems,41 # if a vert has 9 triangles.42 38 43 39 build_quadtree_time = 0 -
anuga_core/source/anuga/fit_interpolate/search_functions.py
r7711 r7713 33 33 # FIXME(Ole): Could we come up with a less confusing structure? 34 34 # FIXME(James): remove this global var 35 LAST_TRIANGLE = [[-10, 35 LAST_TRIANGLE = [[-10, -10, 36 36 (num.array([[max_float, max_float], 37 37 [max_float, max_float], … … 64 64 global search_more_cells_time 65 65 66 # Search the last triangle first 66 if last_triangle[0][1] != -10: 67 # check the last triangle found first 68 element_found, sigma0, sigma1, sigma2, k = \ 69 _search_triangles_of_vertices(last_triangle, x) 70 71 if element_found: 72 return True, sigma0, sigma1, sigma2, k 73 74 # search to bottom of tree from last found leaf 75 branch = last_triangle[0][1] 76 77 if branch == -10: 78 branch = root 79 80 tri_data = branch.search(x[0], x[1]) 81 triangles = _trilist_from_data(mesh, tri_data) 67 82 element_found, sigma0, sigma1, sigma2, k = \ 68 _search_triangles_of_vertices(last_triangle, x) 69 70 if element_found is True: 71 return element_found, sigma0, sigma1, sigma2, k 72 73 74 # Get triangles in the cell that the point is in. 75 tri_indices = root.search(x[0], x[1]) 76 triangles = _trilist_from_indices(mesh, tri_indices) 77 78 element_found, sigma0, sigma1, sigma2, k = \ 79 _search_triangles_of_vertices(triangles, x) 80 81 is_more_elements = True 82 83 # while not element_found and is_more_elements: 84 # triangles, branch = root.expand_search() 85 # if branch == []: 86 # # Searching all the verts from the root cell that haven't 87 # # been searched. This is the last try 88 # element_found, sigma0, sigma1, sigma2, k = \ 89 # _search_triangles_of_vertices(triangles, x) 90 # is_more_elements = False 91 # else: 92 # element_found, sigma0, sigma1, sigma2, k = \ 93 # _search_triangles_of_vertices(triangles, x) 94 95 83 _search_triangles_of_vertices(triangles, x) 84 if element_found: 85 return True, sigma0, sigma1, sigma2, k 86 87 # search rest of tree 88 element_found = False 89 while branch: 90 if not branch.parent: 91 # search from top of tree if we are at root 92 siblings = [root] 93 else: 94 siblings = branch.get_siblings() 95 96 for sibling in siblings: 97 tri_data = sibling.search(x[0], x[1]) 98 triangles = _trilist_from_data(mesh, tri_data) 99 element_found, sigma0, sigma1, sigma2, k = \ 100 _search_triangles_of_vertices(triangles, x) 101 if element_found: 102 return True, sigma0, sigma1, sigma2, k 103 104 branch = branch.parent 105 if branch: 106 tri_data = branch.test_leaves(x[0], x[1]) 107 triangles = _trilist_from_data(mesh, tri_data) 108 element_found, sigma0, sigma1, sigma2, k = \ 109 _search_triangles_of_vertices(triangles, x) 110 if element_found: 111 return True, sigma0, sigma1, sigma2, k 112 96 113 return element_found, sigma0, sigma1, sigma2, k 114 115 116 ## Get triangles in the cell that the point is in. 117 #tri_data = root.search(x[0], x[1]) 118 #triangles = _trilist_from_data(mesh, tri_data) 119 120 #element_found, sigma0, sigma1, sigma2, k = \ 121 #_search_triangles_of_vertices(triangles, x) 122 123 #if element_found: 124 #return True, sigma0, sigma1, sigma2, k 97 125 98 126 … … 116 144 sigma1 = -10.0 117 145 k = -10 146 118 147 # For all vertices in same cell as point x 119 148 element_found = False 120 for k, tri_verts_norms in triangles:149 for k, node, tri_verts_norms in triangles: 121 150 tri = tri_verts_norms[0] 122 151 tri = ensure_numeric(tri) … … 135 164 136 165 # Don't look for any other triangles in the triangle list 137 last_triangle = [[k, tri_verts_norms]]166 last_triangle = [[k, node, tri_verts_norms]] 138 167 break 139 168 … … 141 170 142 171 143 def _trilist_from_ indices(mesh, indices):172 def _trilist_from_data(mesh, indices): 144 173 """return a list of lists. For the inner lists, 145 174 The first element is the triangle index, … … 151 180 152 181 ret_list = [] 153 for i in indices:182 for i, node in indices: 154 183 vertices = mesh.get_vertex_coordinates(triangle_id=i, absolute=True) 155 184 n0 = mesh.get_normal(i, 0) 156 185 n1 = mesh.get_normal(i, 1) 157 186 n2 = mesh.get_normal(i, 2) 158 ret_list.append([i, [vertices, (n0, n1, n2)]])187 ret_list.append([i, node, [vertices, (n0, n1, n2)]]) 159 188 return ret_list 160 189 … … 196 225 class MeshQuadtree(Cell): 197 226 """ A quadtree constructed from the given mesh. 227 This class is actually the root node of a tree, 228 and derives from a Cell. 198 229 """ 199 230 def __init__(self, mesh): … … 205 236 extents = AABB(*mesh.get_extent(absolute=True)) 206 237 extents.grow(1.001) # To avoid round off error 207 Cell.__init__(self, extents )238 Cell.__init__(self, extents, None) # root has no parent 208 239 209 240 N = len(mesh) -
anuga_core/source/anuga/fit_interpolate/test_meshquad.py
r7711 r7713 43 43 44 44 # test a point that falls within a triangle 45 result = Q.search(10, 10 , get_vertices=True)45 result = Q.search(10, 10) 46 46 assert type(result) in [types.ListType,types.TupleType],\ 47 47 'should be a list' 48 pos, index= result[0]48 index, parent = result[0] 49 49 self.assertEqual(index, 1) 50 50 … … 138 138 results = Q.search(4.5, 3) 139 139 assert len(results) == 1 140 self.assertEqual(results[0], 2) 140 idx, _ = results[0] 141 self.assertEqual(idx, 2) 141 142 results = Q.search(5,4.) 142 143 self.assertEqual(len(results),1) 143 self.assertEqual(results[0], 2) 144 idx, _ = results[0] 145 self.assertEqual(idx, 2) 144 146 145 147 def NOtest_get_parent(): -
anuga_core/source/anuga/fit_interpolate/test_search_functions.py
r7711 r7713 5 5 from search_functions import search_tree_of_vertices, set_last_triangle 6 6 from search_functions import _search_triangles_of_vertices 7 from search_functions import _trilist_from_ indices7 from search_functions import _trilist_from_data 8 8 from search_functions import compute_interpolation_values, MeshQuadtree 9 9 … … 156 156 x = ensure_numeric([0.5, 0.5]) 157 157 158 triangles = _trilist_from_ indices(mesh, root.search(x[0], x[1]))158 triangles = _trilist_from_data(mesh, root.search(x[0], x[1])) 159 159 160 160 found, sigma0, sigma1, sigma2, k = \ … … 175 175 [10, 3]]: 176 176 177 triangles = _trilist_from_ indices(mesh, root.search(x[0], x[1]))177 triangles = _trilist_from_data(mesh, root.search(x[0], x[1])) 178 178 179 179 #print x, candidate_vertices -
anuga_core/source/anuga/geometry/quad.py
r7712 r7713 19 19 """ 20 20 21 def __init__(self, extents, 21 def __init__(self, extents, parent, depth = 0, 22 22 name = 'cell'): 23 23 … … 26 26 27 27 self.extents = extents 28 self.parent = parent 29 self.searched = [-10, -15] 30 self.depth = depth 28 31 29 32 # The points in this cell … … 33 36 34 37 def __repr__(self): 35 str = ' %s: leaves: %d' \36 % (self. name , len(self.leaves))38 str = '(%d)%s: leaves: %d' \ 39 % (self.depth, self.name , len(self.leaves)) 37 40 if self.children: 38 41 str += ', children: %d' % (len(self.children)) … … 82 85 subregion1, subregion2 = self.extents.split() 83 86 if subregion1.is_trivial_in(new_region): 84 self.children = [Cell(subregion1 ), Cell(subregion2)]87 self.children = [Cell(subregion1, self, depth=self.depth+1), Cell(subregion2, self, depth=self.depth+1)] 85 88 self.children[0]._insert(new_leaf) 86 89 return 87 90 elif subregion2.is_trivial_in(new_region): 88 self.children = [Cell(subregion1 ), Cell(subregion2)]91 self.children = [Cell(subregion1, self, depth=self.depth+1), Cell(subregion2, self, depth=self.depth+1)] 89 92 self.children[1]._insert(new_leaf) 90 93 return … … 132 135 133 136 134 def search(self, x, y , get_vertices = False):137 def search(self, x, y): 135 138 """return a list of possible intersections with geometry""" 139 140 # if self.searched == [x,y]: 141 # print 'ERROR: already searched at ', [x,y] 142 # self.searched = [x,y] 136 143 144 intersecting_regions = self.test_leaves(x, y) 145 146 # recurse down into nodes that the point passes through 147 if self.children: 148 for child in self.children: 149 if child.extents.contains(x, y): 150 intersecting_regions.extend(child.search(x, y)) 151 152 return intersecting_regions 153 154 155 def test_leaves(self, x, y): 137 156 intersecting_regions = [] 138 157 … … 141 160 aabb, data = leaf 142 161 if aabb.contains(x, y): 143 if get_vertices: 144 intersecting_regions.append(leaf) 145 else: 146 intersecting_regions.append(data) 147 148 # recurse down into nodes that the point passes through 149 if self.children: 150 for child in self.children: 151 if child.extents.contains(x, y): 152 intersecting_regions.extend(child.search(x, y, get_vertices)) 153 154 return intersecting_regions 155 162 intersecting_regions.append((data, self)) 163 164 return intersecting_regions 165 166 167 def get_siblings(self): 168 """ return parent and siblings of this node. 169 """ 170 if not self.parent: 171 return [] 172 173 siblings = list(self.parent.children) 174 siblings.remove(self) 175 return siblings 176 -
anuga_core/source/anuga/geometry/test_geometry.py
r7711 r7713 59 59 60 60 def test_add_data(self): 61 cell = Cell(AABB(0,10, 0,5) )61 cell = Cell(AABB(0,10, 0,5), None) 62 62 cell.insert([(AABB(1,3, 1, 3), 111), (AABB(8,9, 1, 2), 222), \ 63 63 (AABB(7, 8, 3, 4), 333), (AABB(1, 10, 0, 1), 444)]) … … 70 70 71 71 def test_search(self): 72 test_ region = (AABB(8,9, 1, 2), 222)73 cell = Cell(AABB(0,10, 0,5) )74 cell.insert([(AABB(1,3, 1, 3), 111), test_region, \72 test_tag = 222 73 cell = Cell(AABB(0,10, 0,5), None) 74 cell.insert([(AABB(1,3, 1, 3), 111), (AABB(8,9, 1, 2), test_tag), \ 75 75 (AABB(7, 8, 3, 4), 333), (AABB(1, 10, 0, 1), 444)]) 76 76 77 result = cell.search(x = 8.5, y = 1.5, get_vertices=True)77 result = cell.search(x = 8.5, y = 1.5) 78 78 assert type(result) in [types.ListType,types.TupleType],\ 79 79 'should be a list' 80 self.assertEqual(result, [test_region], 'only 1 point should intersect') 80 assert(len(result) == 1) 81 data, node = result[0] 82 self.assertEqual(data, test_tag, 'only 1 point should intersect') 81 83 84 def test_get_siblings(self): 85 """ Make sure children know their parent. """ 86 cell = Cell(AABB(0,10, 0,5), None) 87 cell.insert([(AABB(1,3, 1, 3), 111), (AABB(8,9, 1, 2), 222)]) 88 assert len(cell.children) == 2 89 assert cell.parent == None 90 assert cell.children[0].parent == cell 91 assert cell.children[1].parent == cell 82 92 83 93 def test_clear_1(self): 84 cell = Cell(AABB(0,10, 0,5) )94 cell = Cell(AABB(0,10, 0,5), None) 85 95 cell.insert([(AABB(1,3, 1, 3), 111), (AABB(8,9, 1, 2), 222), \ 86 96 (AABB(7, 8, 3, 4), 333), (AABB(1, 10, 0, 1), 444)]) … … 94 104 95 105 if __name__ == "__main__": 96 mysuite = unittest.makeSuite(Test_ Quad,'test')106 mysuite = unittest.makeSuite(Test_Geometry, 'test') 97 107 runner = unittest.TextTestRunner() 98 108 runner.run(mysuite)
Note: See TracChangeset
for help on using the changeset viewer.
