| 1 | """Create mesh for lwru2 validation of the Oshiri island tsunami |
|---|
| 2 | """ |
|---|
| 3 | |
|---|
| 4 | #!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
|---|
| 5 | # Assume that the root AnuGA dir is included in your pythonpath |
|---|
| 6 | #!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! |
|---|
| 7 | |
|---|
| 8 | from pmesh.mesh import * |
|---|
| 9 | from pyvolution.coordinate_transforms.geo_reference import Geo_reference |
|---|
| 10 | |
|---|
| 11 | #------------------------------------------------------------- |
|---|
| 12 | if __name__ == "__main__": |
|---|
| 13 | |
|---|
| 14 | |
|---|
| 15 | #Basic geometry |
|---|
| 16 | |
|---|
| 17 | xleft = 0 |
|---|
| 18 | xright = 5.448 |
|---|
| 19 | ybottom = 0 |
|---|
| 20 | ytop = 3.402 |
|---|
| 21 | |
|---|
| 22 | #Outline |
|---|
| 23 | point_sw = [xleft, ybottom] |
|---|
| 24 | point_se = [xright, ybottom] |
|---|
| 25 | point_nw = [xleft, ytop] |
|---|
| 26 | point_ne = [xright, ytop] |
|---|
| 27 | |
|---|
| 28 | #Midway points (left) |
|---|
| 29 | point_ml1 = [xleft + (xright-xleft)/3+1, ytop] |
|---|
| 30 | point_ml2 = [xleft + (xright-xleft)/3+1, ybottom] |
|---|
| 31 | |
|---|
| 32 | #Midway points (right) |
|---|
| 33 | point_mr1 = [xleft + 2*(xright-xleft)/3+0.5, ytop] |
|---|
| 34 | point_mr2 = [xleft + 2*(xright-xleft)/3+0.5, ybottom] |
|---|
| 35 | |
|---|
| 36 | |
|---|
| 37 | geo = Geo_reference(xllcorner = xleft, |
|---|
| 38 | yllcorner = ybottom) |
|---|
| 39 | |
|---|
| 40 | |
|---|
| 41 | print "***********************" |
|---|
| 42 | print "geo ref", geo |
|---|
| 43 | print "***********************" |
|---|
| 44 | |
|---|
| 45 | m = Mesh(geo_reference=geo) |
|---|
| 46 | |
|---|
| 47 | #Boundary |
|---|
| 48 | dict = {} |
|---|
| 49 | dict['points'] = [point_se, #se |
|---|
| 50 | point_ne, |
|---|
| 51 | point_mr1, |
|---|
| 52 | point_ml1, |
|---|
| 53 | point_nw, |
|---|
| 54 | point_sw, |
|---|
| 55 | point_ml2, |
|---|
| 56 | point_mr2] |
|---|
| 57 | |
|---|
| 58 | |
|---|
| 59 | dict['segments'] = [[0,1], [1,2], [2,3], [3,4], |
|---|
| 60 | [4,5], [5,6], [6,7], [7,0], #The outer border |
|---|
| 61 | [2,7], [3,6]] #Separators |
|---|
| 62 | |
|---|
| 63 | dict['segment_tags'] = ['wall', |
|---|
| 64 | 'wall', |
|---|
| 65 | 'wall', |
|---|
| 66 | 'wall', |
|---|
| 67 | 'wave', |
|---|
| 68 | 'wall', |
|---|
| 69 | 'wall', |
|---|
| 70 | 'wall', |
|---|
| 71 | '', #Interior |
|---|
| 72 | ''] #Interior |
|---|
| 73 | |
|---|
| 74 | |
|---|
| 75 | m.addVertsSegs(dict) |
|---|
| 76 | |
|---|
| 77 | |
|---|
| 78 | |
|---|
| 79 | #Localised refined area for gulleys |
|---|
| 80 | dict = {} |
|---|
| 81 | xl = 4.8 |
|---|
| 82 | xr = 5.3 |
|---|
| 83 | yb = 1.6 |
|---|
| 84 | yt = 2.3 |
|---|
| 85 | p0 = [xl, yb] |
|---|
| 86 | p1 = [xr, yb] |
|---|
| 87 | p2 = [xr, yt] |
|---|
| 88 | p3 = [xl, yt] |
|---|
| 89 | |
|---|
| 90 | dict['points'] = [p0, p1, p2, p3] |
|---|
| 91 | |
|---|
| 92 | |
|---|
| 93 | dict['segments'] = [[0,1], [1,2], [2,3], [3,0]] |
|---|
| 94 | dict['segment_tags'] = ['', '', '', ''] |
|---|
| 95 | m.addVertsSegs(dict) |
|---|
| 96 | |
|---|
| 97 | |
|---|
| 98 | base_resolution = 100 |
|---|
| 99 | |
|---|
| 100 | ocean = m.addRegionEN(xleft+1, ybottom+1) |
|---|
| 101 | ocean.setMaxArea(0.01*base_resolution) |
|---|
| 102 | |
|---|
| 103 | mid = m.addRegionEN(point_ml2[0]+1, ybottom+1) |
|---|
| 104 | mid.setMaxArea(0.001*base_resolution) |
|---|
| 105 | |
|---|
| 106 | |
|---|
| 107 | inner = m.addRegionEN(point_mr2[0]+1, ybottom+1) |
|---|
| 108 | inner.setMaxArea(0.0001*base_resolution) |
|---|
| 109 | |
|---|
| 110 | |
|---|
| 111 | gulleys = m.addRegionEN(p0[0]+0.1, p0[1]+0.1) |
|---|
| 112 | gulleys.setMaxArea(0.00001*base_resolution) |
|---|
| 113 | |
|---|
| 114 | |
|---|
| 115 | m.generateMesh('pzq28.0za1000000a') |
|---|
| 116 | |
|---|
| 117 | import filenames |
|---|
| 118 | m.export_mesh_file(filenames.mesh_filename) |
|---|
