source: inundation-numpy-branch/debug/mesh_error_reporting/show_bad_mesh.py @ 3514

Last change on this file since 3514 was 3514, checked in by duncan, 18 years ago

Hi all,
I'm doing a change in the anuga structure, moving the code to

\anuga_core\source\anuga

After you have done an svn update, the PYTHONPATH has to be changed to;
PYTHONPATH = anuga_core/source/

This is part of changes required to make installation of anuga quicker and reducing the size of our sandpits.

If any imports are broken, try fixing them. With adding anuga. to them for example. If this seems to have really broken things, email/phone me.

Cheers
Duncan

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 anuga.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.auto_segment(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.auto_segment(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.