Changeset 7348

Timestamp:

Aug 10, 2009, 5:40:57 AM (16 years ago)

Author:

ole

Message:

Elaborated on variable bedslope example for use with AOGS2009 talk

File:

• anuga_core/documentation/user_manual/demos/channel_variable.py (modified) (7 diffs)

anuga_core/documentation/user_manual/demos/channel_variable.py

@@ -16 +16 @@
 # Setup computational domain
 #------------------------------------------------------------------------------
-length = 40.
+length = 24.
 width = 5.
-dx = dy = 1 #.1           # Resolution: Length of subdivisions on both axes
+dx = dy = 0.05 #.1           # Resolution: Length of subdivisions on both axes
 
 points, vertices, boundary = rectangular_cross(int(length/dx), int(width/dy),
                                                len1=length, len2=width)
 domain = Domain(points, vertices, boundary)
-domain.set_name('channel_variable')                  # Output name
+domain.set_name('channel_variable_2')                  # Output name
 print domain.statistics()
 domain.set_quantities_to_be_stored({'elevation': 1, #2,
@@ -35 +35 @@
 
     z = -x/100
-
+    
     N = len(x)
     for i in range(N):
         # Step
-        if 10 < x[i] < 12:
+        if 2 < x[i] < 4:
             z[i] += 0.4 - 0.05*y[i]
-
-        # Constriction
-        #if 27 < x[i] < 29 and y[i] > 3:
-        #    z[i] += 2
-
-        # Pole
-        #if (x[i] - 34)**2 + (y[i] - 2)**2 < 0.4**2:
-        #    z[i] += 2
-
+    
+        # Permanent pole
+        if (x[i] - 8)**2 + (y[i] - 2)**2 < 0.4**2:
+            z[i] += 1
     return z
+    
 
 def pole_increment(x,y):
@@ -61 +57 @@
     N = len(x)
     for i in range(N):
-        # Pole
-        if (x[i] - 20)**2 + (y[i] - 2)**2 < 0.4**2:
-            z[i] += 0.05
+        # Pole 1
+        if (x[i] - 12)**2 + (y[i] - 3)**2 < 0.4**2:
+            z[i] += 0.01
+            
+    for i in range(N):
+        # Pole 2
+        if (x[i] - 14)**2 + (y[i] - 2)**2 < 0.4**2:
+            z[i] += 0.005
+                        
     return z
     
@@ -87 +89 @@
 shrinking = False
 done = False
-for t in domain.evolve(yieldstep=0.1, finaltime=30.0):
+for t in domain.evolve(yieldstep=0.1, finaltime=40.0):
     print domain.timestepping_statistics()
 
@@ -94 +96 @@
     #print 'Level at gauge point = %.2fm' % w
            
-    z = domain.get_quantity('elevation').\
-        get_values(interpolation_points=[[20, 2]])           
-    print 'Elevation at pole location = %.2fm' % z           
+    #z = domain.get_quantity('elevation').\
+    #    get_values(interpolation_points=[[12, 3]])           
+    #print 'Elevation at pole location = %.2fm' % z           
 
     # Start variable elevation after 10 seconds    
@@ -102 +104 @@
         growing = True
            
-    # Turn around when pole has reached a height of 2 m
-    if z >= 2 and growing:
+    # Stop growing when pole has reached a certain height
+    if t > 16 and growing:
         growing = False
+        shrinking = False
+        
+    # Start shrinking
+    if t > 20:
         shrinking = True
+        growing = False
         
-    # Stop changing when pole has shrunk to 1 m
-    if z <= 1 and shrinking:
+    # Stop changing when pole has shrunk to original level
+    if t > 25 and shrinking:
         done = True
         shrinking = growing = False
+        domain.set_quantity('elevation', topography)
 
     # Grow or shrink               
@@ -117 +125 @@
         
     if shrinking:    
-        domain.add_quantity('elevation', lambda x,y: -pole_increment(x,y))    
+        domain.add_quantity('elevation', lambda x,y: -2*pole_increment(x,y))