source: trunk/anuga_core/source/anuga/file/test_read_sww.py @ 7762

Last change on this file since 7762 was 7762, checked in by hudson, 14 years ago

All tests in file module pass.

File size: 6.9 KB
Line 
1#!/usr/bin/env python
2#
3
4# This file was reverted from changeset:5484 to changeset:5470 on 10th July
5# by Ole.
6
7import unittest
8import copy
9import os
10import numpy as num
11
12
13import sww
14               
15
16class Test_read_sww(unittest.TestCase):
17    # Class variable
18    verbose = False
19
20    def set_verbose(self):
21        Test_read_sww.verbose = True
22       
23    def setUp(self):
24        pass
25
26    def tearDown(self):
27        pass
28
29    def test_read_sww(self):
30        """
31        Save to an sww file and then read back the info.
32        Here we store the info "uniquely"
33        """
34
35        #---------------------------------------------------------------------
36        # Import necessary modules
37        #---------------------------------------------------------------------
38        from anuga.abstract_2d_finite_volumes.mesh_factory import \
39            rectangular_cross
40        from anuga.shallow_water.shallow_water_domain import Domain
41        from boundaries import Reflective_boundary
42        from anuga.abstract_2d_finite_volumes.generic_boundary_conditions\
43                            import Dirichlet_boundary, Time_boundary
44
45        #---------------------------------------------------------------------
46        # Setup computational domain
47        #---------------------------------------------------------------------
48        length = 8.
49        width = 4.
50        dx = dy = 2   # Resolution: Length of subdivisions on both axes
51       
52        inc = 0.05 # Elevation increment
53
54        points, vertices, boundary = rectangular_cross(int(length/dx), 
55                                                       int(width/dy),
56                                                       len1=length, 
57                                                       len2=width)
58        domain = Domain(points, vertices, boundary)
59        domain.set_name('read_sww_test'+str(domain.processor))  # Output name
60        domain.set_quantities_to_be_stored({'elevation': 2,
61                                            'stage': 2,
62                                            'xmomentum': 2,
63                                            'ymomentum': 2,
64                                            'friction': 1})
65
66        domain.set_store_vertices_uniquely(True)
67       
68        #---------------------------------------------------------------------
69        # Setup initial conditions
70        #---------------------------------------------------------------------
71        domain.set_quantity('elevation', 0.0)    # Flat bed initially
72        domain.set_quantity('friction', 0.01)    # Constant friction
73        domain.set_quantity('stage', 0.0)        # Dry initial condition
74
75        #------------------------------------------------------------------
76        # Setup boundary conditions
77        #------------------------------------------------------------------
78        Bi = Dirichlet_boundary([0.4, 0, 0])          # Inflow
79        Br = Reflective_boundary(domain)              # Solid reflective wall
80        Bo = Dirichlet_boundary([-5, 0, 0])           # Outflow
81
82        domain.set_boundary({'left': Bi, 'right': Bo, 'top': Br, 'bottom': Br})
83
84        #-------------------------------------------------------------------
85        # Evolve system through time
86        #-------------------------------------------------------------------
87
88        for t in domain.evolve(yieldstep=1, finaltime=4.0):
89            pass
90           
91           
92        # Check that quantities have been stored correctly   
93        source = domain.get_name() + '.sww'
94
95
96        #x = fid.variables['x'][:]
97        #y = fid.variables['y'][:]
98        #stage = fid.variables['stage'][:]
99        #elevation = fid.variables['elevation'][:]
100        #fid.close()
101                   
102        #assert len(stage.shape) == 2
103        #assert len(elevation.shape) == 2       
104       
105        #M, N = stage.shape
106               
107        sww_file = sww.Read_sww(source)
108
109        #print 'last frame number',sww_file.get_last_frame_number()
110
111        assert num.allclose(sww_file.x, domain.get_vertex_coordinates()[:,0])
112        assert num.allclose(sww_file.y, domain.get_vertex_coordinates()[:,1])
113
114
115        assert num.allclose(sww_file.time, [0.0, 1.0, 2.0, 3.0, 4.0])
116       
117        M = domain.get_number_of_triangles()
118       
119        assert num.allclose(num.reshape(num.arange(3*M), (M,3)), sww_file.vertices)
120
121        last_frame_number = sww_file.get_last_frame_number() 
122        assert last_frame_number == 4
123
124        assert num.allclose(sww_file.get_bounds(), [0.0, length, 0.0, width])
125
126        assert 'stage'     in sww_file.quantities.keys()
127        assert 'friction'  in sww_file.quantities.keys()
128        assert 'elevation' in sww_file.quantities.keys()
129        assert 'xmomentum' in sww_file.quantities.keys()
130        assert 'ymomentum' in sww_file.quantities.keys()
131
132
133        for qname, q in sww_file.read_quantities(last_frame_number).items():
134            assert num.allclose(domain.get_quantity(qname).get_values(), q)
135
136        #-----------------------------------------
137        # Start the evolution off again at frame 3
138        #-----------------------------------------
139        sww_file.read_quantities(last_frame_number-1)
140
141        points, vertices, boundary = rectangular_cross(int(length/dx), 
142                                                       int(width/dy),
143                                                       len1=length, 
144                                                       len2=width)
145        new_domain = Domain(points, vertices, boundary)
146        new_domain.set_quantities_to_be_stored(None)
147
148        new_domain.set_store_vertices_uniquely(True)
149
150        for qname, q in sww_file.read_quantities(last_frame_number-1).items():
151            new_domain.set_quantity(qname, q)   
152
153
154        new_domain.set_starttime(sww_file.get_time())
155        #------------------------------------------------------------------
156        # Setup boundary conditions
157        #------------------------------------------------------------------
158        Bi = Dirichlet_boundary([0.4, 0, 0])          # Inflow
159        Br = Reflective_boundary(new_domain)          # Solid reflective wall
160        Bo = Dirichlet_boundary([-5, 0, 0])           # Outflow
161
162        new_domain.set_boundary({'left': Bi, 'right': Bo, 'top': Br, 'bottom': Br})
163
164        #-------------------------------------------------------------------
165        # Evolve system through time
166        #-------------------------------------------------------------------
167
168        for t in new_domain.evolve(yieldstep=1.0, finaltime=1.0):
169             pass
170
171        # Compare  new_domain and domain quantities
172        for quantity in domain.get_quantity_names():
173            dv = domain.get_quantity(quantity).get_values()
174            ndv = new_domain.get_quantity(quantity).get_values()
175
176            assert num.allclose( dv, ndv)
177
178        # Clean up
179        os.remove(source)
180       
181
182if __name__ == "__main__":
183    suite = unittest.makeSuite(Test_read_sww, 'test')
184    runner = unittest.TextTestRunner() #verbosity=2)
185    runner.run(suite)
186   
Note: See TracBrowser for help on using the repository browser.