source: anuga_core/source/anuga/pmesh/test_mesh_interface.py @ 7699

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

Added breaklines to high level interface for task 305.
File size: 34.7 KB
Line 
1#!/usr/bin/env python
2
3import tempfile
4import unittest
5import os
6from anuga.pmesh.mesh import importMeshFromFile
7from anuga.pmesh.mesh_interface import create_mesh_from_regions
8
9from anuga.pmesh.mesh_interface import _create_mesh_from_regions
10
11from load_mesh.loadASCII import *
12from anuga.utilities.polygon import is_inside_polygon
13from anuga.coordinate_transforms.geo_reference import Geo_reference,DEFAULT_ZONE
14
15
16class TestCase(unittest.TestCase):
17    def setUp(self):
18        pass
19
20    def tearDown(self):
21        pass
22
23    def test_create_mesh_from_regions(self):
24        x=-500
25        y=-1000
26        mesh_geo = geo_reference=Geo_reference(56, x, y)
27
28        # These are the absolute values
29        polygon_absolute = [[0,0], [100,0], [100,100], [0,100]]
30
31        x_p = -10
32        y_p = -40
33        geo_ref_poly = Geo_reference(56, x_p, y_p)
34        polygon = geo_ref_poly.change_points_geo_ref(polygon_absolute)
35
36        boundary_tags = {'walls': [0,1], 'bom': [2,3]}
37
38        inner1_polygon_absolute = [[10,10], [20,10], [20,20], [10,20]]
39        inner1_polygon = geo_ref_poly.\
40                            change_points_geo_ref(inner1_polygon_absolute)
41
42        inner2_polygon_absolute = [[30,30], [40,30], [40,40], [30,40]]
43        inner2_polygon = geo_ref_poly.\
44                            change_points_geo_ref(inner2_polygon_absolute)
45
46        interior_regions = [(inner1_polygon, 5), (inner2_polygon, 10)]
47
48        m = create_mesh_from_regions(polygon,
49                                     boundary_tags,
50                                     10000000,
51                                     interior_regions=interior_regions,
52                                     poly_geo_reference=geo_ref_poly,
53                                     mesh_geo_reference=mesh_geo)
54
55        # Test the mesh instance
56        self.failUnless(len(m.regions)==3, 'FAILED!')
57        segs = m.getUserSegments()
58        self.failUnless(len(segs)==12, 'FAILED!')
59        self.failUnless(len(m.userVertices)==12, 'FAILED!')
60        self.failUnless(segs[0].tag=='walls', 'FAILED!')
61        self.failUnless(segs[1].tag=='walls', 'FAILED!')
62        self.failUnless(segs[2].tag=='bom', 'FAILED!')
63        self.failUnless(segs[3].tag=='bom', 'FAILED!')
64
65        # Assuming the order of the region points is known.
66        # (This isn't true, if you consider create_mesh_from_regions
67        # a black box)
68        poly_point = m.getRegions()[0]
69
70        # poly_point values are relative to the mesh geo-ref
71        # make them absolute
72        msg = ('Expected point (%s,%s) to be inside polygon %s'
73               % (str(poly_point.x+x), str(poly_point.y+y),
74                  str(polygon_absolute)))
75        self.failUnless(is_inside_polygon([poly_point.x+x, poly_point.y+y],
76                                          polygon_absolute, closed=False),
77                        msg)
78
79        # Assuming the order of the region points is known.
80        # (This isn't true, if you consider create_mesh_from_regions
81        # a black box)
82        poly_point = m.getRegions()[1]
83
84        # poly_point values are relative to the mesh geo-ref
85        # make them absolute
86        self.failUnless(is_inside_polygon([poly_point.x+x, poly_point.y+y],
87                                          inner1_polygon_absolute,
88                                          closed=False),
89                        'FAILED!')
90
91        # Assuming the order of the region points is known.
92        # (This isn't true, if you consider create_mesh_from_regions
93        # a black box)
94        poly_point = m.getRegions()[2]
95
96        # poly_point values are relative to the mesh geo-ref
97        # make them absolute
98        self.failUnless(is_inside_polygon([poly_point.x+x, poly_point.y+y],
99                                          inner2_polygon_absolute,
100                                          closed=False),
101                        'FAILED!')
102
103    def test_create_mesh_from_regions_with_caching(self):
104        x=-500
105        y=-1000
106        mesh_geo = geo_reference=Geo_reference(56, x, y)
107
108        # These are the absolute values
109        polygon_absolute = [[0,0], [100,0], [100,100], [0,100]]
110
111        x_p = -10
112        y_p = -40
113        geo_ref_poly = Geo_reference(56, x_p, y_p)
114        polygon = geo_ref_poly.change_points_geo_ref(polygon_absolute)
115
116        boundary_tags = {'walls': [0,1], 'bom': [2,3]}
117
118        inner1_polygon_absolute = [[10,10], [20,10], [20,20], [10,20]]
119        inner1_polygon = geo_ref_poly.\
120                            change_points_geo_ref(inner1_polygon_absolute)
121
122        inner2_polygon_absolute = [[30,30], [40,30], [40,40], [30,40]]
123        inner2_polygon = geo_ref_poly.\
124                             change_points_geo_ref(inner2_polygon_absolute)
125
126        interior_regions = [(inner1_polygon, 5), (inner2_polygon, 10)]
127
128        interior_holes = None
129
130        # Clear cache first
131        from anuga.caching import cache
132
133        cache(_create_mesh_from_regions,
134              (polygon, boundary_tags),
135              {'minimum_triangle_angle': 28.0,
136               'maximum_triangle_area': 10000000,
137               'interior_regions': interior_regions,
138               'interior_holes': interior_holes,
139               'poly_geo_reference': geo_ref_poly,
140               'mesh_geo_reference': mesh_geo,
141               'verbose': False},
142              verbose=False,
143              clear=1)
144
145        m = create_mesh_from_regions(polygon,
146                                     boundary_tags,
147                                     maximum_triangle_area=10000000,
148                                     interior_regions=interior_regions,
149                                     poly_geo_reference=geo_ref_poly,
150                                     mesh_geo_reference=mesh_geo,
151                                     verbose=False,
152                                     use_cache=True)
153
154
155        # Test the mesh instance
156        self.failUnless(len(m.regions)==3, 'FAILED!')
157        segs = m.getUserSegments()
158        self.failUnless(len(segs)==12, 'FAILED!')
159        self.failUnless(len(m.userVertices)==12, 'FAILED!')
160        self.failUnless(segs[0].tag=='walls', 'FAILED!')
161        self.failUnless(segs[1].tag=='walls', 'FAILED!')
162        self.failUnless(segs[2].tag=='bom', 'FAILED!')
163        self.failUnless(segs[3].tag=='bom', 'FAILED!')
164
165        # Assuming the order of the region points is known.
166        # (This isn't true, if you consider create_mesh_from_regions
167        # a black box)
168        poly_point = m.getRegions()[0]
169
170        # poly_point values are relative to the mesh geo-ref
171        # make them absolute
172        self.failUnless(is_inside_polygon([poly_point.x+x, poly_point.y+y],
173                                          polygon_absolute,
174                                          closed=False),
175                        'FAILED!')
176
177        # Assuming the order of the region points is known.
178        # (This isn't true, if you consider create_mesh_from_regions
179        # a black box)
180        poly_point = m.getRegions()[1]
181
182        # poly_point values are relative to the mesh geo-ref
183        # make them absolute
184        self.failUnless(is_inside_polygon([poly_point.x+x, poly_point.y+y],
185                                          inner1_polygon_absolute,
186                                          closed=False),
187                        'FAILED!')
188
189        # Assuming the order of the region points is known.
190        # (This isn't true, if you consider create_mesh_from_regions
191        # a black box)
192        poly_point = m.getRegions()[2]
193
194        # poly_point values are relative to the mesh geo-ref
195        # make them absolute
196        self.failUnless(is_inside_polygon([poly_point.x+x, poly_point.y+y],
197                                          inner2_polygon_absolute,
198                                          closed=False),
199                        'FAILED!')
200
201        # Now create m using cached values
202        m_cache = create_mesh_from_regions(polygon,
203                                           boundary_tags,
204                                           10000000,
205                                           interior_regions=interior_regions,
206                                           poly_geo_reference=geo_ref_poly,
207                                           mesh_geo_reference=mesh_geo,
208                                           verbose=False,
209                                           use_cache=True)
210
211    def test_create_mesh_from_regions2(self):
212        # These are the absolute values
213        min_x = -10
214        min_y = -88
215        polygon_absolute = [[min_x,min_y], [1000,100], [1000,1000], [100,1000]]
216
217        x_p = -10
218        y_p = -40
219        zone = 808
220        geo_ref_poly = Geo_reference(zone, x_p, y_p)
221        polygon = geo_ref_poly.change_points_geo_ref(polygon_absolute)
222
223        boundary_tags = {'walls': [0,1], 'bom': [2,3]}
224
225        inner1_polygon_absolute = [[10,10], [20,10], [20,20], [10,20]]
226        inner1_polygon = geo_ref_poly.\
227                            change_points_geo_ref(inner1_polygon_absolute)
228
229        inner2_polygon_absolute = [[30,30], [40,30], [40,40], [30,40]]
230        inner2_polygon = geo_ref_poly.\
231                            change_points_geo_ref(inner2_polygon_absolute)
232
233        interior_regions = [(inner1_polygon, 5), (inner2_polygon, 10)]
234        m = create_mesh_from_regions(polygon,
235                                     boundary_tags,
236                                     10000000,
237                                     interior_regions=interior_regions,
238                                     poly_geo_reference=geo_ref_poly)
239
240        # Test the mesh instance
241        self.failUnless(len(m.regions)==3, 'FAILED!')
242        segs = m.getUserSegments()
243        self.failUnless(len(segs)==12, 'FAILED!')
244        self.failUnless(len(m.userVertices)==12, 'FAILED!')
245        self.failUnless(segs[0].tag=='walls', 'FAILED!')
246        self.failUnless(segs[1].tag=='walls', 'FAILED!')
247        self.failUnless(segs[2].tag=='bom', 'FAILED!')
248        self.failUnless(segs[3].tag=='bom', 'FAILED!')
249        self.failUnless(m.geo_reference.get_zone()==zone, 'FAILED!')
250        self.failUnless(m.geo_reference.get_xllcorner()==min_x, 'FAILED!')
251        self.failUnless(m.geo_reference.get_yllcorner()==min_y, 'FAILED!')             
252               
253    def test_create_mesh_from_regions3(self):
254        # These are the absolute values
255        min_x = -10
256        min_y = -88
257        polygon = [[min_x,min_y], [1000,100], [1000,1000], [100,1000]]
258
259
260        x_p = -10
261        y_p = -40
262        geo_ref_poly = Geo_reference(56, x_p, y_p)
263
264        boundary_tags = {'walls': [0,1], 'bom': [2,3]}
265
266        inner1_polygon_absolute = [[10,10], [20,10], [20,20], [10,20]]
267        inner1_polygon = geo_ref_poly.\
268                            change_points_geo_ref(inner1_polygon_absolute)
269
270        inner2_polygon_absolute = [[30,30], [40,30], [40,40], [30,40]]
271        inner2_polygon = geo_ref_poly.\
272                            change_points_geo_ref(inner2_polygon_absolute)
273
274        interior_regions = [(inner1_polygon, 5), (inner2_polygon, 10)]
275        m = create_mesh_from_regions(polygon,
276                                     boundary_tags,
277                                     10000000,
278                                     interior_regions=interior_regions)
279
280        # Test the mesh instance
281        self.failUnless(len(m.regions) == 3, 'FAILED!')
282        segs = m.getUserSegments()
283        self.failUnless(len(segs) == 12, 'FAILED!')
284        self.failUnless(len(m.userVertices) == 12, 'FAILED!')
285        self.failUnless(segs[0].tag == 'walls', 'FAILED!')
286        self.failUnless(segs[1].tag == 'walls', 'FAILED!')
287        self.failUnless(segs[2].tag == 'bom', 'FAILED!')
288        self.failUnless(segs[3].tag == 'bom', 'FAILED!')
289        self.failUnless(m.geo_reference.get_zone() == DEFAULT_ZONE, 'FAILED!')
290        self.failUnless(m.geo_reference.get_xllcorner() == min_x, 'FAILED!')
291        self.failUnless(m.geo_reference.get_yllcorner() == min_y, 'FAILED!')
292
293    def test_create_mesh_from_regions4(self):
294        file_name = tempfile.mktemp('.tsh')
295
296        # These are the absolute values
297        density_outer = 1000
298        min_outer = 0
299        max_outer = 1000
300        polygon_outer = [[min_outer,min_outer], [max_outer,min_outer],
301                         [max_outer,max_outer], [min_outer,max_outer]]
302
303        density_inner1 = 10000000
304        inner_buffer = 100
305        min_inner1 = min_outer + inner_buffer
306        max_inner1 = max_outer - inner_buffer
307        inner1_polygon = [[min_inner1,min_inner1], [max_inner1,min_inner1],
308                          [max_inner1,max_inner1], [min_inner1,max_inner1]]
309
310        boundary_tags = {'walls': [0,1], 'bom': [2,3]}
311
312        interior_regions = [(inner1_polygon, density_inner1)]
313        create_mesh_from_regions(polygon_outer,
314                                 boundary_tags,
315                                 density_outer,
316                                 interior_regions=interior_regions,
317                                 filename=file_name, verbose=False)
318
319        m = importMeshFromFile(file_name)
320
321        self.failUnless(len(m.getTriangulation()) <= 900,
322                        'Test mesh interface failed!')
323        self.failUnless(len(m.getTriangulation()) >= 200,
324                        'Test mesh interface failed!')
325
326        create_mesh_from_regions(polygon_outer,
327                                 boundary_tags,
328                                 interior_regions=interior_regions,
329                                 filename=file_name,
330                                 verbose=False)
331
332        m = importMeshFromFile(file_name)
333
334        self.failUnless(len(m.getTriangulation()) <= 100,
335                        'Test mesh interface failed!')
336
337        os.remove(file_name)
338
339    def test_create_mesh_from_regions5(self):
340        file_name = tempfile.mktemp('.tsh')
341
342        # These are the absolute values
343        density_outer = 10000000
344        min_outer = 0
345        max_outer = 1000
346        polygon_outer = [[min_outer,min_outer], [max_outer,min_outer],
347                         [max_outer,max_outer], [min_outer,max_outer]]
348
349        density_inner1 = 1000
350        inner_buffer = 100
351        min_inner1 = min_outer + inner_buffer
352        max_inner1 = max_outer - inner_buffer
353        inner1_polygon = [[min_inner1,min_inner1], [max_inner1,min_inner1],
354                          [max_inner1,max_inner1], [min_inner1,max_inner1]]
355
356        boundary_tags = {'walls': [0,1], 'bom': [2,3]}
357
358        interior_regions = [(inner1_polygon, density_inner1)]
359        create_mesh_from_regions(polygon_outer,
360                                 boundary_tags,
361                                 density_outer,
362                                 interior_regions=interior_regions,
363                                 filename=file_name,
364                                 verbose=False)
365
366        m = importMeshFromFile(file_name)
367        self.failUnless(len(m.getTriangulation()) <= 2000,
368                        'Test mesh interface failed!')
369        self.failUnless(len(m.getTriangulation()) >= 900,
370                        'Test mesh interface failed!')
371
372        os.remove(file_name)
373
374    def test_create_mesh_from_regions6(self):
375        file_name = tempfile.mktemp('.tsh')
376
377        # These are the absolute values
378        density_outer = 1000
379        min_outer = 0
380        max_outer = 1000
381        polygon_outer = [[min_outer,min_outer], [max_outer,min_outer],
382                         [max_outer,max_outer], [min_outer,max_outer]]
383
384        delta = 10
385        density_inner1 = 1000
386        min_inner1 = min_outer + delta
387        max_inner1 = max_outer - delta
388        inner1_polygon = [[min_inner1,min_inner1], [max_inner1,min_inner1],
389                          [max_inner1,max_inner1], [min_inner1,max_inner1]]
390
391        density_inner2 = 10000000
392        min_inner2 = min_outer +  2*delta
393        max_inner2 = max_outer -  2*delta
394        inner2_polygon = [[min_inner2,min_inner2], [max_inner2,min_inner2],
395                          [max_inner2,max_inner2], [min_inner2,max_inner2]]
396
397        boundary_tags = {'walls': [0,1], 'bom': [2,3]}
398
399        interior_regions = [(inner1_polygon, density_inner1),
400                            (inner2_polygon, density_inner2)]
401        create_mesh_from_regions(polygon_outer,
402                                 boundary_tags,
403                                 density_outer,
404                                 interior_regions=interior_regions,
405                                 filename=file_name,
406                                 verbose=False)
407
408        m = importMeshFromFile(file_name)
409        self.failUnless(len(m.getTriangulation()) <= 2000,
410                        'Test mesh interface failed!')
411        self.failUnless(len(m.getTriangulation()) >= 900,
412                        'Test mesh interface failed!')
413
414        os.remove(file_name)
415
416    def test_create_mesh_from_regions7(self):
417        file_name = tempfile.mktemp('.tsh')
418
419        # These are the absolute values
420        density_outer = 1001
421        min_outer = 0
422        max_outer = 1000
423        polygon_outer = [[min_outer,min_outer], [max_outer,min_outer],
424                         [max_outer,max_outer], [min_outer,max_outer]]
425
426        delta = 10
427        density_inner1 = 100000000
428        min_inner1 = min_outer + delta
429        max_inner1 = max_outer - delta
430        inner1_polygon = [[min_inner1,min_inner1], [max_inner1,min_inner1],
431                          [max_inner1,max_inner1], [min_inner1,max_inner1]]
432
433        density_inner2 = 1000
434        min_inner2 = min_outer +  2*delta
435        max_inner2 = max_outer -  2*delta
436        inner2_polygon = [[min_inner2,min_inner2], [max_inner2,min_inner2],
437                          [max_inner2,max_inner2], [min_inner2,max_inner2]]
438
439        boundary_tags = {'walls': [0,1], 'bom': [2,3]}
440
441        # Note the list order is important
442        # The last region added will be the region triangle uses,
443        # if two regions points are in the same bounded area.
444        interior_regions = [(inner2_polygon, density_inner2),
445                            (inner1_polygon, density_inner1)]
446        create_mesh_from_regions(polygon_outer,
447                                 boundary_tags,
448                                 density_outer,
449                                 interior_regions=interior_regions,
450                                 filename=file_name,
451                                 verbose=False)
452
453        m = importMeshFromFile(file_name)
454        self.failUnless(len(m.getTriangulation()) <= 3000,
455                        'Test mesh interface failed!')
456        self.failUnless(len(m.getTriangulation()) >= 2000,
457                        'Test mesh interface failed!')
458
459        os.remove(file_name)
460
461    def test_create_mesh_from_regions_interior_regions(self):
462        '''Test that create_mesh_from_regions fails when an interior
463        region is outside bounding polygon.
464        '''
465
466        # These are the absolute values
467        min_x = 10
468        min_y = 88
469        polygon = [[min_x,min_y], [1000,100], [1000,1000], [100,1000]]
470
471        boundary_tags = {'walls': [0,1], 'bom': [2,3]}
472
473        # This one is inside bounding polygon - should pass
474        inner_polygon = [[800,400], [900,500], [800,600]]
475
476        interior_regions = [(inner_polygon, 5)]
477        m = create_mesh_from_regions(polygon,
478                                     boundary_tags,
479                                     10000000,
480                                     interior_regions=interior_regions)
481
482        # This one sticks outside bounding polygon - should fail
483        inner_polygon = [[800,400], [900,500], [800,600], [200, 995]]
484        inner_polygon1 = [[800,400], [1100,500], [800,600]]
485        interior_regions = [[inner_polygon, 50], [inner_polygon1, 50]]
486
487        try:
488            m = create_mesh_from_regions(polygon,
489                                         boundary_tags,
490                                         10000000,
491                                         interior_regions=interior_regions,
492                                         verbose=False)
493        except:
494            pass
495        else:
496            msg = 'Interior polygon sticking outside bounding polygon should '
497            msg += 'cause an Exception to be raised'
498            raise Exception, msg
499
500    def test_create_mesh_from_regions_interior_regions1(self):
501        '''Test that create_mesh_from_regions fails
502        when an interior region is outside bounding polygon.
503        '''
504
505        # These are the values
506        d0 = [310000, 7690000]
507        d1 = [280000, 7690000]
508        d2 = [270000, 7645000]
509        d3 = [240000, 7625000]
510        d4 = [270000, 7580000]
511        d5 = [300000, 7590000]
512        d6 = [340000, 7610000]
513        poly_all = [d0, d1, d2, d3, d4, d5, d6]
514
515        i0 = [304000, 7607000]
516        i1 = [302000, 7605000]
517        i2 = [304000, 7603000]
518        i3 = [307000, 7602000]
519        i4 = [309000, 7603000]
520#        i4 = [310000, 7580000]
521        i5 = [307000, 7606000]
522        poly_onslow = [i0, i1, i2, i3, i4, i5]
523
524        # Thevenard Island
525        j0 = [294000, 7629000]
526        j1 = [285000, 7625000]
527        j2 = [294000, 7621000]
528        j3 = [299000, 7625000]
529        poly_thevenard = [j0, j1, j2, j3]
530
531        # med res around onslow
532        l0 = [300000, 7610000]
533        l1 = [285000, 7600000]
534        l2 = [300000, 7597500]
535        l3 = [310000, 7770000] # this one is outside
536#        l3 = [310000, 7630000] # this one is NOT outside
537        l4 = [315000, 7610000]
538        poly_coast = [l0, l1, l2, l3, l4]
539
540        # general coast and local area to onslow region
541        m0 = [270000, 7581000]
542        m1 = [300000, 7591000]
543        m2 = [339000, 7610000]
544        m3 = [330000, 7630000]
545        m4 = [290000, 7640000]
546        m5 = [260000, 7600000]
547        poly_region = [m0, m1, m2, m3, m4, m5]
548
549        # This one sticks outside bounding polygon - should fail
550        interior_regions = [[poly_onslow, 50000], [poly_region, 50000],
551                            [poly_coast, 100000], [poly_thevenard, 100000]]
552        boundary_tags = {'walls': [0,1], 'bom': [2]}
553
554        try:
555            m = create_mesh_from_regions(poly_all,
556                                         boundary_tags,
557                                         10000000,
558                                         interior_regions=interior_regions,
559                                         verbose=False)
560        except:
561            pass
562        else:
563            msg = 'Interior polygon sticking outside bounding polygon should '
564            msg += 'cause an Exception to be raised'
565            raise Exception, msg
566
567    def FIXME_test_create_mesh_with_multiply_tagged_segments(self):
568        '''Test that create_mesh_from_regions fails when
569        segments are listed repeatedly in boundary_tags.
570        '''
571
572        # These are the absolute values
573        min_x = 10
574        min_y = 88
575        polygon = [[min_x,min_y], [1000,100], [1000,1000], [100,1000]]
576        boundary_tags = {'walls': [0,1], 'bom': [1,2]}
577
578        # This one is inside bounding polygon - should pass
579        inner_polygon = [[800,400], [900,500], [800,600]]
580        interior_regions = [(inner_polygon, 5)]
581        m = create_mesh_from_regions(polygon,
582                                     boundary_tags,
583                                     10000000,
584                                     interior_regions=interior_regions,
585                                     verbose=False)
586
587        # This one sticks outside bounding polygon - should fail
588        inner_polygon = [[800,400], [900,500], [800,600]]
589        interior_regions = [(inner_polygon, 5)]
590
591        m = create_mesh_from_regions(polygon,
592                                     boundary_tags,
593                                     10000000,
594                                     interior_regions=interior_regions)
595        try:
596            m = create_mesh_from_regions(polygon,
597                                         boundary_tags,
598                                         10000000,
599                                         interior_regions=interior_regions)
600        except:
601            pass
602        else:
603            msg = 'Tags are listed repeatedly, but create mesh from regions '
604            msg += 'does not cause an Exception to be raised'
605            raise Exception, msg
606
607           
608    def test_create_mesh_with_segments_out_of_bounds(self):
609        """Test that create_mesh_from_regions fails when a segment is out of bounds.
610        """
611       
612        # These are the absolute values
613        min_x = 10 
614        min_y = 88
615        polygon = [[min_x,min_y],[1000,100],[1000,1000],[100,1000]]
616
617       
618        boundary_tags = {'walls': [0,1], 'bom':[2,3], 'out': [5]}
619
620        # This one is inside bounding polygon - should pass
621        inner_polygon = [[800,400],[900,500],[800,600]]
622       
623        interior_regions = [(inner_polygon, 5)]
624
625
626        try:
627            m = create_mesh_from_regions(polygon,
628                                         boundary_tags,
629                                         10000000,
630                                         interior_regions=interior_regions)
631        except:
632            pass
633        else:
634            msg = 'Tags are listed repeatedly, but create mesh from regions '
635            msg += 'does not cause an Exception to be raised'
636            raise Exception, msg
637           
638    def test_create_mesh_with_breaklines(self):
639         # These are the absolute values
640        polygon = [[100,100], [1000,100], [1000,1000], [100,1000]]
641
642        boundary_tags = {'walls': [0,1], 'bom': [2,3]}
643
644        m = create_mesh_from_regions(polygon,
645                                     boundary_tags,
646                                     10000000,
647                                     breaklines=[[[50,50],[2000,2000]]])
648                                                                         
649        self.failUnless(len(m.regions) == 1, 'FAILED!')
650        segs = m.getUserSegments()
651        self.failUnless(len(segs) == 5, 'FAILED!')
652        self.failUnless(len(m.userVertices) == 6, 'FAILED!')
653               
654    def test_create_mesh_from_regions_with_duplicate_verts(self):
655        # These are the absolute values
656        polygon_absolute = [[0.0, 0.0],
657                            [0, 4.0],
658                            [4.0, 4.0],
659                            [4.0, 0.0],
660                            [4.0, 0.0]]
661        x_p = -10
662        y_p = -40
663        zone = 808
664        geo_ref_poly = Geo_reference(zone, x_p, y_p)
665        polygon = geo_ref_poly.change_points_geo_ref(polygon_absolute)
666        boundary_tags = {'50': [0],
667                         '40': [1],
668                         '30': [2],
669                         'no where seg': [3],
670                         '20': [4]}
671        m = create_mesh_from_regions(polygon,
672                                     boundary_tags,
673                                     10000000,
674                                     poly_geo_reference=geo_ref_poly,
675                                     verbose=False)
676
677        fileName = 'badmesh.tsh'
678        #m.export_mesh_file(fileName)
679
680    def concept_create_mesh_from_regions_with_ungenerate(self):
681        x=0
682        y=0
683        mesh_geo = geo_reference=Geo_reference(56, x, y)
684
685        # These are the absolute values
686        polygon_absolute = [[0,0], [100,0], [100,100], [0,100]]
687        x_p = -10
688        y_p = -40
689        geo_ref_poly = Geo_reference(56, x_p, y_p)
690        polygon = geo_ref_poly.change_points_geo_ref(polygon_absolute)
691
692        boundary_tags = {'walls': [0,1], 'bom': [2]}
693
694        inner1_polygon_absolute = [[10,10], [20,10], [20,20], [10,20]]
695        inner1_polygon = geo_ref_poly.\
696                            change_points_geo_ref(inner1_polygon_absolute)
697
698        inner2_polygon_absolute = [[30,30], [40,30], [40,40], [30,40]]
699        inner2_polygon = geo_ref_poly.\
700                            change_points_geo_ref(inner2_polygon_absolute)
701
702        max_area = 10000000
703        interior_regions = [(inner1_polygon, 5), (inner2_polygon, 10)]
704        m = create_mesh_from_regions(polygon,
705                                     boundary_tags,
706                                     max_area,
707                                     interior_regions=interior_regions,
708                                     poly_geo_reference=geo_ref_poly,
709                                     mesh_geo_reference=mesh_geo)
710
711        m.export_mesh_file('a_test_mesh_iknterface.tsh')
712
713        fileName = tempfile.mktemp('.txt')
714        file = open(fileName, 'w')
715        file.write('         1       ??      ??\n\
716       90.0       90.0\n\
717       81.0       90.0\n\
718       81.0       81.0\n\
719       90.0       81.0\n\
720       90.0       90.0\n\
721END\n\
722         2      ?? ??\n\
723       10.0       80.0\n\
724       10.0       90.0\n\
725       20.0       90.0\n\
726       10.0       80.0\n\
727END\n\
728END\n')
729        file.close()
730
731        m.import_ungenerate_file(fileName, tag='wall')
732        os.remove(fileName)
733        m.generate_mesh(maximum_triangle_area=max_area, verbose=False)
734        m.export_mesh_file('b_test_mesh_iknterface.tsh')
735
736    def concept_ungenerateII(self):
737        from anuga.shallow_water import Domain, Reflective_boundary, \
738                            Dirichlet_boundary
739
740        x=0
741        y=0
742        mesh_geo = geo_reference=Geo_reference(56, x, y)
743
744        # These are the absolute values
745        polygon_absolute = [[0,0], [100,0], [100,100], [0,100]]
746        x_p = -10
747        y_p = -40
748        geo_ref_poly = Geo_reference(56, x_p, y_p)
749        polygon = geo_ref_poly.change_points_geo_ref(polygon_absolute)
750
751        boundary_tags = {'wall': [0,1,3], 'wave': [2]}
752
753        inner1_polygon_absolute = [[10,10], [20,10], [20,20], [10,20]]
754        inner1_polygon = geo_ref_poly.\
755                            change_points_geo_ref(inner1_polygon_absolute)
756
757        inner2_polygon_absolute = [[30,30], [40,30], [40,40], [30,40]]
758        inner2_polygon = geo_ref_poly.\
759                            change_points_geo_ref(inner2_polygon_absolute)
760
761        max_area = 1
762        interior_regions = [(inner1_polygon, 5), (inner2_polygon, 10)]
763        m = create_mesh_from_regions(polygon,
764                                     boundary_tags,
765                                     max_area,
766                                     interior_regions=interior_regions,
767                                     poly_geo_reference=geo_ref_poly,
768                                     mesh_geo_reference=mesh_geo)
769
770        m.export_mesh_file('a_test_mesh_iknterface.tsh')
771
772        fileName = tempfile.mktemp('.txt')
773        file = open(fileName, 'w')
774        file.write('         1       ??      ??\n\
775       90.0       90.0\n\
776       81.0       90.0\n\
777       81.0       81.0\n\
778       90.0       81.0\n\
779       90.0       90.0\n\
780END\n\
781         2      ?? ??\n\
782       10.0       80.0\n\
783       10.0       90.0\n\
784       20.0       90.0\n\
785       10.0       80.0\n\
786END\n\
787END\n')
788        file.close()
789
790        m.import_ungenerate_file(fileName)      #, tag='wall')
791        os.remove(fileName)
792        m.generate_mesh(maximum_triangle_area=max_area, verbose=False)
793        mesh_filename = 'bento_b.tsh'
794        m.export_mesh_file(mesh_filename)
795
796        domain = Domain(mesh_filename, use_cache = False)
797
798        Br = Reflective_boundary(domain)
799        Bd = Dirichlet_boundary([3, 0, 0])
800        domain.set_boundary({'wall': Br, 'wave': Bd})
801        yieldstep = 0.1
802        finaltime = 10
803        for t in domain.evolve(yieldstep, finaltime):
804            domain.write_time()
805
806    def concept_ungenerateIII(self):
807        from anuga.shallow_water import Domain, Reflective_boundary, \
808                            Dirichlet_boundary
809        from anuga.pmesh.mesh_interface import create_mesh_from_regions
810
811        # These are the absolute values
812        polygon = [[0,0], [100,0], [100,100], [0,100]]
813
814        boundary_tags = {'wall': [0,1,3], 'wave': [2]}
815        inner1_polygon = [[10,10], [20,10], [20,20], [10,20]]
816        inner2_polygon = [[30,30], [40,30], [40,40], [30,40]]
817
818        max_area = 1
819        interior_regions = [(inner1_polygon, 5), (inner2_polygon, 10)]
820        m = create_mesh_from_regions(polygon,
821                                     boundary_tags,
822                                     max_area,
823                                     interior_regions=interior_regions)
824
825        fileName = tempfile.mktemp('.txt')
826        file = open(fileName, 'w')
827        file.write('         1       ??      ??\n\
828       90.0       90.0\n\
829       81.0       90.0\n\
830       81.0       81.0\n\
831       90.0       81.0\n\
832       90.0       90.0\n\
833END\n\
834         2      ?? ??\n\
835       10.0       80.0\n\
836       10.0       90.0\n\
837       20.0       90.0\n\
838       10.0       80.0\n\
839END\n\
840END\n')
841        file.close()
842
843        m.import_ungenerate_file(fileName)
844        os.remove(fileName)
845        m.generate_mesh(maximum_triangle_area=max_area, verbose=False)
846        mesh_filename = 'mesh.tsh'
847        m.export_mesh_file(mesh_filename)
848
849        domain = Domain(mesh_filename, use_cache=False)
850
851        Br = Reflective_boundary(domain)
852        Bd = Dirichlet_boundary([3, 0, 0])
853        domain.set_boundary({'wall': Br, 'wave': Bd})
854        yieldstep = 0.1
855        finaltime = 10
856        for t in domain.evolve(yieldstep, finaltime):
857            domain.write_time()
858
859    def test_create_mesh_from_regions_check_segs(self):
860        '''Test that create_mesh_from_regions fails
861        when an interior region is outside bounding polygon.
862        '''
863
864        # These are the absolute values
865        min_x = 10
866        min_y = 88
867        polygon = [[min_x,min_y], [1000,100], [1000,1000], [100,1000]]
868        boundary_tags = {'walls': [0,1,3], 'bom': [2]}
869
870        # This one is inside bounding polygon - should pass
871        inner_polygon = [[800,400], [900,500], [800,600]]
872        interior_regions = [(inner_polygon, 5)]
873        m = create_mesh_from_regions(polygon,
874                                     boundary_tags,
875                                     10000000,
876                                     interior_regions=interior_regions)
877
878        boundary_tags = {'walls': [0,1,3,4], 'bom': [2]}
879        try:
880            m = create_mesh_from_regions(polygon,
881                                         boundary_tags,
882                                         10000000,
883                                         interior_regions=interior_regions)
884        except:
885            pass
886        else:
887            msg = 'Segment out of bounds not caught '
888            raise Exception, msg
889
890################################################################################
891
892if __name__ == "__main__":
893    suite = unittest.makeSuite(TestCase,'test')
894    runner = unittest.TextTestRunner() #verbosity=2)
895    runner.run(suite)
896
Note: See TracBrowser for help on using the repository browser.