Changeset 5289

Timestamp:

May 7, 2008, 5:02:52 PM (16 years ago)

Author:

kristy

Message:

cone is created, working on boundary and wave dimension

File:

• anuga_validation/circular_island/simple_circular/circular.py (modified) (7 diffs)

anuga_validation/circular_island/simple_circular/circular.py

@@ -13 +13 @@
 from anuga.shallow_water import Time_boundary
 from anuga.pmesh.mesh_interface import create_mesh_from_regions
+from anuga.shallow_water import Transmissive_Momentum_Set_Stage_boundary
+from math import tan, sqrt, sin, pi
 
 #------------------------------------------------------------------------------
@@ -19 +21 @@
 length = 30.
 width = 25.
+Cx = 12.96                  # centre of island on the x axis
+Cy = 13.8                   # centre of island on the y axis
 dx = dy = 1     # Resolution: Length of subdivisions on both axes
+water_depth = 0.32
 
-#Define polygon so that 0,0 is the centre - unstructured mesh
-##poly_domain = [[-width/2,length/2],[width/2,length/2],
-##               [width/2,-length/2],[-width/2,-length/2]]
+#boundary
 poly_domain = [[0,0],[length,0],[length,width],[0,width]]
 meshname = 'test.msh'
@@ -29 +32 @@
 
 #Create interior region
-#Dome = [[-4,4],[4,4],[4,-4],[-4,-4]]
-Dome = [[(length/2)-4,(width/2)-4],[(length/2)+4,(width/2)-4,],
-        [(length/2)+4,(width/2)+4],[(length/2)-4,(width/2)+4]]
-remainder_res=5
+Dome = [[(Cx)-4,(Cy)-4],[(Cx)+4,(Cy)-4,],
+        [(Cx)+4,(Cy)+4],[(Cx)-4,(Cy)+4]]
+
+remainder_res=1
 Dome_res = .01
+
 interior_dome = [[Dome, Dome_res]]
 
@@ -43 +47 @@
                          filename=meshname,
                          interior_regions = interior_dome,
-                         use_cache=False,
+                         use_cache=True,
                          verbose=True)
 
-domain = Domain(meshname, use_cache=False, verbose = True)
+domain = Domain(meshname, use_cache=True, verbose = True)
 domain.set_name('circular')                  # Output name
 print domain.statistics()
@@ -57 +61 @@
     """
     
-    water_depth = 0.32                                
-    z= 0*x
+                                    
+    z= 0*x                      # defining z for all values other than the if statements
+    r= 3.6                      # radius, provided in document
+    angle = 14                  # angle, provided in document
+    h= r*tan(angle/57.2957795)  # finding height of cone if not truncated
+    
 
     N = len(x)
-##    print 'Hello',N
-##    from Numeric import zeros
-##    z = zeros(N)
     for i in range(N):
-
-        """# Square
-        if (-2<x[i]-length/2<2) and (-2 <y[i]-width/2< 2):
-            z[i] += 2"""
-        #print 'Hello',i,x[i],y[i],z[i]
-        # dome
-        #if (x[i]-length/2)**2 + (y[i]-width/2)**2 <1.1**2:
-        if (x[i]-length/2)**2 + (y[i]-width/2)**2 <1.1**2:
+       
+        #truncated top
+        if (x[i]-Cx)**2 + (y[i]-Cy)**2 <1.1**2:
             z[i] += 0.625
 
-##        if (x[i]-length/2)**2 + (y[i]-width/2)**2 >3.6**2:
-##            z[i] += 0
-
-        if (x[i]-length/2)**2 + (y[i]-width/2)**2 <3.6**2 and (x[i]-length/2)**2 + (y[i]-width/2)**2 >1.1**2:
-           z[i] = (x[i])/4 + 15
-
+        # cone
+        if (x[i]-Cx)**2 + (y[i]-Cy)**2 <r**2 and (x[i]-Cx)**2 + (y[i]-Cy)**2 >1.1**2:
+           z[i] = -(sqrt(((x[i]-Cx)**2+(y[i]-Cy)**2)/((r/h)**2))-h)
     return z   
-
-#The cone has a 1/4 slope and a 7.2 diameter this make the height 0.9m
-#    to get the stage at 0 elevation + cone height and - water depth
-
-            
-            
-
-##        if z <0-water_depth:
-##            z[i] =0-water_depth
-##        if z > .625-water_depth:
-##            z[i] =0.625-water_depth
-##        z.append(z)
-    
-       
-
 
 domain.set_quantity('elevation', topography, verbose=True)  # Use function for elevation
 domain.set_quantity('friction', 0.01)         # Constant friction 
-domain.set_quantity('stage',
-                    expression='elevation')   # Dry initial condition
+domain.set_quantity('stage',water_depth)   # Dry initial condition
 
 
@@ -107 +88 @@
 # Setup boundary conditions
 #------------------------------------------------------------------------------
-Bi = Dirichlet_boundary([0.4, 0, 0])          # Inflow
-Br = Reflective_boundary(domain)              # Solid reflective wall
-Bo = Dirichlet_boundary([-5, 0, 0])           # Outflow
+# Create boundary function from timeseries provided in file
 
-domain.set_boundary({'left': Bi, 'right': Bo, 'top': Br, 'wavemaker': Br})
+##boundary_filename="ts2cnew1_input_20_80sec_new"
+##prepare_timeboundary(boundary_filename+'.txt')
+##
+##function = file_function(boundary_filename+'.tms',
+##                         domain, verbose=True)
+def wave_form(t):
+    return 0.1*sin(2*pi*t/50.)
+
+# Create and assign boundary objects
+Bw = Dirichlet_boundary([water_depth, 0, 0])          #wall
+Bt = Transmissive_Momentum_Set_Stage_boundary(domain, wave_form) #wavemaker
+domain.set_boundary({'left': Bw, 'right': Bw, 'top': Bw, 'wavemaker': Bt})
+
+
 
 
@@ -117 +109 @@
 # Evolve system through time
 #------------------------------------------------------------------------------
-for t in domain.evolve(yieldstep = 0.2, finaltime = 16.0):
+for t in domain.evolve(yieldstep = 0.2, finaltime = 100.0):
     print domain.timestepping_statistics()
 
 
-    if domain.get_quantity('stage').\
-           get_values(interpolation_points=[[10, 2.5]]) > 0:        
-        print 'Stage > 0: Changing to outflow boundary'
-        domain.set_boundary({'right': Bo})
+