Changeset 9673


Ignore:
Timestamp:
Feb 13, 2015, 4:58:07 PM (9 years ago)
Author:
davies
Message:

New pumping station function + tests, and a slight update to a bridge validation report

Location:
trunk
Files:
4 edited

Legend:

Unmodified
Added
Removed
  • trunk/anuga_core/anuga/structures/internal_boundary_functions.py

    r9657 r9673  
    401401    """
    402402
    403     def __init__(self, pump_capacity, hw_to_start_pumping, verbose=True):
    404         """
    405             @param pump_capacity m^3/s
    406 
    407         """
     403    def __init__(self, domain, pump_capacity, hw_to_start_pumping, hw_to_stop_pumping,
     404                 initial_pump_rate=0., pump_rate_of_increase = 1.0e+100,
     405                 pump_rate_of_decrease = 1.0e+100, verbose=True):
     406        """
     407            @param domain ANUGA domain
     408            @param pump_capacity (m^3/s)
     409            @param hw_to_start_pumping Turn pumps on if hw exceeds this (m)
     410            @param hw_to_stop_pumping Turn pumps off if hw below this (m)
     411            @param initial_pump_rate rate of pumps at start of simulation  (m^3/s)
     412            @param pump_rate_of_increase Accelleration of pump rate when turning on (m^3/s/s)
     413            @param pump_rate_of_decrease Decelleration of pump rate when turning off (m^3/s/s)
     414            @param verbose
     415        """
     416
    408417        self.pump_capacity = pump_capacity
    409418        self.hw_to_start_pumping = hw_to_start_pumping
     419        self.hw_to_stop_pumping = hw_to_stop_pumping
     420
     421        self.pump_rate_of_increase = pump_rate_of_increase
     422        self.pump_rate_of_decrease = pump_rate_of_decrease
     423
     424        self.domain=domain
     425        self.last_time_called = domain.get_time()
     426        self.time = domain.get_time()
     427
     428        if hw_to_start_pumping < hw_to_stop_pumping:
     429            raise Exception('hw_to_start_pumping should be >= hw_to_stop_pumping')
     430
     431        if initial_pump_rate > pump_capacity:
     432            raise Exception('Initial pump rate is > pump capacity')
     433       
     434        if ((self.pump_rate_of_increase < 0.) | (self.pump_rate_of_decrease < 0.)):
     435            raise Exception('Pump rates of increase / decrease MUST be non-negative')
     436
     437        if ( (pump_capacity < 0.) | (initial_pump_rate < 0.) ):
     438            raise Exception('Pump rates cannot be negative')
     439
     440        self.pump_rate = initial_pump_rate
    410441
    411442        if verbose:
     
    419450    def __call__(self, hw_in, tw_in):
    420451        """
    421             @param hw_in The stage (or energy) at the headwater site
    422             @param tw_in The stage (or energy) at the tailwater site
    423         """
    424 
    425         if(hw_in > self.hw_to_start_pumping):
    426             Q = self.pump_capacity
     452            @param hw_in The stage (or energy) at the headwater site (m)
     453            @param tw_in The stage (or energy) at the tailwater site (m)
     454        """
     455
     456        # Compute the time since last called, so we can increase / decrease the pump rate if needed
     457        self.time = self.domain.get_time()
     458        if self.time > self.last_time_called:
     459            dt = self.time - self.last_time_called
     460            self.last_time_called = self.time
    427461        else:
    428             Q = 0.
    429 
    430         return Q
    431 
     462            dt = 0.
     463            if self.time != self.last_time_called:
     464                print self.time
     465                print self.last_time_called
     466                print self.time - self.last_time_called
     467                raise Exception('Impossible timestepping, ask Gareth')
     468         
     469        # Increase / decrease the pump rate if needed
     470        if hw_in < self.hw_to_stop_pumping:
     471            self.pump_rate = max(0., self.pump_rate - dt*self.pump_rate_of_decrease)
     472        elif hw_in > self.hw_to_start_pumping:
     473            self.pump_rate = min(self.pump_capacity, self.pump_rate + dt*self.pump_rate_of_increase)
     474
     475        return self.pump_rate
     476
  • trunk/anuga_core/anuga/structures/tests/test_internal_boundary_functions.py

    r9657 r9673  
    33
    44import numpy
     5import anuga
    56from anuga.structures import internal_boundary_functions
    67from anuga.structures.internal_boundary_functions import hecras_internal_boundary_function
     8from anuga.structures.internal_boundary_functions import pumping_station_function
    79from anuga.utilities.system_tools import get_pathname_from_package
    810
     11## Data for the domain used to test the pumping station function
     12boundaryPolygon = [ [0., 0.], [0., 100.], [100.0, 100.0], [100.0, 0.0]]
     13wallLoc = 50.
     14# The boundary polygon + riverwall breaks the mesh into multiple regions
     15# Must define the resolution in these areas with an xy point + maximum area
     16# Otherwise triangle.c gets confused
     17regionPtAreas = [ [99., 99., 20.0*20.0*0.5],
     18                  [1., 1., 20.0*20.0*0.5] ]
    919
    1020class Test_internal_boundary_functions(unittest.TestCase):
     
    1323
    1424        path = get_pathname_from_package('anuga.structures')
    15         self.input_file = filename1 = os.path.join(path, 'tests', 'data', 'hecras_bridge_table.csv')
    16         self.hb = hecras_internal_boundary_function(self.input_file, verbose=False)
     25        self.input_hecras_file = os.path.join(path, 'tests', 'data', 'hecras_bridge_table.csv')
    1726
    1827        return
     
    2332        return
    2433
     34    def create_domain(self, wallHeight, InitialOceanStage, InitialLandStage,
     35                      riverWall=None, riverWall_Par=None):
     36        # Riverwall = list of lists, each with a set of x,y,z (and optional QFactor) values
     37
     38        if(riverWall is None):
     39            riverWall = { 'centralWall':
     40                           [ [wallLoc, 0.0, wallHeight],
     41                             [wallLoc, 100.0, wallHeight]]
     42                        }
     43        if(riverWall_Par is None):
     44            riverWall_Par = {'centralWall':{'Qfactor':1.0}}
     45        # Make the domain
     46        anuga.create_mesh_from_regions(boundaryPolygon,
     47                                 boundary_tags={'left': [0],
     48                                                'top': [1],
     49                                                'right': [2],
     50                                                'bottom': [3]},
     51                                   maximum_triangle_area = 200.,
     52                                   minimum_triangle_angle = 28.0,
     53                                   filename = 'testRiverwall.msh',
     54                                   interior_regions =[ ], #[ [higherResPolygon, 1.*1.*0.5],
     55                                                          #  [midResPolygon, 3.0*3.0*0.5]],
     56                                   breaklines=riverWall.values(),
     57                                   use_cache=False,
     58                                   verbose=False,
     59                                   regionPtArea=regionPtAreas)
     60
     61        domain = anuga.create_domain_from_file('testRiverwall.msh')
     62
     63        # 05/05/2014 -- riverwalls only work with DE0 and DE1
     64        domain.set_flow_algorithm('DE0')
     65        domain.set_name('test_riverwall')
     66
     67        def topography(x,y):
     68            return -x/150.
     69
     70        def stagefun(x,y):
     71            stg = InitialOceanStage*(x>=wallLoc) + InitialLandStage*(x<wallLoc)
     72            return stg
     73
     74        # NOTE: Setting quantities at centroids is important for exactness of tests
     75        domain.set_quantity('elevation',topography,location='centroids')     
     76        domain.set_quantity('stage', stagefun,location='centroids')           
     77       
     78        domain.riverwallData.create_riverwalls(riverWall,riverWall_Par,verbose=False)
     79
     80        # Boundary conditions
     81        Br = anuga.Reflective_boundary(domain)
     82        domain.set_boundary({'left': Br, 'right': Br, 'top': Br, 'bottom':Br})
     83
     84        return domain
     85
    2586    def test_basic_properties(self):
    2687
    27         assert self.hb.name == self.input_file
     88        self.hb = hecras_internal_boundary_function(self.input_hecras_file, verbose=False)
     89        assert self.hb.name == self.input_hecras_file
    2890        assert len(self.hb.hw_max_given_tw) == len(self.hb.nonfree_flow_tw)
    2991        assert len(self.hb.nonfree_flow_curves) == len(self.hb.nonfree_flow_tw)
     
    3193        return
    3294
    33     def test_function_values(self):
     95    def test_hecras_internal_boundary_function_values(self):
     96
     97        self.hb = hecras_internal_boundary_function(self.input_hecras_file, verbose=False)
    3498
    3599        # Stationary states #
     
    90154        return
    91155
     156    def test_pumping_station_function(self):
     157
     158        domain = self.create_domain(wallHeight=10.,
     159            InitialOceanStage=4.,
     160            InitialLandStage=0.)
     161
     162        ps_function = pumping_station_function(
     163            domain=domain,
     164            pump_capacity=1.0,
     165            hw_to_start_pumping=0.0,
     166            hw_to_stop_pumping=-1.0,
     167            initial_pump_rate=0.,
     168            pump_rate_of_increase = 1.0,
     169            pump_rate_of_decrease = 1.0,
     170            verbose=False)
     171
     172
     173        # Pump should be starting at zero       
     174        assert(ps_function(1., 1.) == 0.)
     175   
     176        # As time increases, so will the pump rate
     177        domain.time = 0.5
     178        assert(numpy.allclose(ps_function(1., 1.), 0.5))
     179        assert(numpy.allclose(ps_function(1., 1.), 0.5))
     180        domain.time = 0.75
     181        assert(numpy.allclose(ps_function(1., 1.), 0.75))
     182
     183        # Pump rate maximum = 1.0
     184        domain.time = 1.5
     185        assert(numpy.allclose(ps_function(1., 1.), 1.0))
     186
     187        # Force the pump rate to zero, and it should start increasing as time
     188        # increases
     189        ps_function.pump_rate = 0.
     190        assert (ps_function(1., 1.) == 0.)
     191        domain.time = 1.75
     192        assert (ps_function(1., 1.) == 0.25)
     193        domain.time = 2.5
     194        assert (ps_function(1., 1.) == 1.0)
     195        # Can't increase any more
     196        domain.time = 3.0
     197        assert (ps_function(1., 1.) == 1.0)
     198       
     199        # Let's try to decrease the pump
     200        assert (ps_function(-1.5, 1.) == 1.0)
     201        domain.time = 3.5
     202        assert (ps_function(-1.5, 1.) == 0.5)
     203        domain.time = 3.9
     204
     205        assert (numpy.allclose(ps_function(-1.5, 1.), 0.1))
     206        domain.time = 4.0
     207        assert (numpy.allclose(ps_function(-1.5, 1.), 0.0))
     208        domain.time = 5.0
     209        assert (numpy.allclose(ps_function(-1.5, 1.), 0.0))
     210
     211       
     212        return
     213
    92214
    93215if __name__ == "__main__":
  • trunk/anuga_core/validation_tests/behaviour_only/bridge_hecras2/results.tex

    r9672 r9673  
    6161diffusion causes additional drag for overbank flows.
    6262
    63 Several additional factors will contribute to deviations between the models. ANUGA
    64 models cross-channel variations in the water surface elevation, which is
     63Several additional factors will contribute to deviations between the models.
     64ANUGA models cross-channel variations in the water surface elevation, which is
    6565assumed constant in HECRAS. The under-bridge flux in ANUGA is based on the
    6666water elevations in the central channel up-and-down stream of the bridge, thus
    67 the 2D representation will have an effect on the under-bridge flow. The models
     67the 2D representation will have an effect on the under-bridge flow. Further,
     68since the ANUGA model here is fairly coarsely resolved, the flow details near
     69the bridge are noticably affected by the choice of flow algorithm. The models
    6870also flux the momentum under the bridge in different ways. ANUGA's method is to
    6971compute the average momentum in each direction along the upstream bridge
  • trunk/anuga_work/development/gareth/tests/ras_bridge/channel_floodplain1.py

    r9672 r9673  
    187187    hecras_discharge_function,
    188188    exchange_lines=[bridge_in, bridge_out],
    189     enquiry_gap=5.0,
     189    enquiry_gap=0.01,
    190190    zero_outflow_momentum=False,
    191     smoothing_timescale=0.0,
     191    smoothing_timescale=30.0,
    192192    logging=True)
    193193   
Note: See TracChangeset for help on using the changeset viewer.