source: branches/numpy/anuga/abstract_2d_finite_volumes/test_pmesh2domain.py @ 6883

Last change on this file since 6883 was 6410, checked in by rwilson, 16 years ago

numpy changes.

File size: 9.5 KB
Line 
1#!/usr/bin/env python
2#
3
4import unittest
5
6#from anuga.pyvolution.pmesh2domain import *
7from pmesh2domain import *
8
9from anuga.shallow_water import Domain,\
10     Reflective_boundary, Dirichlet_boundary,\
11     Transmissive_boundary
12
13from anuga.coordinate_transforms.geo_reference import Geo_reference
14from anuga.pmesh.mesh import importMeshFromFile
15
16import numpy as num
17
18
19class Test_pmesh2domain(unittest.TestCase):
20
21    def setUp(self):
22        pass
23
24    def tearDown(self):
25        pass
26
27    def test_pmesh2Domain(self):
28         import os
29         import tempfile
30
31         fileName = tempfile.mktemp(".tsh")
32         file = open(fileName,"w")
33         file.write("4 3 # <vertex #> <x> <y> [attributes]\n \
340 0.0 0.0 0.0 0.0 0.01 \n \
351 1.0 0.0 10.0 10.0 0.02  \n \
362 0.0 1.0 0.0 10.0 0.03  \n \
373 0.5 0.25 8.0 12.0 0.04  \n \
38# Vert att title  \n \
39elevation  \n \
40stage  \n \
41friction  \n \
422 # <triangle #> [<vertex #>] [<neigbouring triangle #>]  \n\
430 0 3 2 -1  -1  1 dsg\n\
441 0 1 3 -1  0 -1   ole nielsen\n\
454 # <segment #> <vertex #>  <vertex #> [boundary tag] \n\
460 1 0 2 \n\
471 0 2 3 \n\
482 2 3 \n\
493 3 1 1 \n\
503 0 # <x> <y> [attributes] ...Mesh Vertices... \n \
510 216.0 -86.0 \n \
521 160.0 -167.0 \n \
532 114.0 -91.0 \n \
543 # <vertex #>  <vertex #> [boundary tag] ...Mesh Segments... \n \
550 0 1 0 \n \
561 1 2 0 \n \
572 2 0 0 \n \
580 # <x> <y> ...Mesh Holes... \n \
590 # <x> <y> <attribute>...Mesh Regions... \n \
600 # <x> <y> <attribute>...Mesh Regions, area... \n\
61#Geo reference \n \
6256 \n \
63140 \n \
64120 \n")
65         file.close()
66
67         tags = {}
68         b1 =  Dirichlet_boundary(conserved_quantities = num.array([0.0]))
69         b2 =  Dirichlet_boundary(conserved_quantities = num.array([1.0]))
70         b3 =  Dirichlet_boundary(conserved_quantities = num.array([2.0]))
71         tags["1"] = b1
72         tags["2"] = b2
73         tags["3"] = b3
74
75         domain = pmesh_to_domain_instance(fileName, Domain)
76         os.remove(fileName)
77         #print "domain.tagged_elements", domain.tagged_elements
78         ## check the quantities
79         #print domain.quantities['elevation'].vertex_values
80         answer = [[0., 8., 0.],
81                   [0., 10., 8.]]
82         assert num.allclose(domain.quantities['elevation'].vertex_values,
83                             answer)
84
85         #print domain.quantities['stage'].vertex_values
86         answer = [[0., 12., 10.],
87                   [0., 10., 12.]]
88         assert num.allclose(domain.quantities['stage'].vertex_values,
89                             answer)
90
91         #print domain.quantities['friction'].vertex_values
92         answer = [[0.01, 0.04, 0.03],
93                   [0.01, 0.02, 0.04]]
94         assert num.allclose(domain.quantities['friction'].vertex_values,
95                             answer)
96
97         #print domain.quantities['friction'].vertex_values
98         tagged_elements = domain.get_tagged_elements()
99         assert num.allclose(tagged_elements['dsg'][0],0)
100         assert num.allclose(tagged_elements['ole nielsen'][0],1)
101
102         self.failUnless( domain.boundary[(1, 0)]  == '1',
103                          "test_tags_to_boundaries  failed. Single boundary wasn't added.")
104         self.failUnless( domain.boundary[(1, 2)]  == '2',
105                          "test_tags_to_boundaries  failed. Single boundary wasn't added.")
106         self.failUnless( domain.boundary[(0, 1)]  == '3',
107                          "test_tags_to_boundaries  failed. Single boundary wasn't added.")
108         self.failUnless( domain.boundary[(0, 0)]  == 'exterior',
109                          "test_tags_to_boundaries  failed. Single boundary wasn't added.")
110         #print "domain.boundary",domain.boundary
111         self.failUnless( len(domain.boundary)  == 4,
112                          "test_pmesh2Domain Too many boundaries")
113         #FIXME change to use get_xllcorner
114         #print "d.geo_reference.xllcorner",domain.geo_reference.xllcorner
115         self.failUnless(domain.geo_reference.xllcorner  == 140.0,
116                          "bad geo_referece")
117    #************
118   
119    def test_pmesh2Domain_instance(self):
120         import os
121         import tempfile
122
123         fileName = tempfile.mktemp(".tsh")
124         file = open(fileName,"w")
125         file.write("4 3 # <vertex #> <x> <y> [attributes]\n \
1260 0.0 0.0 0.0 0.0 0.01 \n \
1271 1.0 0.0 10.0 10.0 0.02  \n \
1282 0.0 1.0 0.0 10.0 0.03  \n \
1293 0.5 0.25 8.0 12.0 0.04  \n \
130# Vert att title  \n \
131elevation  \n \
132stage  \n \
133friction  \n \
1342 # <triangle #> [<vertex #>] [<neigbouring triangle #>]  \n\
1350 0 3 2 -1  -1  1 dsg\n\
1361 0 1 3 -1  0 -1   ole nielsen\n\
1374 # <segment #> <vertex #>  <vertex #> [boundary tag] \n\
1380 1 0 2 \n\
1391 0 2 3 \n\
1402 2 3 \n\
1413 3 1 1 \n\
1423 0 # <x> <y> [attributes] ...Mesh Vertices... \n \
1430 216.0 -86.0 \n \
1441 160.0 -167.0 \n \
1452 114.0 -91.0 \n \
1463 # <vertex #>  <vertex #> [boundary tag] ...Mesh Segments... \n \
1470 0 1 0 \n \
1481 1 2 0 \n \
1492 2 0 0 \n \
1500 # <x> <y> ...Mesh Holes... \n \
1510 # <x> <y> <attribute>...Mesh Regions... \n \
1520 # <x> <y> <attribute>...Mesh Regions, area... \n\
153#Geo reference \n \
15456 \n \
155140 \n \
156120 \n")
157         file.close()
158
159         mesh_instance = importMeshFromFile(fileName)
160       
161         tags = {}
162         b1 =  Dirichlet_boundary(conserved_quantities = num.array([0.0]))
163         b2 =  Dirichlet_boundary(conserved_quantities = num.array([1.0]))
164         b3 =  Dirichlet_boundary(conserved_quantities = num.array([2.0]))
165         tags["1"] = b1
166         tags["2"] = b2
167         tags["3"] = b3
168
169         domain = pmesh_instance_to_domain_instance(mesh_instance, Domain)
170
171         os.remove(fileName)
172         #print "domain.tagged_elements", domain.tagged_elements
173         ## check the quantities
174         #print domain.quantities['elevation'].vertex_values
175         answer = [[0., 8., 0.],
176                   [0., 10., 8.]]
177         assert num.allclose(domain.quantities['elevation'].vertex_values,
178                             answer)
179
180         #print domain.quantities['stage'].vertex_values
181         answer = [[0., 12., 10.],
182                   [0., 10., 12.]]
183         assert num.allclose(domain.quantities['stage'].vertex_values,
184                             answer)
185
186         #print domain.quantities['friction'].vertex_values
187         answer = [[0.01, 0.04, 0.03],
188                   [0.01, 0.02, 0.04]]
189         assert num.allclose(domain.quantities['friction'].vertex_values,
190                             answer)
191
192         #print domain.quantities['friction'].vertex_values
193         tagged_elements = domain.get_tagged_elements()         
194         assert num.allclose(tagged_elements['dsg'][0],0)
195         assert num.allclose(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 = num.array([0.0]))
233         b2 =  Dirichlet_boundary(conserved_quantities = num.array([1.0]))
234         b3 =  Dirichlet_boundary(conserved_quantities = num.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#-------------------------------------------------------------
257
258if __name__ == "__main__":
259    suite = unittest.makeSuite(Test_pmesh2domain, 'test')
260    runner = unittest.TextTestRunner()
261    runner.run(suite)
Note: See TracBrowser for help on using the repository browser.