Changeset 8050


Ignore:
Timestamp:
Oct 22, 2010, 5:10:19 PM (8 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:
3 added
3 edited

Legend:

Unmodified
Added
Removed
  • trunk/anuga_core/source/anuga/structures/inlet.py

    r8049 r8050  
    99    """
    1010
    11     def __init__(self, domain, polyline, enquiry_pt,  outward_culvert_vector=None, verbose=False):
     11    def __init__(self, domain, polyline, verbose=False):
    1212
    1313        self.domain = domain
    1414        self.domain_bounding_polygon = self.domain.get_boundary_polygon()
    1515        self.polyline = polyline
    16         self.enquiry_pt = enquiry_pt
    17         self.outward_culvert_vector = outward_culvert_vector
    1816        self.verbose = verbose
    1917
    20         self.compute_indices()
     18        self.compute_triangle_indices()
    2119        self.compute_area()
    2220
    2321
    24     def compute_indices(self):
     22    def compute_triangle_indices(self):
    2523
    2624        # Get boundary (in absolute coordinates)
     
    3533                assert is_inside_polygon(point, bounding_polygon), msg
    3634               
    37         point = self.enquiry_pt
    38         msg = 'Enquiry Point %s ' %  str(point)
    39         msg += ' did not fall within the domain boundary.'
    40         assert is_inside_polygon(point, bounding_polygon), msg
     35
    4136
    4237        self.triangle_indices = polyline_overlap(vertex_coordinates, self.polyline)
     
    4641            msg += 'No triangles intersecting polyline '
    4742            raise Exception, msg
    48            
    49         self.enquiry_index = self.domain.get_triangle_containing_point(self.enquiry_pt)
    5043
    51         if self.enquiry_index in self.triangle_indices:
    52             msg = 'Enquiry point %s' % (self.enquiry_pt)
    53             msg += 'is in an inlet triangle'
    54             raise Exception, msg
     44
    5545
    5646    def compute_area(self):
     
    180170
    181171
    182     def get_enquiry_stage(self):
    183 
    184         return self.domain.quantities['stage'].centroid_values[self.enquiry_index]
    185 
    186 
    187     def get_enquiry_xmom(self):
    188 
    189         return self.domain.quantities['xmomentum'].centroid_values[self.enquiry_index]
    190 
    191     def get_enquiry_ymom(self):
    192 
    193         return self.domain.quantities['ymomentum'].centroid_values[self.enquiry_index]
    194 
    195 
    196     def get_enquiry_elevation(self):
    197 
    198         return self.domain.quantities['elevation'].centroid_values[self.enquiry_index]
    199 
    200     def get_enquiry_height(self):
    201 
    202         return self.get_enquiry_stage() - self.get_enquiry_elevation()
    203 
    204 
    205     def get_enquiry_velocity(self):
    206 
    207             height = self.get_enquiry_height()
    208             u = self.get_enquiry_xmom()/(height + velocity_protection/height)
    209             v = self.get_enquiry_ymom()/(height + velocity_protection/height)
    210 
    211             return u, v
    212 
    213 
    214     def get_enquiry_xvelocity(self):
    215 
    216             height = self.get_enquiry_height()
    217             return self.get_enquiry_xmom()/(height + velocity_protection/height)
    218 
    219     def get_enquiry_yvelocity(self):
    220 
    221             height = self.get_enquiry_height()
    222             return self.get_enquiry_ymom()/(height + velocity_protection/height)
    223 
    224 
    225     def get_enquiry_speed(self):
    226 
    227             u, v = self.get_enquiry_velocity()
    228 
    229             return math.sqrt(u**2 + v**2)
    230 
    231 
    232     def get_enquiry_velocity_head(self):
    233 
    234         return 0.5*self.get_enquiry_speed()**2/g
    235 
    236 
    237     def get_enquiry_total_energy(self):
    238 
    239         return self.get_enquiry_velocity_head() + self.get_enquiry_stage()
    240 
    241 
    242     def get_enquiry_specific_energy(self):
    243 
    244         return self.get_enquiry_velocity_head() + self.get_enquiry_height()
    245 
    246172
    247173    def set_heights(self,height):
  • trunk/anuga_core/source/anuga/structures/structure_operator.py

    r8049 r8050  
    22import numpy as num
    33import math
    4 import inlet
     4import inlet_enquiry
    55
    66from anuga.utilities.system_tools import log_to_file
     
    8989        enquiry_point0 = self.inlet_equiry_points[0]
    9090        outward_vector0 = self.culvert_vector
    91         self.inlets.append(inlet.Inlet(self.domain, polyline0, enquiry_point0, outward_vector0))
     91        self.inlets.append(inlet_enquiry.Inlet_enquiry(self.domain, polyline0,
     92                           enquiry_point0, outward_vector0, self.verbose))
    9293
    9394        polyline1 = self.inlet_polylines[1]
    9495        enquiry_point1 = self.inlet_equiry_points[1]
    9596        outward_vector1  = - self.culvert_vector
    96         self.inlets.append(inlet.Inlet(self.domain, polyline1, enquiry_point1, outward_vector1))
     97        self.inlets.append(inlet_enquiry.Inlet_enquiry(self.domain, polyline1,
     98                           enquiry_point1, outward_vector1, self.verbose))
    9799
    98100        self.set_logging(logging)
  • trunk/anuga_core/source/anuga/structures/testing_culvert.py

    r8007 r8050  
    1111
    1212from anuga.structures.boyd_pipe_operator import Boyd_pipe_operator
     13from anuga.structures.inlet_operator import Inlet_operator
    1314                           
    1415#from anuga.culvert_flows.culvert_routines import boyd_generalised_culvert_model
     
    101102
    102103
    103 #culvert2 = Generic_box_culvert(domain,
    104 #                              end_point0=[19.0, 2.5],
    105 #                              end_point1=[25.0, 2.5],
    106 #                              width=1.00,
    107 #                              verbose=False)
     104line = [[0.0, 5.0], [0.0, 10.0]]
     105Q = 5.0
     106#Inlet_operator(domain, line, Q)
    108107
    109108
    110109
    111 
    112 #print domain.fractional_step_operators
    113 
    114 #domain.apply_fractional_steps()
    115110
    116111##-----------------------------------------------------------------------
     
    143138    domain.write_time()
    144139
    145     if domain.get_time() > 150.5 and domain.get_time() < 151.5 :
    146         Bi = anuga.Dirichlet_boundary([0.0, 0.0, 0.0])
    147         domain.set_boundary({'left': Bi, 'right': Br, 'top': Br, 'bottom': Br})
     140    #if domain.get_time() > 150.5 and domain.get_time() < 151.5 :
     141        #Bi = anuga.Dirichlet_boundary([0.0, 0.0, 0.0])
     142        #domain.set_boundary({'left': Bi, 'right': Br, 'top': Br, 'bottom': Br})
    148143
    149144    #delta_w = culvert.inlet.stage - culvert.outlet.stage
Note: See TracChangeset for help on using the changeset viewer.