Changeset 7348
 Timestamp:
 Aug 10, 2009, 5:40:57 AM (14 years ago)
 File:

 1 edited
Legend:
 Unmodified
 Added
 Removed

anuga_core/documentation/user_manual/demos/channel_variable.py
r7346 r7348 16 16 # Setup computational domain 17 17 # 18 length = 40.18 length = 24. 19 19 width = 5. 20 dx = dy = 1#.1 # Resolution: Length of subdivisions on both axes20 dx = dy = 0.05 #.1 # Resolution: Length of subdivisions on both axes 21 21 22 22 points, vertices, boundary = rectangular_cross(int(length/dx), int(width/dy), 23 23 len1=length, len2=width) 24 24 domain = Domain(points, vertices, boundary) 25 domain.set_name('channel_variable ') # Output name25 domain.set_name('channel_variable_2') # Output name 26 26 print domain.statistics() 27 27 domain.set_quantities_to_be_stored({'elevation': 1, #2, … … 35 35 36 36 z = x/100 37 37 38 38 N = len(x) 39 39 for i in range(N): 40 40 # Step 41 if 10 < x[i] < 12:41 if 2 < x[i] < 4: 42 42 z[i] += 0.4  0.05*y[i] 43 44 # Constriction 45 #if 27 < x[i] < 29 and y[i] > 3: 46 # z[i] += 2 47 48 # Pole 49 #if (x[i]  34)**2 + (y[i]  2)**2 < 0.4**2: 50 # z[i] += 2 51 43 44 # Permanent pole 45 if (x[i]  8)**2 + (y[i]  2)**2 < 0.4**2: 46 z[i] += 1 52 47 return z 48 53 49 54 50 def pole_increment(x,y): … … 61 57 N = len(x) 62 58 for i in range(N): 63 # Pole 64 if (x[i]  20)**2 + (y[i]  2)**2 < 0.4**2: 65 z[i] += 0.05 59 # Pole 1 60 if (x[i]  12)**2 + (y[i]  3)**2 < 0.4**2: 61 z[i] += 0.01 62 63 for i in range(N): 64 # Pole 2 65 if (x[i]  14)**2 + (y[i]  2)**2 < 0.4**2: 66 z[i] += 0.005 67 66 68 return z 67 69 … … 87 89 shrinking = False 88 90 done = False 89 for t in domain.evolve(yieldstep=0.1, finaltime= 30.0):91 for t in domain.evolve(yieldstep=0.1, finaltime=40.0): 90 92 print domain.timestepping_statistics() 91 93 … … 94 96 #print 'Level at gauge point = %.2fm' % w 95 97 96 z = domain.get_quantity('elevation').\97 get_values(interpolation_points=[[20, 2]])98 print 'Elevation at pole location = %.2fm' % z98 #z = domain.get_quantity('elevation').\ 99 # get_values(interpolation_points=[[12, 3]]) 100 #print 'Elevation at pole location = %.2fm' % z 99 101 100 102 # Start variable elevation after 10 seconds … … 102 104 growing = True 103 105 104 # Turn around when pole has reached a height of 2 m105 if z >= 2and growing:106 # Stop growing when pole has reached a certain height 107 if t > 16 and growing: 106 108 growing = False 109 shrinking = False 110 111 # Start shrinking 112 if t > 20: 107 113 shrinking = True 114 growing = False 108 115 109 # Stop changing when pole has shrunk to 1 m110 if z <= 1and shrinking:116 # Stop changing when pole has shrunk to original level 117 if t > 25 and shrinking: 111 118 done = True 112 119 shrinking = growing = False 120 domain.set_quantity('elevation', topography) 113 121 114 122 # Grow or shrink … … 117 125 118 126 if shrinking: 119 domain.add_quantity('elevation', lambda x,y:  pole_increment(x,y))127 domain.add_quantity('elevation', lambda x,y: 2*pole_increment(x,y)) 120 128
Note: See TracChangeset
for help on using the changeset viewer.