source: inundation/pyvolution/test_pmesh2domain.py @ 2689

Last change on this file since 2689 was 2533, checked in by ole, 18 years ago

Work on area histogram

File size: 9.7 KB
Line 
1#!/usr/bin/env python
2#
3
4import unittest
5
6from Numeric import allclose, array
7
8from shallow_water import Domain
9
10from pmesh2domain import *
11
12from shallow_water import Domain, Reflective_boundary, Dirichlet_boundary,\
13     Transmissive_boundary
14
15from coordinate_transforms.geo_reference import Geo_reference
16
17#This is making pyvolution dependent on pmesh.
18# not good.  this should be in a seperate package.(directory)
19from pmesh.mesh import importMeshFromFile
20
21class Test_pmesh2domain(unittest.TestCase):
22
23    def setUp(self):
24        pass
25
26    def tearDown(self):
27        pass
28
29    def test_pmesh2Domain(self):
30         import os
31         import tempfile
32
33         fileName = tempfile.mktemp(".tsh")
34         file = open(fileName,"w")
35         file.write("4 3 # <vertex #> <x> <y> [attributes]\n \
360 0.0 0.0 0.0 0.0 0.01 \n \
371 1.0 0.0 10.0 10.0 0.02  \n \
382 0.0 1.0 0.0 10.0 0.03  \n \
393 0.5 0.25 8.0 12.0 0.04  \n \
40# Vert att title  \n \
41elevation  \n \
42stage  \n \
43friction  \n \
442 # <triangle #> [<vertex #>] [<neigbouring triangle #>]  \n\
450 0 3 2 -1  -1  1 dsg\n\
461 0 1 3 -1  0 -1   ole nielsen\n\
474 # <segment #> <vertex #>  <vertex #> [boundary tag] \n\
480 1 0 2 \n\
491 0 2 3 \n\
502 2 3 \n\
513 3 1 1 \n\
523 0 # <x> <y> [attributes] ...Mesh Vertices... \n \
530 216.0 -86.0 \n \
541 160.0 -167.0 \n \
552 114.0 -91.0 \n \
563 # <vertex #>  <vertex #> [boundary tag] ...Mesh Segments... \n \
570 0 1 0 \n \
581 1 2 0 \n \
592 2 0 0 \n \
600 # <x> <y> ...Mesh Holes... \n \
610 # <x> <y> <attribute>...Mesh Regions... \n \
620 # <x> <y> <attribute>...Mesh Regions, area... \n\
63#Geo reference \n \
6456 \n \
65140 \n \
66120 \n")
67         file.close()
68
69         tags = {}
70         b1 =  Dirichlet_boundary(conserved_quantities = array([0.0]))
71         b2 =  Dirichlet_boundary(conserved_quantities = array([1.0]))
72         b3 =  Dirichlet_boundary(conserved_quantities = array([2.0]))
73         tags["1"] = b1
74         tags["2"] = b2
75         tags["3"] = b3
76
77         domain = pmesh_to_domain_instance(fileName, Domain)
78         os.remove(fileName)
79         #print "domain.tagged_elements", domain.tagged_elements
80         ## check the quantities
81         #print domain.quantities['elevation'].vertex_values
82         answer = [[0., 8., 0.],
83                   [0., 10., 8.]]
84         assert allclose(domain.quantities['elevation'].vertex_values,
85                        answer)
86
87         #print domain.quantities['stage'].vertex_values
88         answer = [[0., 12., 10.],
89                   [0., 10., 12.]]
90         assert allclose(domain.quantities['stage'].vertex_values,
91                        answer)
92
93         #print domain.quantities['friction'].vertex_values
94         answer = [[0.01, 0.04, 0.03],
95                   [0.01, 0.02, 0.04]]
96         assert allclose(domain.quantities['friction'].vertex_values,
97                        answer)
98
99         #print domain.quantities['friction'].vertex_values
100         assert allclose(domain.tagged_elements['dsg'][0],0)
101         assert allclose(domain.tagged_elements['ole nielsen'][0],1)
102
103         self.failUnless( domain.boundary[(1, 0)]  == '1',
104                          "test_tags_to_boundaries  failed. Single boundary wasn't added.")
105         self.failUnless( domain.boundary[(1, 2)]  == '2',
106                          "test_tags_to_boundaries  failed. Single boundary wasn't added.")
107         self.failUnless( domain.boundary[(0, 1)]  == '3',
108                          "test_tags_to_boundaries  failed. Single boundary wasn't added.")
109         self.failUnless( domain.boundary[(0, 0)]  == 'exterior',
110                          "test_tags_to_boundaries  failed. Single boundary wasn't added.")
111         #print "domain.boundary",domain.boundary
112         self.failUnless( len(domain.boundary)  == 4,
113                          "test_pmesh2Domain Too many boundaries")
114         #FIXME change to use get_xllcorner
115         #print "d.geo_reference.xllcorner",domain.geo_reference.xllcorner
116         self.failUnless(domain.geo_reference.xllcorner  == 140.0,
117                          "bad geo_referece")
118    #************
119   
120    def test_pmesh2Domain_instance(self):
121         import os
122         import tempfile
123
124         fileName = tempfile.mktemp(".tsh")
125         file = open(fileName,"w")
126         file.write("4 3 # <vertex #> <x> <y> [attributes]\n \
1270 0.0 0.0 0.0 0.0 0.01 \n \
1281 1.0 0.0 10.0 10.0 0.02  \n \
1292 0.0 1.0 0.0 10.0 0.03  \n \
1303 0.5 0.25 8.0 12.0 0.04  \n \
131# Vert att title  \n \
132elevation  \n \
133stage  \n \
134friction  \n \
1352 # <triangle #> [<vertex #>] [<neigbouring triangle #>]  \n\
1360 0 3 2 -1  -1  1 dsg\n\
1371 0 1 3 -1  0 -1   ole nielsen\n\
1384 # <segment #> <vertex #>  <vertex #> [boundary tag] \n\
1390 1 0 2 \n\
1401 0 2 3 \n\
1412 2 3 \n\
1423 3 1 1 \n\
1433 0 # <x> <y> [attributes] ...Mesh Vertices... \n \
1440 216.0 -86.0 \n \
1451 160.0 -167.0 \n \
1462 114.0 -91.0 \n \
1473 # <vertex #>  <vertex #> [boundary tag] ...Mesh Segments... \n \
1480 0 1 0 \n \
1491 1 2 0 \n \
1502 2 0 0 \n \
1510 # <x> <y> ...Mesh Holes... \n \
1520 # <x> <y> <attribute>...Mesh Regions... \n \
1530 # <x> <y> <attribute>...Mesh Regions, area... \n\
154#Geo reference \n \
15556 \n \
156140 \n \
157120 \n")
158         file.close()
159
160         mesh_instance = importMeshFromFile(fileName)
161       
162         tags = {}
163         b1 =  Dirichlet_boundary(conserved_quantities = array([0.0]))
164         b2 =  Dirichlet_boundary(conserved_quantities = array([1.0]))
165         b3 =  Dirichlet_boundary(conserved_quantities = array([2.0]))
166         tags["1"] = b1
167         tags["2"] = b2
168         tags["3"] = b3
169
170         domain = pmesh_instance_to_domain_instance(mesh_instance, Domain)
171
172         os.remove(fileName)
173         #print "domain.tagged_elements", domain.tagged_elements
174         ## check the quantities
175         #print domain.quantities['elevation'].vertex_values
176         answer = [[0., 8., 0.],
177                   [0., 10., 8.]]
178         assert allclose(domain.quantities['elevation'].vertex_values,
179                        answer)
180
181         #print domain.quantities['stage'].vertex_values
182         answer = [[0., 12., 10.],
183                   [0., 10., 12.]]
184         assert allclose(domain.quantities['stage'].vertex_values,
185                        answer)
186
187         #print domain.quantities['friction'].vertex_values
188         answer = [[0.01, 0.04, 0.03],
189                   [0.01, 0.02, 0.04]]
190         assert allclose(domain.quantities['friction'].vertex_values,
191                        answer)
192
193         #print domain.quantities['friction'].vertex_values
194         assert allclose(domain.tagged_elements['dsg'][0],0)
195         assert allclose(domain.tagged_elements['ole nielsen'][0],1)
196
197         self.failUnless( domain.boundary[(1, 0)]  == '1',
198                          "test_tags_to_boundaries  failed. Single boundary wasn't added.")
199         self.failUnless( domain.boundary[(1, 2)]  == '2',
200                          "test_tags_to_boundaries  failed. Single boundary wasn't added.")
201         self.failUnless( domain.boundary[(0, 1)]  == '3',
202                          "test_tags_to_boundaries  failed. Single boundary wasn't added.")
203         self.failUnless( domain.boundary[(0, 0)]  == 'exterior',
204                          "test_tags_to_boundaries  failed. Single boundary wasn't added.")
205         #print "domain.boundary",domain.boundary
206         self.failUnless( len(domain.boundary)  == 4,
207                          "test_pmesh2Domain Too many boundaries")
208         #FIXME change to use get_xllcorner
209         #print "d.geo_reference.xllcorner",domain.geo_reference.xllcorner
210         self.failUnless(domain.geo_reference.xllcorner  == 140.0,
211                          "bad geo_referece")
212         
213    #***********
214    def old_test_tags_to_boundaries (self):
215         meshDict = {}
216         p0 = [0.0, 0.0]
217         p1 = [1.0, 0.0]
218         p2 = [0.0, 1.0]
219         p3 = [0.646446609407, 0.353553390593]
220         meshDict['vertices'] = [p0,p1,p2,p3]
221         meshDict['vertex_attributes'] = [[0.0, 0.0,7.0],[10.0, 0.0,7.0],[0.0, 10.0,7.0],[6.46446609407, 3.53553390593,7.0]]
222         meshDict['triangles'] = [[0,3,2],[0,1,3]]
223         meshDict['triangle_tags'] = [6.6,6.6]
224         meshDict['triangle_neighbors'] = [[-1,-1,1],[-1,0,-1]]
225         meshDict['segments'] = [[1,0],[0,2],[2,3],[3,1]]
226         meshDict['segment_tags'] = [2,3,1,1]
227
228         domain = Domain.pmesh_dictionary_to_domain(meshDict)
229
230         #domain.set_tag_dict(tag_dict)
231         #Boundary tests
232         b1 =  Dirichlet_boundary(conserved_quantities = array([0.0]))
233         b2 =  Dirichlet_boundary(conserved_quantities = array([1.0]))
234         b3 =  Dirichlet_boundary(conserved_quantities = array([1.0]))
235         #test adding a boundary
236         tags = {}
237         tags[1] = b1
238         tags[2] = b2
239         tags[3] = b3
240         domain.set_boundary(tags)
241         inverted_id = Volume.instances[0].neighbours[0]
242         id = -(inverted_id+1)
243         boundary_value = Boundary_value.instances[id]
244         boundary_obj = boundary_value.boundary_object
245
246         self.failUnless( boundary_obj  == b1,
247                          "test_tags_to_boundaries  failed. Single boundary wasn't added.")
248
249         inverted_id = Volume.instances[0].neighbours[1]
250         id = -(inverted_id+1)
251         boundary_value = Boundary_value.instances[id]
252         boundary_obj = boundary_value.boundary_object
253         self.failUnless( boundary_obj  == b3,
254                          "test_tags_to_boundaries  failed. Single boundary wasn't added.")
255
256#-------------------------------------------------------------
257if __name__ == "__main__":
258    suite = unittest.makeSuite(Test_pmesh2domain, 'test')
259    runner = unittest.TextTestRunner()
260    runner.run(suite)
261
262
263
264
Note: See TracBrowser for help on using the repository browser.