Changeset 8048

Timestamp:

Oct 21, 2010, 3:37:57 PM (14 years ago)

Author:

habili

Message:

Changed from using polygons to polylines. The inlets are now represented by triangles that intersect or contain a polyline.

Location:

trunk/anuga_core/source/anuga/structures

Files:

• inlet.py (modified) (5 diffs)
• structure_operator.py (modified) (5 diffs)

trunk/anuga_core/source/anuga/structures/inlet.py

@@ -1 @@
-from anuga.geometry.polygon import inside_polygon, is_inside_polygon, polygon_overlap
+from anuga.geometry.polygon import inside_polygon, is_inside_polygon, polygon_overlap, polyline_overlap
 from anuga.config import velocity_protection, g
 import math
@@ -9 +9 @@
     """
 
-    def __init__(self, domain, polygon, enquiry_pt,  outward_culvert_vector=None, verbose=False):
+    def __init__(self, domain, polyline, enquiry_pt,  outward_culvert_vector=None, verbose=False):
 
         self.domain = domain
         self.domain_bounding_polygon = self.domain.get_boundary_polygon()
-        self.polygon = polygon
+        self.polyline = polyline
         self.enquiry_pt = enquiry_pt
         self.outward_culvert_vector = outward_culvert_vector
@@ -29 +29 @@
         vertex_coordinates = self.domain.get_vertex_coordinates(absolute=True)
 
-        # Check that polygon lies within the mesh.
-        for point in self.polygon:  
+        # Check that polyline lies within the mesh.
+        for point in self.polyline:  
                 msg = 'Point %s ' %  str(point)
                 msg += ' did not fall within the domain boundary.'
@@ -40 +40 @@
         assert is_inside_polygon(point, bounding_polygon), msg
 
-        self.triangle_indices = polygon_overlap(vertex_coordinates, self.polygon)
-        #self.triangle_indices_cen = inside_polygon(domain_centroids, self.polygon, verbose=self.verbose)
+        self.triangle_indices = polyline_overlap(vertex_coordinates, self.polyline)
 
         if len(self.triangle_indices) == 0:
-            region = 'Inlet polygon=%s' % (self.polygon)
+            region = 'Inlet polyline=%s' % (self.polyline)
             msg = 'No triangles have been identified in region '
             raise Exception, msg
@@ -54 +53 @@
         
         # Compute inlet area as the sum of areas of triangles identified
-        # by polygon. Must be called after compute_inlet_triangle_indices().
+        # by polyline. Must be called after compute_inlet_triangle_indices().
         if len(self.triangle_indices) == 0:
-            region = 'Inlet polygon=%s' % (self.inlet_polygon)
+            region = 'Inlet polyline=%s' % (self.inlet_polyline)
             msg = 'No triangles have been identified in region '
             raise Exception, msg

trunk/anuga_core/source/anuga/structures/structure_operator.py

@@ -86 +86 @@
 
         self.inlets = []
-        polygon0 = self.inlet_polygons[0]
+        polyline0 = self.inlet_polylines[0]
         enquiry_point0 = self.inlet_equiry_points[0]
         outward_vector0 = self.culvert_vector
-        self.inlets.append(inlet.Inlet(self.domain, polygon0, enquiry_point0, outward_vector0))
-
-        polygon1 = self.inlet_polygons[1]
-        exchange_polygon1 = self.inlet_equiry_points[1]
+        self.inlets.append(inlet.Inlet(self.domain, polyline0, enquiry_point0, outward_vector0))
+
+        polyline1 = self.inlet_polylines[1]
+        enquiry_point1 = self.inlet_equiry_points[1]
         outward_vector1  = - self.culvert_vector
-        self.inlets.append(inlet.Inlet(self.domain, polygon1, exchange_polygon1, outward_vector1))
+        self.inlets.append(inlet.Inlet(self.domain, polyline1, enquiry_point1, outward_vector1))
 
         self.set_logging(logging)
@@ -196 +196 @@
     def __create_exchange_polygons(self):
 
-        """Create polygons at the end of a culvert inlet and outlet.
-        At either end two polygons will be created; one for the actual flow to pass through and one a little further away
+        """Create polylines at the end of a culvert inlet and outlet.
+        At either end two polylines will be created; one for the actual flow to pass through and one a little further away
         for enquiring the total energy at both ends of the culvert and transferring flow.
         """
@@ -218 +218 @@
         # Short hands
         w = 0.5*self.width*self.culvert_normal # Perpendicular vector of 1/2 width
-        h = self.apron*self.culvert_vector    # Vector of length=height in the
+        #h = self.apron*self.culvert_vector    # Vector of length=height in the
                              # direction of the culvert
 
-        gap = 1.5*h + self.enquiry_gap
-
-        self.inlet_polygons = []
+        #gap = 1.5*h + self.enquiry_gap
+        gap = 1.5*self.culvert_vector + self.enquiry_gap
+
+        self.inlet_polylines = []
         self.inlet_equiry_points = []
 
@@ -231 +232 @@
             p0 = self.end_points[i] + w
             p1 = self.end_points[i] - w
-            p2 = p1 + i0*h
-            p3 = p0 + i0*h
             ep = self.end_points[i] + i0*gap
 
-            self.inlet_polygons.append(num.array([p0, p1, p2, p3]))
-            self.inlet_equiry_points.append(ep)
-
-        # Check that enquiry points are outside inlet polygons
-        for i in [0,1]:
-            polygon = self.inlet_polygons[i]
-            ep = self.inlet_equiry_points[i]
-           
-            area = anuga.polygon_area(polygon)
-            
-            msg = 'Polygon %s ' %(polygon)
-            msg += ' has area = %f' % area
-            assert area > 0.0, msg
-
-            msg = 'Enquiry point falls inside an exchange polygon.'
-            assert not anuga.inside_polygon(ep, polygon), msg
-            
+            self.inlet_polylines.append(num.array([p0, p1]))
+            self.inlet_equiry_points.append(ep)            
     
     def discharge_routine(self):
@@ -284 +268 @@
 
             message += 'polygon\n'
-            message += '%s' % inlet.polygon
+            message += '%s' % inlet.polyline
             message += '\n'