source: inundation/debug/mesh_error_reporting/show_bad_mesh.py @ 2378

Last change on this file since 2378 was 1954, checked in by duncan, 19 years ago

show generation of a bad mesh.

File size: 2.9 KB
Line 
1"""
2This will produce a mesh which has vertex 10 that has no triangle.
3
4Why?
5   
6   Ole Nielsen, Stephen Roberts, Duncan Gray, Christopher Zoppou
7   Geoscience Australia, 2004.
8"""
9
10
11import unittest
12from math import sqrt
13import time
14
15
16from pyvolution.least_squares import Interpolation
17from Numeric import allclose, array, transpose
18
19from coordinate_transforms.geo_reference import Geo_reference
20
21def distance(x, y):
22    return sqrt( sum( (array(x)-array(y))**2 ))
23
24def linear_function(point):
25    point = array(point)
26    return point[:,0]+point[:,1]
27
28
29class Test_Least_Squares(unittest.TestCase):
30
31    def setUp(self):
32        pass
33
34    def tearDown(self):
35        pass
36
37    def test_expand_search3(self):
38        from random import seed, random
39        from pmesh.mesh import Mesh
40
41        num_of_points = 100
42       
43        seed(2)
44        m = Mesh()
45        m.addUserVertex(0,0)
46        m.addUserVertex(100,0)
47        m.addUserVertex(0,100)
48        m.addUserVertex(100,100)
49       
50        m.autoSegment(alpha = 100 )
51       
52        dict = {}
53        dict['points'] = [[10,10],[90,20]]
54        dict['segments'] = [[0,1]] 
55        dict['segment_tags'] = ['wall1']   
56        m.addVertsSegs(dict)
57   
58        dict = {}
59        dict['points'] = [[10,90],[40,20]]
60        dict['segments'] = [[0,1]] 
61        dict['segment_tags'] = ['wall2']   
62        m.addVertsSegs(dict)
63       
64        dict = {}
65        dict['points'] = [[20,90],[60,60]]
66        dict['segments'] = [[0,1]] 
67        dict['segment_tags'] = ['wall3'] 
68        m.addVertsSegs(dict)
69       
70        dict = {}
71        dict['points'] = [[60,20],[90,90]]
72        dict['segments'] = [[0,1]] 
73        dict['segment_tags'] = ['wall4']   
74        m.addVertsSegs(dict)
75       
76        m.addVertsSegs(dict)
77        m.generateMesh(mode = "Q")
78        m.export_mesh_file("aaaa.tsh")
79        mesh_dict =  m.Mesh2IOTriangulationDict()
80
81        points = []
82        point_atts = []
83        for point in range(num_of_points):
84            points.append([random()*100, random()*100])
85            point_atts.append(10.0)
86        #m.autoSegment(alpha = 100 )
87
88        #print 'points', points
89        t0 = time.time()
90        interp = Interpolation(mesh_dict['vertices'],
91                               mesh_dict['triangles'], points,
92                               alpha=0.2, expand_search=True,
93                               verbose = True,
94                               max_points_per_cell = 4)
95        calc = interp.fit_points(point_atts )
96        print "interp.expanded_quad_searches", interp.expanded_quad_searches
97        #print "calc", calc
98        print 'That took %.2f seconds' %(time.time()-t0)
99       
100#-------------------------------------------------------------
101if __name__ == "__main__":
102    suite = unittest.makeSuite(Test_Least_Squares,'test')
103    runner = unittest.TextTestRunner(verbosity=1)
104    runner.run(suite)
Note: See TracBrowser for help on using the repository browser.