Context Navigation

test_mesh_interface.py @ 7248

Last change on this file since 7248 was 7207, checked in by rwilson, 16 years ago
Back-merge from Numeric trunk. Final one?
File size: 34.1 KB

Rev	Line
[2276]	1	#!/usr/bin/env python
[6304]	2
[2276]	3	import tempfile
	4	import unittest
[3193]	5	import os
[3514]	6	from anuga.pmesh.mesh import importMeshFromFile
[3715]	7	from anuga.pmesh.mesh_interface import create_mesh_from_regions
[3741]	8
[3742]	9	from anuga.pmesh.mesh_interface import _create_mesh_from_regions
	10
[2276]	11	from load_mesh.loadASCII import *
[3514]	12	from anuga.utilities.polygon import is_inside_polygon
	13	from anuga.coordinate_transforms.geo_reference import Geo_reference,DEFAULT_ZONE
[2276]	14
[6304]	15
[2276]	16	class TestCase(unittest.TestCase):
	17	def setUp(self):
	18	pass
[6304]	19
[2276]	20	def tearDown(self):
	21	pass
	22
	23	def test_create_mesh_from_regions(self):
	24	x=-500
	25	y=-1000
[6304]	26	mesh_geo = geo_reference=Geo_reference(56, x, y)
	27
[2276]	28	# These are the absolute values
[6304]	29	polygon_absolute = [[0,0], [100,0], [100,100], [0,100]]
	30
[2276]	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
[6304]	36	boundary_tags = {'walls': [0,1], 'bom': [2,3]}
[2276]	37
[6304]	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
[2276]	48	m = create_mesh_from_regions(polygon,
	49	boundary_tags,
	50	10000000,
	51	interior_regions=interior_regions,
[2295]	52	poly_geo_reference=geo_ref_poly,
	53	mesh_geo_reference=mesh_geo)
[2276]	54
[2282]	55	# Test the mesh instance
[6304]	56	self.failUnless(len(m.regions)==3, 'FAILED!')
[2282]	57	segs = m.getUserSegments()
[6304]	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!')
[2282]	64
[2295]	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]
[6304]	69
[2295]	70	# poly_point values are relative to the mesh geo-ref
	71	# make them absolute
[6360]	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)))
[6304]	75	self.failUnless(is_inside_polygon([poly_point.x+x, poly_point.y+y],
	76	polygon_absolute, closed=False),
[6360]	77	msg)
[6304]	78
[2295]	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]
[6304]	83
[2295]	84	# poly_point values are relative to the mesh geo-ref
	85	# make them absolute
[6304]	86	self.failUnless(is_inside_polygon([poly_point.x+x, poly_point.y+y],
[3052]	87	inner1_polygon_absolute,
[6304]	88	closed=False),
[2295]	89	'FAILED!')
[6304]	90
[2295]	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]
[6304]	95
[2295]	96	# poly_point values are relative to the mesh geo-ref
	97	# make them absolute
[6304]	98	self.failUnless(is_inside_polygon([poly_point.x+x, poly_point.y+y],
[3052]	99	inner2_polygon_absolute,
[6304]	100	closed=False),
[2295]	101	'FAILED!')
[2402]	102
[3715]	103	def test_create_mesh_from_regions_with_caching(self):
	104	x=-500
	105	y=-1000
[6304]	106	mesh_geo = geo_reference=Geo_reference(56, x, y)
	107
[3715]	108	# These are the absolute values
[6304]	109	polygon_absolute = [[0,0], [100,0], [100,100], [0,100]]
	110
[3715]	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
[6304]	116	boundary_tags = {'walls': [0,1], 'bom': [2,3]}
[3715]	117
[6304]	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)
[3715]	121
[6304]	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
[3719]	128	interior_holes = None
[3715]	129
	130	# Clear cache first
	131	from anuga.caching import cache
[6304]	132
[3715]	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,
[3719]	138	'interior_holes': interior_holes,
[3715]	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
[6304]	156	self.failUnless(len(m.regions)==3, 'FAILED!')
[3715]	157	segs = m.getUserSegments()
[6304]	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!')
[3715]	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]
[6304]	169
[3715]	170	# poly_point values are relative to the mesh geo-ref
	171	# make them absolute
[6304]	172	self.failUnless(is_inside_polygon([poly_point.x+x, poly_point.y+y],
	173	polygon_absolute,
	174	closed=False),
[3715]	175	'FAILED!')
[6304]	176
[3715]	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]
[6304]	181
[3715]	182	# poly_point values are relative to the mesh geo-ref
	183	# make them absolute
[6304]	184	self.failUnless(is_inside_polygon([poly_point.x+x, poly_point.y+y],
[3715]	185	inner1_polygon_absolute,
[6304]	186	closed=False),
[3715]	187	'FAILED!')
[6304]	188
[3715]	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]
[6304]	193
[3715]	194	# poly_point values are relative to the mesh geo-ref
	195	# make them absolute
[6304]	196	self.failUnless(is_inside_polygon([poly_point.x+x, poly_point.y+y],
[3715]	197	inner2_polygon_absolute,
[6304]	198	closed=False),
[3715]	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
[2402]	211	def test_create_mesh_from_regions2(self):
	212	# These are the absolute values
[6304]	213	min_x = -10
[3056]	214	min_y = -88
[6304]	215	polygon_absolute = [[min_x,min_y], [1000,100], [1000,1000], [100,1000]]
	216
[2402]	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
[6304]	223	boundary_tags = {'walls': [0,1], 'bom': [2,3]}
[2402]	224
[6304]	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)]
[2402]	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
[6304]	241	self.failUnless(len(m.regions)==3, 'FAILED!')
[2402]	242	segs = m.getUserSegments()
[6304]	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!')
[2402]	252
	253	def test_create_mesh_from_regions3(self):
	254	# These are the absolute values
[6304]	255	min_x = -10
[3056]	256	min_y = -88
[6304]	257	polygon = [[min_x,min_y], [1000,100], [1000,1000], [100,1000]]
[2402]	258
[6304]	259
[2402]	260	x_p = -10
	261	y_p = -40
	262	geo_ref_poly = Geo_reference(56, x_p, y_p)
	263
[6304]	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)]
[2402]	275	m = create_mesh_from_regions(polygon,
	276	boundary_tags,
	277	10000000,
	278	interior_regions=interior_regions)
	279
	280	# Test the mesh instance
[6304]	281	self.failUnless(len(m.regions) == 3, 'FAILED!')
[2402]	282	segs = m.getUserSegments()
[6304]	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!')
[2402]	292
[3193]	293	def test_create_mesh_from_regions4(self):
[6304]	294	file_name = tempfile.mktemp('.tsh')
[3193]	295
	296	# These are the absolute values
	297	density_outer = 1000
[6304]	298	min_outer = 0
[3193]	299	max_outer = 1000
[6304]	300	polygon_outer = [[min_outer,min_outer], [max_outer,min_outer],
	301	[max_outer,max_outer], [min_outer,max_outer]]
	302
[3193]	303	density_inner1 = 10000000
	304	inner_buffer = 100
	305	min_inner1 = min_outer + inner_buffer
	306	max_inner1 = max_outer - inner_buffer
[6304]	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
[3193]	312	interior_regions = [(inner1_polygon, density_inner1)]
[6304]	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
[3193]	319	m = importMeshFromFile(file_name)
[6304]	320
[4902]	321	self.failUnless(len(m.getTriangulation()) <= 900,
[3193]	322	'Test mesh interface failed!')
[4902]	323	self.failUnless(len(m.getTriangulation()) >= 200,
[3193]	324	'Test mesh interface failed!')
[6304]	325
	326	create_mesh_from_regions(polygon_outer,
	327	boundary_tags,
	328	interior_regions=interior_regions,
	329	filename=file_name,
	330	verbose=False)
	331
[3193]	332	m = importMeshFromFile(file_name)
[6304]	333
[4902]	334	self.failUnless(len(m.getTriangulation()) <= 100,
[3193]	335	'Test mesh interface failed!')
	336
	337	os.remove(file_name)
[6304]	338
[3193]	339	def test_create_mesh_from_regions5(self):
[6304]	340	file_name = tempfile.mktemp('.tsh')
[3193]	341
	342	# These are the absolute values
[6304]	343	density_outer = 10000000
	344	min_outer = 0
[3193]	345	max_outer = 1000
[6304]	346	polygon_outer = [[min_outer,min_outer], [max_outer,min_outer],
	347	[max_outer,max_outer], [min_outer,max_outer]]
	348
[3193]	349	density_inner1 = 1000
	350	inner_buffer = 100
	351	min_inner1 = min_outer + inner_buffer
	352	max_inner1 = max_outer - inner_buffer
[6304]	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
[3193]	358	interior_regions = [(inner1_polygon, density_inner1)]
[6304]	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
[3193]	366	m = importMeshFromFile(file_name)
[6304]	367	self.failUnless(len(m.getTriangulation()) <= 2000,
[3193]	368	'Test mesh interface failed!')
[4902]	369	self.failUnless(len(m.getTriangulation()) >= 900,
[3193]	370	'Test mesh interface failed!')
	371
	372	os.remove(file_name)
[6304]	373
[3195]	374	def test_create_mesh_from_regions6(self):
[6304]	375	file_name = tempfile.mktemp('.tsh')
[3195]	376
	377	# These are the absolute values
	378	density_outer = 1000
[6304]	379	min_outer = 0
[3195]	380	max_outer = 1000
[6304]	381	polygon_outer = [[min_outer,min_outer], [max_outer,min_outer],
	382	[max_outer,max_outer], [min_outer,max_outer]]
[3195]	383
	384	delta = 10
	385	density_inner1 = 1000
	386	min_inner1 = min_outer + delta
	387	max_inner1 = max_outer - delta
[6304]	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
[3195]	392	min_inner2 = min_outer + 2*delta
	393	max_inner2 = max_outer - 2*delta
[6304]	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
[3715]	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)
[6304]	407
[3195]	408	m = importMeshFromFile(file_name)
[6304]	409	self.failUnless(len(m.getTriangulation()) <= 2000,
[3195]	410	'Test mesh interface failed!')
[4902]	411	self.failUnless(len(m.getTriangulation()) >= 900,
[3195]	412	'Test mesh interface failed!')
	413
	414	os.remove(file_name)
[6304]	415
[3195]	416	def test_create_mesh_from_regions7(self):
[6304]	417	file_name = tempfile.mktemp('.tsh')
[3195]	418
	419	# These are the absolute values
[3196]	420	density_outer = 1001
[6304]	421	min_outer = 0
[3195]	422	max_outer = 1000
[6304]	423	polygon_outer = [[min_outer,min_outer], [max_outer,min_outer],
	424	[max_outer,max_outer], [min_outer,max_outer]]
[3195]	425
	426	delta = 10
	427	density_inner1 = 100000000
	428	min_inner1 = min_outer + delta
	429	max_inner1 = max_outer - delta
[6304]	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
[3195]	434	min_inner2 = min_outer + 2*delta
	435	max_inner2 = max_outer - 2*delta
[6304]	436	inner2_polygon = [[min_inner2,min_inner2], [max_inner2,min_inner2],
	437	[max_inner2,max_inner2], [min_inner2,max_inner2]]
[3195]	438
[6304]	439	boundary_tags = {'walls': [0,1], 'bom': [2,3]}
	440
	441	# Note the list order is important
[3195]	442	# The last region added will be the region triangle uses,
	443	# if two regions points are in the same bounded area.
[6304]	444	interior_regions = [(inner2_polygon, density_inner2),
	445	(inner1_polygon, density_inner1)]
[3715]	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)
[6304]	452
[3195]	453	m = importMeshFromFile(file_name)
[6304]	454	self.failUnless(len(m.getTriangulation()) <= 3000,
[3195]	455	'Test mesh interface failed!')
[4902]	456	self.failUnless(len(m.getTriangulation()) >= 2000,
[3195]	457	'Test mesh interface failed!')
	458
	459	os.remove(file_name)
	460
[3056]	461	def test_create_mesh_from_regions_interior_regions(self):
[6304]	462	'''Test that create_mesh_from_regions fails when an interior
	463	region is outside bounding polygon.
	464	'''
[3013]	465
[3038]	466	# These are the absolute values
[6304]	467	min_x = 10
[3038]	468	min_y = 88
[6304]	469	polygon = [[min_x,min_y], [1000,100], [1000,1000], [100,1000]]
	470
	471	boundary_tags = {'walls': [0,1], 'bom': [2,3]}
	472
[3038]	473	# This one is inside bounding polygon - should pass
[6304]	474	inner_polygon = [[800,400], [900,500], [800,600]]
[3868]	475
[3038]	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
[6304]	483	inner_polygon = [[800,400], [900,500], [800,600], [200, 995]]
	484	inner_polygon1 = [[800,400], [1100,500], [800,600]]
[3868]	485	interior_regions = [[inner_polygon, 50], [inner_polygon1, 50]]
[3038]	486
	487	try:
	488	m = create_mesh_from_regions(polygon,
	489	boundary_tags,
	490	10000000,
[3715]	491	interior_regions=interior_regions,
	492	verbose=False)
[3038]	493	except:
	494	pass
	495	else:
	496	msg = 'Interior polygon sticking outside bounding polygon should '
	497	msg += 'cause an Exception to be raised'
[6304]	498	raise Exception, msg
[3038]	499
[3868]	500	def test_create_mesh_from_regions_interior_regions1(self):
[6304]	501	'''Test that create_mesh_from_regions fails
	502	when an interior region is outside bounding polygon.
	503	'''
[3113]	504
[3868]	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]
[6304]	513	poly_all = [d0, d1, d2, d3, d4, d5, d6]
[3868]	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
[6304]	524	# Thevenard Island
[3868]	525	j0 = [294000, 7629000]
	526	j1 = [285000, 7625000]
	527	j2 = [294000, 7621000]
	528	j3 = [299000, 7625000]
	529	poly_thevenard = [j0, j1, j2, j3]
	530
[6304]	531	# med res around onslow
[3868]	532	l0 = [300000, 7610000]
	533	l1 = [285000, 7600000]
	534	l2 = [300000, 7597500]
[6304]	535	l3 = [310000, 7770000] # this one is outside
	536	# l3 = [310000, 7630000] # this one is NOT outside
[3868]	537	l4 = [315000, 7610000]
	538	poly_coast = [l0, l1, l2, l3, l4]
	539
[6304]	540	# general coast and local area to onslow region
[3868]	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
[6304]	550	interior_regions = [[poly_onslow, 50000], [poly_region, 50000],
	551	[poly_coast, 100000], [poly_thevenard, 100000]]
	552	boundary_tags = {'walls': [0,1], 'bom': [2]}
[3868]	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'
[6304]	565	raise Exception, msg
[3868]	566
[3126]	567	def FIXME_test_create_mesh_with_multiply_tagged_segments(self):
[6304]	568	'''Test that create_mesh_from_regions fails when
[3113]	569	segments are listed repeatedly in boundary_tags.
[6304]	570	'''
[3038]	571
[3113]	572	# These are the absolute values
[6304]	573	min_x = 10
[3113]	574	min_y = 88
[6304]	575	polygon = [[min_x,min_y], [1000,100], [1000,1000], [100,1000]]
	576	boundary_tags = {'walls': [0,1], 'bom': [1,2]}
[3038]	577
[3113]	578	# This one is inside bounding polygon - should pass
[6304]	579	inner_polygon = [[800,400], [900,500], [800,600]]
[3113]	580	interior_regions = [(inner_polygon, 5)]
	581	m = create_mesh_from_regions(polygon,
	582	boundary_tags,
	583	10000000,
[6304]	584	interior_regions=interior_regions,
	585	verbose=False)
[3113]	586
	587	# This one sticks outside bounding polygon - should fail
[6304]	588	inner_polygon = [[800,400], [900,500], [800,600]]
[3113]	589	interior_regions = [(inner_polygon, 5)]
[7035]	590
	591	m = create_mesh_from_regions(polygon,
	592	boundary_tags,
	593	10000000,
	594	interior_regions=interior_regions)
[3113]	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'
[6304]	605	raise Exception, msg
[3113]	606
[7035]	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
	639
	640
[3013]	641	def test_create_mesh_from_regions_with_duplicate_verts(self):
	642	# These are the absolute values
[6304]	643	polygon_absolute = [[0.0, 0.0],
	644	[0, 4.0],
	645	[4.0, 4.0],
	646	[4.0, 0.0],
	647	[4.0, 0.0]]
[3013]	648	x_p = -10
	649	y_p = -40
	650	zone = 808
	651	geo_ref_poly = Geo_reference(zone, x_p, y_p)
	652	polygon = geo_ref_poly.change_points_geo_ref(polygon_absolute)
[6304]	653	boundary_tags = {'50': [0],
	654	'40': [1],
	655	'30': [2],
	656	'no where seg': [3],
	657	'20': [4]}
[3013]	658	m = create_mesh_from_regions(polygon,
	659	boundary_tags,
	660	10000000,
[6304]	661	poly_geo_reference=geo_ref_poly,
	662	verbose=False)
[3013]	663
[6304]	664	fileName = 'badmesh.tsh'
[3013]	665	#m.export_mesh_file(fileName)
[6304]	666
[5193]	667	def concept_create_mesh_from_regions_with_ungenerate(self):
	668	x=0
	669	y=0
[6304]	670	mesh_geo = geo_reference=Geo_reference(56, x, y)
	671
[5193]	672	# These are the absolute values
[6304]	673	polygon_absolute = [[0,0], [100,0], [100,100], [0,100]]
[5193]	674	x_p = -10
	675	y_p = -40
	676	geo_ref_poly = Geo_reference(56, x_p, y_p)
	677	polygon = geo_ref_poly.change_points_geo_ref(polygon_absolute)
	678
[6304]	679	boundary_tags = {'walls': [0,1], 'bom': [2]}
[5193]	680
[6304]	681	inner1_polygon_absolute = [[10,10], [20,10], [20,20], [10,20]]
	682	inner1_polygon = geo_ref_poly.\
	683	change_points_geo_ref(inner1_polygon_absolute)
	684
	685	inner2_polygon_absolute = [[30,30], [40,30], [40,40], [30,40]]
	686	inner2_polygon = geo_ref_poly.\
	687	change_points_geo_ref(inner2_polygon_absolute)
	688
[5193]	689	max_area = 10000000
[6304]	690	interior_regions = [(inner1_polygon, 5), (inner2_polygon, 10)]
[5193]	691	m = create_mesh_from_regions(polygon,
	692	boundary_tags,
	693	max_area,
	694	interior_regions=interior_regions,
	695	poly_geo_reference=geo_ref_poly,
	696	mesh_geo_reference=mesh_geo)
[6304]	697
	698	m.export_mesh_file('a_test_mesh_iknterface.tsh')
	699
	700	fileName = tempfile.mktemp('.txt')
	701	file = open(fileName, 'w')
	702	file.write(' 1 ?? ??\n\
[5193]	703	90.0 90.0\n\
	704	81.0 90.0\n\
	705	81.0 81.0\n\
	706	90.0 81.0\n\
	707	90.0 90.0\n\
	708	END\n\
	709	2 ?? ??\n\
	710	10.0 80.0\n\
	711	10.0 90.0\n\
	712	20.0 90.0\n\
	713	10.0 80.0\n\
	714	END\n\
[6304]	715	END\n')
	716	file.close()
	717
[5193]	718	m.import_ungenerate_file(fileName, tag='wall')
	719	os.remove(fileName)
[6304]	720	m.generate_mesh(maximum_triangle_area=max_area, verbose=False)
[5193]	721	m.export_mesh_file('b_test_mesh_iknterface.tsh')
[5638]	722
	723	def concept_ungenerateII(self):
	724	from anuga.shallow_water import Domain, Reflective_boundary, \
	725	Dirichlet_boundary
[6304]	726
[5638]	727	x=0
	728	y=0
[6304]	729	mesh_geo = geo_reference=Geo_reference(56, x, y)
	730
[5638]	731	# These are the absolute values
[6304]	732	polygon_absolute = [[0,0], [100,0], [100,100], [0,100]]
[5638]	733	x_p = -10
	734	y_p = -40
	735	geo_ref_poly = Geo_reference(56, x_p, y_p)
	736	polygon = geo_ref_poly.change_points_geo_ref(polygon_absolute)
	737
[6304]	738	boundary_tags = {'wall': [0,1,3], 'wave': [2]}
[5638]	739
[6304]	740	inner1_polygon_absolute = [[10,10], [20,10], [20,20], [10,20]]
	741	inner1_polygon = geo_ref_poly.\
	742	change_points_geo_ref(inner1_polygon_absolute)
	743
	744	inner2_polygon_absolute = [[30,30], [40,30], [40,40], [30,40]]
	745	inner2_polygon = geo_ref_poly.\
	746	change_points_geo_ref(inner2_polygon_absolute)
	747
[5638]	748	max_area = 1
[6304]	749	interior_regions = [(inner1_polygon, 5), (inner2_polygon, 10)]
[5638]	750	m = create_mesh_from_regions(polygon,
	751	boundary_tags,
	752	max_area,
	753	interior_regions=interior_regions,
	754	poly_geo_reference=geo_ref_poly,
	755	mesh_geo_reference=mesh_geo)
[6304]	756
	757	m.export_mesh_file('a_test_mesh_iknterface.tsh')
	758
	759	fileName = tempfile.mktemp('.txt')
	760	file = open(fileName, 'w')
	761	file.write(' 1 ?? ??\n\
[5638]	762	90.0 90.0\n\
	763	81.0 90.0\n\
	764	81.0 81.0\n\
	765	90.0 81.0\n\
	766	90.0 90.0\n\
	767	END\n\
	768	2 ?? ??\n\
	769	10.0 80.0\n\
	770	10.0 90.0\n\
	771	20.0 90.0\n\
	772	10.0 80.0\n\
	773	END\n\
[6304]	774	END\n')
	775	file.close()
	776
	777	m.import_ungenerate_file(fileName) #, tag='wall')
[5638]	778	os.remove(fileName)
[6304]	779	m.generate_mesh(maximum_triangle_area=max_area, verbose=False)
	780	mesh_filename = 'bento_b.tsh'
[5638]	781	m.export_mesh_file(mesh_filename)
	782
	783	domain = Domain(mesh_filename, use_cache = False)
[6304]	784
[5638]	785	Br = Reflective_boundary(domain)
[6304]	786	Bd = Dirichlet_boundary([3, 0, 0])
	787	domain.set_boundary({'wall': Br, 'wave': Bd})
[5638]	788	yieldstep = 0.1
	789	finaltime = 10
[6304]	790	for t in domain.evolve(yieldstep, finaltime):
[5638]	791	domain.write_time()
[6304]	792
[5638]	793	def concept_ungenerateIII(self):
	794	from anuga.shallow_water import Domain, Reflective_boundary, \
	795	Dirichlet_boundary
	796	from anuga.pmesh.mesh_interface import create_mesh_from_regions
[6304]	797
[5638]	798	# These are the absolute values
[6304]	799	polygon = [[0,0], [100,0], [100,100], [0,100]]
[5638]	800
[6304]	801	boundary_tags = {'wall': [0,1,3], 'wave': [2]}
	802	inner1_polygon = [[10,10], [20,10], [20,20], [10,20]]
	803	inner2_polygon = [[30,30], [40,30], [40,40], [30,40]]
	804
[5638]	805	max_area = 1
[6304]	806	interior_regions = [(inner1_polygon, 5), (inner2_polygon, 10)]
[5638]	807	m = create_mesh_from_regions(polygon,
	808	boundary_tags,
	809	max_area,
	810	interior_regions=interior_regions)
[6304]	811
	812	fileName = tempfile.mktemp('.txt')
	813	file = open(fileName, 'w')
	814	file.write(' 1 ?? ??\n\
[5638]	815	90.0 90.0\n\
	816	81.0 90.0\n\
	817	81.0 81.0\n\
	818	90.0 81.0\n\
	819	90.0 90.0\n\
	820	END\n\
	821	2 ?? ??\n\
	822	10.0 80.0\n\
	823	10.0 90.0\n\
	824	20.0 90.0\n\
	825	10.0 80.0\n\
	826	END\n\
[6304]	827	END\n')
	828	file.close()
	829
	830	m.import_ungenerate_file(fileName)
[5638]	831	os.remove(fileName)
[6304]	832	m.generate_mesh(maximum_triangle_area=max_area, verbose=False)
	833	mesh_filename = 'mesh.tsh'
[5638]	834	m.export_mesh_file(mesh_filename)
	835
[6304]	836	domain = Domain(mesh_filename, use_cache=False)
	837
[5638]	838	Br = Reflective_boundary(domain)
[6304]	839	Bd = Dirichlet_boundary([3, 0, 0])
	840	domain.set_boundary({'wall': Br, 'wave': Bd})
[5638]	841	yieldstep = 0.1
	842	finaltime = 10
[6304]	843	for t in domain.evolve(yieldstep, finaltime):
[5638]	844	domain.write_time()
[5717]	845
	846	def test_create_mesh_from_regions_check_segs(self):
[6304]	847	'''Test that create_mesh_from_regions fails
	848	when an interior region is outside bounding polygon.
	849	'''
[5717]	850
	851	# These are the absolute values
[6304]	852	min_x = 10
[5717]	853	min_y = 88
[6304]	854	polygon = [[min_x,min_y], [1000,100], [1000,1000], [100,1000]]
	855	boundary_tags = {'walls': [0,1,3], 'bom': [2]}
	856
[5717]	857	# This one is inside bounding polygon - should pass
[6304]	858	inner_polygon = [[800,400], [900,500], [800,600]]
[5717]	859	interior_regions = [(inner_polygon, 5)]
	860	m = create_mesh_from_regions(polygon,
	861	boundary_tags,
	862	10000000,
	863	interior_regions=interior_regions)
	864
[6304]	865	boundary_tags = {'walls': [0,1,3,4], 'bom': [2]}
[5717]	866	try:
	867	m = create_mesh_from_regions(polygon,
	868	boundary_tags,
	869	10000000,
	870	interior_regions=interior_regions)
	871	except:
	872	pass
	873	else:
[6304]	874	msg = 'Segment out of bounds not caught '
	875	raise Exception, msg
[5717]	876
[6360]	877	################################################################################
[6304]	878
[2276]	879	if __name__ == "__main__":
	880	suite = unittest.makeSuite(TestCase,'test')
[3056]	881	runner = unittest.TextTestRunner() #verbosity=2)
[2276]	882	runner.run(suite)
[6304]	883

Note: See TracBrowser for help on using the repository browser.

Context Navigation

source: branches/numpy/anuga/pmesh/test_mesh_interface.py @ 7248

Download in other formats: