| 1 | |
|---|
| 2 | import Numeric |
|---|
| 3 | import math |
|---|
| 4 | import random |
|---|
| 5 | # add inundation dir to your pythonpath |
|---|
| 6 | from anuga.pmesh.mesh import Mesh |
|---|
| 7 | from anuga.coordinate_transforms.geo_reference import Geo_reference |
|---|
| 8 | |
|---|
| 9 | WidtH = 200 # width of boudary in metres |
|---|
| 10 | |
|---|
| 11 | #random.uniform(0,1) #proper implemantation of random generator |
|---|
| 12 | |
|---|
| 13 | def create_mesh(maximum_triangle_area, depth, |
|---|
| 14 | mesh_file=None, |
|---|
| 15 | triangles_in_name = False): |
|---|
| 16 | """ |
|---|
| 17 | triangles_in_name, if True is used to append the number of |
|---|
| 18 | triangles in the mesh to the mesh file name. |
|---|
| 19 | """ |
|---|
| 20 | |
|---|
| 21 | WidtH = 200 # width of boudary in metres |
|---|
| 22 | breadth = depth |
|---|
| 23 | |
|---|
| 24 | print "building footprint" |
|---|
| 25 | print depth * breadth , "m^2" |
|---|
| 26 | block = 625 |
|---|
| 27 | BL = block**0.5 |
|---|
| 28 | |
|---|
| 29 | porosity = breadth/BL |
|---|
| 30 | print porosity, " Building porosity" |
|---|
| 31 | # create a mesh instance of class Mesh |
|---|
| 32 | m = Mesh() |
|---|
| 33 | |
|---|
| 34 | # Boundary of problem |
|---|
| 35 | outer_polygon = [[0,0],[5*WidtH,0],[5*WidtH,WidtH],[0,WidtH]] |
|---|
| 36 | m.add_region_from_polygon(outer_polygon, tags={'wall':[0,2], 'wave':[3], 'back':[1]}) |
|---|
| 37 | |
|---|
| 38 | # inner polygons => building boundaries |
|---|
| 39 | |
|---|
| 40 | whs_front = 15/2 |
|---|
| 41 | lhs_front = 15/2 |
|---|
| 42 | whs_rear = 5/2 |
|---|
| 43 | lhs_rear = 5/2 |
|---|
| 44 | #Th = (45 *(3.14159/180)) # sets an initial rotation |
|---|
| 45 | Th = 0 |
|---|
| 46 | ForDep = (0.2*WidtH) + (BL-whs_front) |
|---|
| 47 | RearDep = 1.2*WidtH |
|---|
| 48 | |
|---|
| 49 | Breadths = Numeric.arrayrange( -(BL/2), WidtH+(BL/2), (BL)) |
|---|
| 50 | #Breadths = Numeric.arrayrange( (BL/2), 5*WidtH-(BL/2), (BL)) # For ortho-offset |
|---|
| 51 | #print Breadths, "Breadths" |
|---|
| 52 | Depths_first = Numeric.arrayrange( ForDep, RearDep, BL ) |
|---|
| 53 | Depths_second = Numeric.arrayrange( RearDep+BL, 2.2*WidtH, BL ) |
|---|
| 54 | #print Depths, "Depths" |
|---|
| 55 | |
|---|
| 56 | for i,D in enumerate(Depths_first): |
|---|
| 57 | Breadths = Breadths + ((-1)**i)*(0.5*BL/2) #Used to offset buildings |
|---|
| 58 | # First block of buildings |
|---|
| 59 | for B in Breadths: |
|---|
| 60 | wh1 = (-whs_front) * math.cos(Th) + (-lhs_front) * math.sin(Th) |
|---|
| 61 | lh1 = (-lhs_front) * math.cos(Th) - (-whs_front) * math.sin(Th) |
|---|
| 62 | wh2 = (+whs_front) * math.cos(Th) + (-lhs_front) * math.sin(Th) |
|---|
| 63 | lh2 = (-lhs_front) * math.cos(Th) - (+whs_front) * math.sin(Th) |
|---|
| 64 | wh3 = (+whs_front) * math.cos(Th) + (+lhs_front) * math.sin(Th) |
|---|
| 65 | lh3 = (+lhs_front) * math.cos(Th) - (+whs_front) * math.sin(Th) |
|---|
| 66 | wh4 = (-whs_front) * math.cos(Th) + (+lhs_front) * math.sin(Th) |
|---|
| 67 | lh4 = (+lhs_front) * math.cos(Th) - (-whs_front) * math.sin(Th) |
|---|
| 68 | polygon = [[D+wh1,B+lh1],[D+wh2,B+lh2],[D+wh3,B+lh3],[D+wh4,B+lh4]] |
|---|
| 69 | m.add_hole_from_polygon(polygon, tags={'wall':[0,1,2,3]}) |
|---|
| 70 | # Th = Th + (37.3 *(3.14159/180)) # keeps rotating individual buildings. |
|---|
| 71 | |
|---|
| 72 | for i,D in enumerate(Depths_second): |
|---|
| 73 | Breadths = Breadths + ((-1)**i)*(0.5*BL/2) #Used to offset buildings |
|---|
| 74 | # Second block of buildings |
|---|
| 75 | for B in Breadths: |
|---|
| 76 | wh1 = (-whs_rear) * math.cos(Th) + (-lhs_rear) * math.sin(Th) |
|---|
| 77 | lh1 = (-lhs_rear) * math.cos(Th) - (-whs_rear) * math.sin(Th) |
|---|
| 78 | wh2 = (+whs_rear) * math.cos(Th) + (-lhs_rear) * math.sin(Th) |
|---|
| 79 | lh2 = (-lhs_rear) * math.cos(Th) - (+whs_rear) * math.sin(Th) |
|---|
| 80 | wh3 = (+whs_rear) * math.cos(Th) + (+lhs_rear) * math.sin(Th) |
|---|
| 81 | lh3 = (+lhs_rear) * math.cos(Th) - (+whs_rear) * math.sin(Th) |
|---|
| 82 | wh4 = (-whs_rear) * math.cos(Th) + (+lhs_rear) * math.sin(Th) |
|---|
| 83 | lh4 = (+lhs_rear) * math.cos(Th) - (-whs_rear) * math.sin(Th) |
|---|
| 84 | polygon = [[D+wh1,B+lh1],[D+wh2,B+lh2],[D+wh3,B+lh3],[D+wh4,B+lh4]] |
|---|
| 85 | m.add_hole_from_polygon(polygon, tags={'wall':[0,1,2,3]}) |
|---|
| 86 | # measures apparent depth to inundating wave |
|---|
| 87 | # App_breadth = depth*(math.sin(math.pi()*Th/180) + math.cos(math.pi()*Th/180))/BL |
|---|
| 88 | |
|---|
| 89 | |
|---|
| 90 | m.generate_mesh(maximum_triangle_area=maximum_triangle_area) |
|---|
| 91 | triangle_count = m.get_triangle_count() |
|---|
| 92 | |
|---|
| 93 | if mesh_file is None: |
|---|
| 94 | return m, triangle_count |
|---|
| 95 | else: |
|---|
| 96 | if triangles_in_name is True: |
|---|
| 97 | mesh_file = mesh_file[:-4] + '_Str_Off(hlf)_D=' + str(whs_front*2) + '_' + str(whs_rear*2) + '_' + str(triangle_count) \ |
|---|
| 98 | + mesh_file[-4:] |
|---|
| 99 | m.export_mesh_file(mesh_file) |
|---|
| 100 | return mesh_file, triangle_count |
|---|
| 101 | |
|---|
| 102 | #------------------------------------------------------------- |
|---|
| 103 | if __name__ == "__main__": |
|---|
| 104 | _, triangle_count = create_mesh(100,15,mesh_file="test.tsh") |
|---|
| 105 | print "triangle_count",triangle_count |
|---|
