Changeset 5031

Timestamp:

Feb 14, 2008, 1:50:57 PM (15 years ago)

Author:

nick

Message:

update circular island

File:

• anuga_validation/conical_island/run_circular.py (modified) (18 diffs)

anuga_validation/conical_island/run_circular.py

@@ -10 +10 @@
 from anuga.shallow_water import Transmissive_Momentum_Set_Stage_boundary
 from anuga.shallow_water import Dirichlet_boundary, Time_boundary
-from anuga.shallow_water.data_manager import get_maximum_inundation_data
+from anuga.shallow_water.data_manager import get_maximum_inundation_data, start_screen_catcher
 from anuga.abstract_2d_finite_volumes.util import file_function, sww2csv_gauges
 from anuga.pmesh.mesh_interface import create_mesh_from_regions
@@ -23 +23 @@
 #-------------------------
 
+#start_screen_catcher()
+
 def get_xy(x=0,y=0,r=1,angle=1):
+    #return x and y co-ordinates 
     x1 = x + cos(angle)*r
     y1 = y + sin(-angle)*r
@@ -91 +94 @@
           247.5, 270.0, 292.5, 315.0, 337.5)
 
-x = 12.96
-y = 13.80
+#x = 12.96
+#y = 13.80
+center_x = 12.96
+center_y = 13.80
+
+#wavelenght = 
+
 r1 = 1.5
 r2 = 2.5
@@ -98 +106 @@
 res=.05
 #res=1.
+L=4.11
+
 
 #create polygons (circles) to have higher resolution
@@ -107 +117 @@
     angle = ((angle-180)/180.0)*pi
 
-    p = get_xy(x,y,3.6,angle)
+    p = get_xy(center_x,center_y,3.6,angle)
     poly1.append([p[0],p[1]])
 
@@ -117 +127 @@
     #convert Degs to Rads
     angle = ((angle-180)/180.0)*pi
-    p = get_xy(x,y,2.4,angle)
+    p = get_xy(center_x,center_y,2.4,angle)
     poly2.append([p[0],p[1]])
     
@@ -127 +137 @@
     #convert Degs to Rads
     angle = ((angle-180)/180.0)*pi
-    p = get_xy(x,y,1.6,angle)
+    p = get_xy(center_x,center_y,1.6,angle)
     poly3.append([p[0],p[1]])
     
@@ -137 +147 @@
     #convert Degs to Rads
     angle = ((angle-180)/180.0)*pi
-    p = get_xy(x,y,1.1,angle)
+    p = get_xy(center_x,center_y,1.1,angle)
     poly4.append([p[0],p[1]])
     
@@ -145 +155 @@
 #plot_polygons(poly,'poly.png')
 
-poly_all= [[0,7],[0,25],[30,25],[30,7]]
+#poly_all= [[0,7],[0,25],[30,25],[30,7]]
+poly_all= [[0,L],[0,25],[30,25],[30,L]]
 create_mesh_from_regions(poly_all,
                          boundary_tags={'wall': [0,2],'wave': [3], 'buffer':[1]},
@@ -164 +175 @@
 domain.set_quantity('stage', water_height)
 
-#attribute_dic, title_index_dic = csv2dict('sqrt_table.csv')
-attribute_dic, title_index_dic = csv2dict('sqrt_table_new.csv')
+
+attribute_dic, title_index_dic = csv2dict('sqrt_table.csv')
+#print attribute_dic
 print attribute_dic['number'][10]
 
 def test(x,y):
-    center_x = 12.96
+#    center_x = 12.96
 #    center_y = 13.80
-    center_y = 6.80
-
-    du = (center_x-x)**2+(center_y-y)**2
+#    center_y = 6.80
+
+    du = (center_x-x)**2+(center_y-L-y)**2
 #    print 'x',min(x),max(x)
 #    print 'y',min(y),max(y)
@@ -182 +194 @@
         zz=-float(attribute_dic['sqrt'][nn])
 #        print nn,zz
-        # this is the elevation for a cone with 1/4 slope and is ?? below the water .
+
+#    z = -square_root((center_x-x)**2+(center_y-y)**2)
         z = zz*.25+0.8975
 
@@ -189 +202 @@
         if z > .625:
             z=0.625
-        print 'hello',d,nn,zz,z
-#        print 'hello',d,n1,n2,n3,nn,zz,z
+#        print 'hello',z
         z1.append(z)
     return z1
 
-domain.set_quantity('elevation',test, alpha=1., verbose=True)
+domain.set_quantity('elevation',test, verbose=True)
 
 
@@ -200 +212 @@
 # Set simulation parameters
 #-------------------------
-domain.set_name('test1')  # Name of output sww file
+model_name='wavetank_model'
+domain.set_name(model_name)  # Name of output sww file
 domain.set_default_order(1)               # Apply second order scheme 
 #domain.set_all_limiters(0.9)              # Max second order scheme (old lim)
@@ -208 +221 @@
 domain.beta_h = 0.0
 
-#Timings on AMD64-242 (beta_h=0)
-#  tight_slope_limiters = 0:
-#    3035s - 3110s
-#  tight_slope_limiters = 1:
-#    3000s - 3008s
-#
-# beta_h>0: In the order of 3200s
 
 #-------------------------
@@ -223 +229 @@
 #function = file_function(project.boundary_filename,
 #                         domain, verbose=True)
-def waveform(t): 
-        return (sin(t*2*pi))
 
 # Create and assign boundary objects
@@ -236 +240 @@
 Bts = Transmissive_Momentum_Set_Stage_boundary(domain, function)
 
-#Bts = Transmissive_Momentum_Set_Stage_boundary(domain, waveform)
 Bw = Time_boundary(domain=domain,
                    f=lambda t: [(0.1*sin(t*2*pi)), 0.0, 0.0])
@@ -281 +284 @@
     angle = ((angle-180)/180.0)*pi
 #    angle = angle1*3.14157
-    p1 = get_xy(x,y,r1,angle-0.01745*d)
-    p2 = get_xy(x,y,r2,angle-0.01745*d)
-    p3 = get_xy(x,y,r2,angle+0.01745*d)
-    p4 = get_xy(x,y,r1,angle+0.01745*d)
+    p1 = get_xy(center_x,center_y,r1,angle-0.01745*d)
+    p2 = get_xy(center_x,center_y,r2,angle-0.01745*d)
+    p3 = get_xy(center_x,center_y,r2,angle+0.01745*d)
+    p4 = get_xy(center_x,center_y,r1,angle+0.01745*d)
 #    print p1,p2,p3,p4,angle
     poly_segment =[[p1[0],p1[1]],[p2[0],p2[1]],[p3[0],p3[1]],[p4[0],p4[1]]]
@@ -290 +293 @@
     #print i,poly[i]
     
-    run_up, x_y = get_maximum_inundation_data(filename='test.sww',polygon=poly_segment, verbose=False)
-    
-    print 'maximum_inundation_data',((angle/pi)*180)+180, run_up, x_y 
-
-
-
-
-
-
-
+    run_up, x_y = get_maximum_inundation_data(filename=model_name+'.sww',polygon=poly_segment, verbose=False)
+    
+    print 'maxd %.3f, %.6f, %.6f, %.5s, %.5s '%(((angle/pi)*180)+180, run_up*100, run_up-water_height*100, x_y[0],x_y[1]) 
+
+
+
+
+
+
+