Changeset 8050 for trunk/anuga_core/source/anuga

Timestamp:

Oct 22, 2010, 5:10:19 PM (13 years ago)

Author:

steve

Message:

Added in an inlet_operator class to implement inflow "boundary condition" by setting up a inlet defined by a line close to the boundary. Look at testing_culvert_inlet.py as an example.

Location:

trunk/anuga_core/source/anuga/structures

Files:

• inlet.py (modified) (4 diffs)
• inlet_enquiry.py (added)
• inlet_operator.py (added)
• structure_operator.py (modified) (2 diffs)
• testing_culvert.py (modified) (3 diffs)
• testing_culvert_inlet.py (added)

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

@@ -9 +9 @@
     """
 
-    def __init__(self, domain, polyline, enquiry_pt,  outward_culvert_vector=None, verbose=False):
+    def __init__(self, domain, polyline, verbose=False):
 
         self.domain = domain
         self.domain_bounding_polygon = self.domain.get_boundary_polygon()
         self.polyline = polyline
-        self.enquiry_pt = enquiry_pt
-        self.outward_culvert_vector = outward_culvert_vector
         self.verbose = verbose
 
-        self.compute_indices()
+        self.compute_triangle_indices()
         self.compute_area()
 
 
-    def compute_indices(self):
+    def compute_triangle_indices(self):
 
         # Get boundary (in absolute coordinates)
@@ -35 +33 @@
                 assert is_inside_polygon(point, bounding_polygon), msg
                 
-        point = self.enquiry_pt
-        msg = 'Enquiry Point %s ' %  str(point)
-        msg += ' did not fall within the domain boundary.'
-        assert is_inside_polygon(point, bounding_polygon), msg
+
 
         self.triangle_indices = polyline_overlap(vertex_coordinates, self.polyline)
@@ -46 +41 @@
             msg += 'No triangles intersecting polyline '
             raise Exception, msg
-            
-        self.enquiry_index = self.domain.get_triangle_containing_point(self.enquiry_pt)
 
-        if self.enquiry_index in self.triangle_indices:
-            msg = 'Enquiry point %s' % (self.enquiry_pt)
-            msg += 'is in an inlet triangle'
-            raise Exception, msg
+
 
     def compute_area(self):
@@ -180 +170 @@
 
 
-    def get_enquiry_stage(self):
-
-        return self.domain.quantities['stage'].centroid_values[self.enquiry_index]
-
-
-    def get_enquiry_xmom(self):
-
-        return self.domain.quantities['xmomentum'].centroid_values[self.enquiry_index]
-
-    def get_enquiry_ymom(self):
-
-        return self.domain.quantities['ymomentum'].centroid_values[self.enquiry_index]
-
-
-    def get_enquiry_elevation(self):
-
-        return self.domain.quantities['elevation'].centroid_values[self.enquiry_index]
-
-    def get_enquiry_height(self):
-
-        return self.get_enquiry_stage() - self.get_enquiry_elevation()
-
-
-    def get_enquiry_velocity(self):
-
-            height = self.get_enquiry_height()
-            u = self.get_enquiry_xmom()/(height + velocity_protection/height)
-            v = self.get_enquiry_ymom()/(height + velocity_protection/height)
-
-            return u, v
-
-
-    def get_enquiry_xvelocity(self):
-
-            height = self.get_enquiry_height()
-            return self.get_enquiry_xmom()/(height + velocity_protection/height)
-
-    def get_enquiry_yvelocity(self):
-
-            height = self.get_enquiry_height()
-            return self.get_enquiry_ymom()/(height + velocity_protection/height)
-
-
-    def get_enquiry_speed(self):
-
-            u, v = self.get_enquiry_velocity()
-
-            return math.sqrt(u**2 + v**2)
-
-
-    def get_enquiry_velocity_head(self):
-
-        return 0.5*self.get_enquiry_speed()**2/g
-
-
-    def get_enquiry_total_energy(self):
-
-        return self.get_enquiry_velocity_head() + self.get_enquiry_stage()
-
-
-    def get_enquiry_specific_energy(self):
-
-        return self.get_enquiry_velocity_head() + self.get_enquiry_height()
-
 
     def set_heights(self,height):

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

@@ -2 +2 @@
 import numpy as num
 import math
-import inlet
+import inlet_enquiry
 
 from anuga.utilities.system_tools import log_to_file
@@ -89 +89 @@
         enquiry_point0 = self.inlet_equiry_points[0]
         outward_vector0 = self.culvert_vector
-        self.inlets.append(inlet.Inlet(self.domain, polyline0, enquiry_point0, outward_vector0))
+        self.inlets.append(inlet_enquiry.Inlet_enquiry(self.domain, polyline0,
+                           enquiry_point0, outward_vector0, self.verbose))
 
         polyline1 = self.inlet_polylines[1]
         enquiry_point1 = self.inlet_equiry_points[1]
         outward_vector1  = - self.culvert_vector
-        self.inlets.append(inlet.Inlet(self.domain, polyline1, enquiry_point1, outward_vector1))
+        self.inlets.append(inlet_enquiry.Inlet_enquiry(self.domain, polyline1,
+                           enquiry_point1, outward_vector1, self.verbose))
 
         self.set_logging(logging)

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

@@ -11 +11 @@
 
 from anuga.structures.boyd_pipe_operator import Boyd_pipe_operator
+from anuga.structures.inlet_operator import Inlet_operator
                             
 #from anuga.culvert_flows.culvert_routines import boyd_generalised_culvert_model
@@ -101 +102 @@
 
 
-#culvert2 = Generic_box_culvert(domain,
-#                              end_point0=[19.0, 2.5],
-#                              end_point1=[25.0, 2.5],
-#                              width=1.00,
-#                              verbose=False)
+line = [[0.0, 5.0], [0.0, 10.0]]
+Q = 5.0
+#Inlet_operator(domain, line, Q)
 
 
 
-
-#print domain.fractional_step_operators
-
-#domain.apply_fractional_steps()
 
 ##-----------------------------------------------------------------------
@@ -143 +138 @@
     domain.write_time()
 
-    if domain.get_time() > 150.5 and domain.get_time() < 151.5 :
-        Bi = anuga.Dirichlet_boundary([0.0, 0.0, 0.0])
-        domain.set_boundary({'left': Bi, 'right': Br, 'top': Br, 'bottom': Br})
+    #if domain.get_time() > 150.5 and domain.get_time() < 151.5 :
+        #Bi = anuga.Dirichlet_boundary([0.0, 0.0, 0.0])
+        #domain.set_boundary({'left': Bi, 'right': Br, 'top': Br, 'bottom': Br})
 
     #delta_w = culvert.inlet.stage - culvert.outlet.stage