Changeset 7496


Ignore:
Timestamp:
Sep 7, 2009, 6:51:23 PM (15 years ago)
Author:
ole
Message:

Fixed conflict

File:
1 edited

Legend:

Unmodified
Added
Removed
  • anuga_core/source/anuga/culvert_flows/test_culvert_class.py

    r7492 r7496  
    2222
    2323import numpy as num
     24
     25
     26# Helper functions
     27def run_culvert_flow_problem(depth):
     28    """Run flow with culvert given depth
     29    """
     30
     31    length = 40.
     32    width = 5.
     33
     34    dx = dy = 1           # Resolution: Length of subdivisions on both axes
     35
     36    points, vertices, boundary = rectangular_cross(int(length/dx),
     37                                                   int(width/dy),
     38                                                   len1=length,
     39                                                   len2=width)
     40    domain = Domain(points, vertices, boundary)   
     41    domain.set_name('Test_culvert_shallow') # Output name
     42    domain.set_default_order(2)
     43
     44
     45    #----------------------------------------------------------------------
     46    # Setup initial conditions
     47    #----------------------------------------------------------------------
     48
     49    def topography(x, y):
     50        """Set up a weir
     51       
     52        A culvert will connect either side
     53        """
     54        # General Slope of Topography
     55        z=-x/1000
     56       
     57        N = len(x)
     58        for i in range(N):
     59
     60           # Sloping Embankment Across Channel
     61            if 5.0 < x[i] < 10.1:
     62                # Cut Out Segment for Culvert face               
     63                if  1.0+(x[i]-5.0)/5.0 <  y[i]  < 4.0 - (x[i]-5.0)/5.0:
     64                   z[i]=z[i]
     65                else:
     66                   z[i] +=  0.5*(x[i] -5.0)    # Sloping Segment  U/S Face
     67            if 10.0 < x[i] < 12.1:
     68               z[i] +=  2.5                    # Flat Crest of Embankment
     69            if 12.0 < x[i] < 14.5:
     70                # Cut Out Segment for Culvert face               
     71                if  2.0-(x[i]-12.0)/2.5 <  y[i]  < 3.0 + (x[i]-12.0)/2.5:
     72                   z[i]=z[i]
     73                else:
     74                   z[i] +=  2.5-1.0*(x[i] -12.0) # Sloping D/S Face
     75                           
     76               
     77        return z
     78
     79
     80    domain.set_quantity('elevation', topography)
     81    domain.set_quantity('friction', 0.01)         # Constant friction
     82    domain.set_quantity('stage',
     83                        expression='elevation + %f' % depth) # Shallow initial condition
     84
     85    # Boyd culvert
     86    culvert = Culvert_flow(domain,
     87                           label='Culvert No. 1',
     88                           description='This culvert is a test unit 1.2m Wide by 0.75m High',   
     89                           end_point0=[9.0, 2.5],
     90                           end_point1=[13.0, 2.5],
     91                           width=1.20, height=0.75,
     92                           culvert_routine=boyd_generalised_culvert_model,
     93                           number_of_barrels=1,
     94                           update_interval=2,
     95                           verbose=False)
     96   
     97
     98    domain.forcing_terms.append(culvert)
     99   
     100
     101    #-----------------------------------------------------------------------
     102    # Setup boundary conditions
     103    #-----------------------------------------------------------------------
     104
     105    # Inflow based on Flow Depth and Approaching Momentum
     106
     107    Br = Reflective_boundary(domain)              # Solid reflective wall
     108    domain.set_boundary({'left': Br, 'right': Br, 'top': Br, 'bottom': Br})
     109
     110
     111   
     112    #-----------------------------------------------------------------------
     113    # Evolve system through time
     114    #-----------------------------------------------------------------------
     115
     116    #print 'depth', depth
     117    ref_volume = domain.get_quantity('stage').get_integral()
     118    for t in domain.evolve(yieldstep = 0.1, finaltime = 25):
     119        new_volume = domain.get_quantity('stage').get_integral()
     120       
     121        msg = ('Total volume has changed: Is %.8f m^3 should have been %.8f m^3'
     122               % (new_volume, ref_volume))
     123        assert num.allclose(new_volume, ref_volume), msg       
     124       
    24125
    25126
     
    100201                               end_point0=[9.0, 2.5],
    101202                               end_point1=[13.0, 2.5],
     203                               width=1.00,
    102204                               use_velocity_head=True,
    103205                               verbose=False)
     
    214316                               end_point0=[9.0, 2.5],
    215317                               end_point1=[13.0, 2.5],
     318                               height=0.75,
    216319                               verbose=False)
    217320
     
    243346
    244347   
    245     def test_that_culvert_rating_limits_flow_in_shallow_inlet_condition(self):
     348
     349       
     350    def test_that_culvert_flows_conserves_volume(self):
     351        """test_that_culvert_flows_conserves_volume
     352
     353        Test that culvert on a sloping dry bed limits flows when very little water
     354        is present at inlet.
     355
     356        Uses helper function: run_culvert_flow_problem(depth):
     357
     358        """
     359
     360        # Try this for a range of depths
     361        for depth in [0.1, 0.2, 0.5, 1.0]:
     362            run_culvert_flow_problem(depth)
     363
     364
     365    def OBSOLETE_XXXtest_that_culvert_rating_limits_flow_in_shallow_inlet_condition(self):
    246366        """test_that_culvert_rating_limits_flow_in_shallow_inlet_condition
    247367       
     
    251371        This one is using the rating curve variant
    252372        """
     373
     374       
    253375
    254376        path = get_pathname_from_package('anuga.culvert_flows')   
     
    307429        domain.set_quantity('stage',
    308430                            expression='elevation + 0.1') # Shallow initial condition
    309                            
    310 
    311         filename = os.path.join(path, 'example_rating_curve.csv')
     431
     432        # Boyd culvert
    312433        culvert = Culvert_flow(domain,
    313                                culvert_description_filename=filename,
     434                               label='Culvert No. 1',
     435                               description='This culvert is a test unit 1.2m Wide by 0.75m High',   
    314436                               end_point0=[9.0, 2.5],
    315437                               end_point1=[13.0, 2.5],
    316                                trigger_depth=0.01,
     438                               width=1.20, height=0.75,
     439                               culvert_routine=boyd_generalised_culvert_model,
     440                               number_of_barrels=1,
     441                               update_interval=2,
    317442                               verbose=False)
     443       
     444        # Rating curve
     445        #filename = os.path.join(path, 'example_rating_curve.csv')
     446        #culvert = Culvert_flow(domain,
     447        #                       culvert_description_filename=filename,
     448        #                       end_point0=[9.0, 2.5],
     449        #                       end_point1=[13.0, 2.5],
     450        #                       trigger_depth=0.01,
     451        #                       verbose=False)
    318452
    319453        domain.forcing_terms.append(culvert)
     
    335469        #-----------------------------------------------------------------------
    336470
     471        print 'depth', 0.1
    337472        ref_volume = domain.get_quantity('stage').get_integral()
    338473        for t in domain.evolve(yieldstep = 0.1, finaltime = 25):
    339474            new_volume = domain.get_quantity('stage').get_integral()
    340475
    341 
    342476            msg = ('Total volume has changed: Is %.8f m^3 should have been %.8f m^3'
    343477                   % (new_volume, ref_volume))
     
    346480       
    347481        return
    348         # Now try this for a range of other depths
    349         for depth in [1.0, 0.5, 0.2, 0.1, 0.05]:
     482        # Now try this again for a depth of 10 cm and for a range of other depths
     483        for depth in [0.1, 0.2, 0.5, 1.0]:
     484            print 'depth', depth
    350485            domain.set_time(0.0)
    351486           
     487            domain.set_quantity('elevation', topography)
     488            domain.set_quantity('friction', 0.01)         # Constant friction
    352489            domain.set_quantity('stage',
    353490                                expression='elevation + %f' % depth)
     
    358495                new_volume = domain.get_quantity('stage').get_integral()
    359496           
    360                 msg = 'Total volume has changed: Is %.2f m^3 should have been %.2f m^3'\
     497                msg = 'Total volume has changed: Is %.8f m^3 should have been %.8f m^3'\
    361498                    % (new_volume, ref_volume)
    362499
     
    439576                               end_point0=[9.0, 2.5],
    440577                               end_point1=[13.0, 2.5],
    441                                width=1.20,height=0.75,
     578                               width=1.20, height=0.75,
    442579                               culvert_routine=boyd_generalised_culvert_model,
    443580                               number_of_barrels=1,
    444581                               update_interval=2,
    445                                verbose=True)
     582                               verbose=False)
    446583       
    447584        domain.forcing_terms.append(culvert)
     
    543680                               culvert_routine=boyd_generalised_culvert_model,       
    544681                               number_of_barrels=1,
    545                                verbose=True)
     682                               verbose=False)
    546683                               
    547684       
     
    551688        culvert(domain)
    552689   
     690
     691   
     692
     693    def test_momentum_jet(self):
     694        """test_momentum_jet
     695       
     696        Test that culvert_class can accept keyword use_momentum_jet
     697        This does not yet imply that the values have been tested. FIXME
     698        """
     699
     700
     701        length = 40.
     702        width = 5.
     703
     704        dx = dy = 1           # Resolution: Length of subdivisions on both axes
     705
     706        points, vertices, boundary = rectangular_cross(int(length/dx),
     707                                                       int(width/dy),
     708                                                       len1=length,
     709                                                       len2=width)
     710        domain = Domain(points, vertices, boundary)   
     711        domain.set_name('Test_culvert_shallow') # Output name
     712        domain.set_default_order(2)
     713
     714
     715        #----------------------------------------------------------------------
     716        # Setup initial conditions
     717        #----------------------------------------------------------------------
     718
     719        def topography(x, y):
     720            """Set up a weir
     721           
     722            A culvert will connect either side
     723            """
     724            # General Slope of Topography
     725            z=-x/1000
     726           
     727            N = len(x)
     728            for i in range(N):
     729
     730               # Sloping Embankment Across Channel
     731                if 5.0 < x[i] < 10.1:
     732                    # Cut Out Segment for Culvert face               
     733                    if  1.0+(x[i]-5.0)/5.0 <  y[i]  < 4.0 - (x[i]-5.0)/5.0:
     734                       z[i]=z[i]
     735                    else:
     736                       z[i] +=  0.5*(x[i] -5.0)    # Sloping Segment  U/S Face
     737                if 10.0 < x[i] < 12.1:
     738                   z[i] +=  2.5                    # Flat Crest of Embankment
     739                if 12.0 < x[i] < 14.5:
     740                    # Cut Out Segment for Culvert face               
     741                    if  2.0-(x[i]-12.0)/2.5 <  y[i]  < 3.0 + (x[i]-12.0)/2.5:
     742                       z[i]=z[i]
     743                    else:
     744                       z[i] +=  2.5-1.0*(x[i] -12.0) # Sloping D/S Face
     745                                   
     746                       
     747            return z
     748
     749
     750        domain.set_quantity('elevation', topography)
     751        domain.set_quantity('friction', 0.01)         # Constant friction
     752        domain.set_quantity('stage',
     753                            expression='elevation + 1.0') # Shallow initial condition
     754
     755        # Boyd culvert
     756        culvert = Culvert_flow(domain,
     757                               label='Culvert No. 1',
     758                               description='This culvert is a test unit 1.2m Wide by 0.75m High',   
     759                               end_point0=[9.0, 2.5],
     760                               end_point1=[13.0, 2.5],
     761                               width=1.20, height=0.75,
     762                               culvert_routine=boyd_generalised_culvert_model,
     763                               number_of_barrels=1,
     764                               use_momentum_jet=True,
     765                               update_interval=2,
     766                               verbose=False)
     767   
     768
     769        domain.forcing_terms.append(culvert)
     770
     771       
     772        # Call
     773        culvert(domain)
     774   
     775
     776        #-----------------------------------------------------------------------
     777        # Setup boundary conditions
     778        #-----------------------------------------------------------------------
     779
     780
     781        Br = Reflective_boundary(domain)              # Solid reflective wall
     782        domain.set_boundary({'left': Br, 'right': Br, 'top': Br, 'bottom': Br})
     783
     784        #-----------------------------------------------------------------------
     785        # Evolve system through time
     786        #-----------------------------------------------------------------------
     787
     788        ref_volume = domain.get_quantity('stage').get_integral()
     789        for t in domain.evolve(yieldstep = 0.1, finaltime = 25):
     790            new_volume = domain.get_quantity('stage').get_integral()
     791           
     792            msg = ('Total volume has changed: Is %.8f m^3 should have been %.8f m^3'
     793                   % (new_volume, ref_volume))
     794            assert num.allclose(new_volume, ref_volume), msg       
     795       
     796
     797       
     798
    553799               
    554800#-------------------------------------------------------------
Note: See TracChangeset for help on using the changeset viewer.