source: anuga_core/source/anuga/shallow_water_balanced/test_swb_reporting.py @ 7573

Last change on this file since 7573 was 7573, checked in by steve, 14 years ago

Committing a version of shallow_water_balanced which passes it unit test
using a version of edge limiting which doesn't limit boundary edges. THis
is useful to allow linear functions to be reconstructed.

Had to play with the transmissive BC and use centroid values which is
set via domain set_centroid_transmissive_bc

File size: 5.5 KB
Line 
1#!/usr/bin/env python
2
3import unittest, os
4import os.path
5from math import pi, sqrt
6import tempfile
7
8from anuga.config import g, epsilon
9from anuga.config import netcdf_mode_r, netcdf_mode_w, netcdf_mode_a
10from anuga.utilities.numerical_tools import mean
11from anuga.utilities.polygon import is_inside_polygon
12from anuga.coordinate_transforms.geo_reference import Geo_reference
13from anuga.abstract_2d_finite_volumes.quantity import Quantity
14from anuga.geospatial_data.geospatial_data import Geospatial_data
15from anuga.abstract_2d_finite_volumes.mesh_factory import rectangular_cross
16
17from anuga.utilities.system_tools import get_pathname_from_package
18from swb_domain import *
19
20import numpy as num
21
22# Get gateway to C implementation of flux function for direct testing
23from shallow_water_ext import flux_function_central as flux_function
24
25
26
27
28class Test_swb_clean(unittest.TestCase):
29    def setUp(self):
30        pass
31
32    def tearDown(self):
33        pass
34
35
36    def test_extrema(self):
37        """Test that extrema of quantities are computed correctly
38        Extrema are updated at every *internal* timestep
39        """
40
41        from anuga.abstract_2d_finite_volumes.mesh_factory \
42                import rectangular_cross
43
44        initial_runup_height = -0.4
45        final_runup_height = -0.3
46
47        #--------------------------------------------------------------
48        # Setup computational domain
49        #--------------------------------------------------------------
50        N = 5
51        points, vertices, boundary = rectangular_cross(N, N)
52        domain = Domain(points, vertices, boundary)
53        domain.set_name('extrema_test')
54
55        #--------------------------------------------------------------
56        # Setup initial conditions
57        #--------------------------------------------------------------
58        def topography(x,y):
59            return -x/2                             # linear bed slope
60
61        domain.set_quantity('elevation', topography)    # function for elevation
62        domain.set_quantity('friction', 0.)             # Zero friction
63        # Constant negative initial stage
64        domain.set_quantity('stage', initial_runup_height)
65        domain.set_quantities_to_be_monitored(['stage', 'height'],
66                                              time_interval=[0.5, 2.7],
67                                              polygon=[[0,0], [0,1],
68                                                       [1,1], [1,0]])
69
70        assert len(domain.quantities_to_be_monitored) == 2
71        assert domain.quantities_to_be_monitored.has_key('stage')
72        assert domain.quantities_to_be_monitored.has_key('height')
73        for key in domain.quantities_to_be_monitored['stage'].keys():
74            assert domain.quantities_to_be_monitored['stage'][key] is None
75
76        #--------------------------------------------------------------
77        # Setup boundary conditions
78        #--------------------------------------------------------------
79        Br = Reflective_boundary(domain)              # Reflective wall
80        # Constant inflow
81        Bd = Dirichlet_boundary([final_runup_height, 0, 0])
82
83        # All reflective to begin with (still water)
84        domain.set_boundary({'left': Br, 'right': Br, 'top': Br, 'bottom': Br})
85
86        #--------------------------------------------------------------
87        # Let triangles adjust and check extrema
88        #--------------------------------------------------------------
89        for t in domain.evolve(yieldstep=0.1, finaltime=1.0):
90            domain.quantity_statistics() # Run it silently
91
92        #--------------------------------------------------------------
93        # Test extrema
94        #--------------------------------------------------------------
95        stage = domain.quantities_to_be_monitored['stage']
96        assert stage['min'] <= stage['max']
97
98        assert num.allclose(stage['min'], initial_runup_height,
99                            rtol=1.0/N)    # First order accuracy
100
101        depth = domain.quantities_to_be_monitored['height']
102
103        assert depth['min'] <= depth['max']
104        assert depth['min'] >= 0.0
105        assert depth['max'] >= 0.0
106
107        #--------------------------------------------------------------
108        # Update boundary to allow inflow
109        #--------------------------------------------------------------
110        domain.set_boundary({'right': Bd})
111
112        #--------------------------------------------------------------
113        # Evolve system through time
114        #--------------------------------------------------------------
115        for t in domain.evolve(yieldstep=0.1, finaltime=3.0):
116            domain.quantity_statistics()    # Run it silently
117
118        #--------------------------------------------------------------
119        # Test extrema again
120        #--------------------------------------------------------------
121        stage = domain.quantities_to_be_monitored['stage']
122        assert stage['min'] <= stage['max']
123
124        assert num.allclose(stage['min'], initial_runup_height,
125                            rtol = 1.0/N) # First order accuracy
126
127        depth = domain.quantities_to_be_monitored['height']
128        assert depth['min'] <= depth['max']
129        assert depth['min'] >= 0.0
130        assert depth['max'] >= 0.0
131
132        # Cleanup
133        os.remove(domain.get_name() + '.sww')
134
135
136
137
138
139#################################################################################
140
141if __name__ == "__main__":
142    suite = unittest.makeSuite(Test_swb_clean, 'test')
143    runner = unittest.TextTestRunner(verbosity=1)
144    runner.run(suite)
Note: See TracBrowser for help on using the repository browser.