Changeset 8234
- Timestamp:
- Oct 24, 2011, 5:37:21 PM (12 years ago)
- Location:
- trunk
- Files:
-
- 5 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/anuga_core/source/anuga/abstract_2d_finite_volumes/generic_domain.py
r8207 r8234 1317 1317 self.update_boundary() 1318 1318 1319 1319 1320 # Update extrema if necessary (for reporting) 1320 1321 self.update_extrema() … … 1453 1454 self.update_timestep(yieldstep, finaltime) 1454 1455 1455 # Update c onserved quantities1456 # Update centroid values of conserved quantities 1456 1457 self.update_conserved_quantities() 1457 1458 -
trunk/anuga_work/development/2010-projects/anuga_1d/channel/channel_domain.py
r8209 r8234 142 142 #distribute_to_vertices_and_edges_limit_s_v_h(self) 143 143 distribute_to_vertices_and_edges_limit_s_v(self) 144 144 145 def update_derived_quantites(self): 146 pass 145 147 146 148 #=============== End of Channel Domain =============================== … … 200 202 201 203 w_C[:] = h_C + z_C 204 205 print w_C 202 206 203 207 # Extrapolate velocity and stage as well as channel geometry. … … 217 221 b_V = width.vertex_values 218 222 219 # These quantites need to update vertex_values223 # Need to update these vertex_values 220 224 a_V = area.vertex_values 221 225 h_V = height.vertex_values -
trunk/anuga_work/development/2010-projects/anuga_1d/channel/channel_flow_1_padarn.py
r7884 r8234 6 6 from anuga_1d.channel.channel_domain import * 7 7 from anuga_1d.config import g, epsilon 8 from anuga_1d.base.generic_mesh import interval_mesh8 from anuga_1d.base.generic_mesh import uniform_mesh 9 9 10 10 … … 42 42 43 43 # Create domain with centroid points as defined above 44 domain = Domain(* interval_mesh(N))44 domain = Domain(*uniform_mesh(N)) 45 45 46 46 … … 69 69 70 70 71 exit()71 #exit() 72 72 73 73 HeightC = domain.quantities['height'].centroid_values -
trunk/anuga_work/development/2010-projects/anuga_1d/channel/channel_flow_2_padarn.py
r8204 r8234 36 36 import time 37 37 38 # Set final time and yield time for simulation 39 finaltime = 10.0 40 yieldstep = finaltime 38 41 39 42 40 # Length of channel (m) … … 76 74 t0 = time.time() 77 75 78 finaltime= 0.0 76 # Set final time and yield time for simulation 77 finaltime = 10.0 78 yieldstep = finaltime 79 79 80 80 #=================================================================== -
trunk/anuga_work/development/2010-projects/anuga_1d/channel/var_depth_only.py
r7884 r8234 2 2 import random 3 3 from math import sqrt, pow, pi 4 from channel_domain _Abimport *5 from Numeric import allclose, array, zeros, ones, Float, take, sqrt6 from config import g, epsilon4 from channel_domain import * 5 from numpy import allclose, array, zeros, ones, take, sqrt 6 from anuga_1d.config import g, epsilon 7 7 8 8 … … 10 10 11 11 # Define functions for initial quantities 12 def initial_area(x): 13 y = zeros(len(x),Float) 14 for i in range(len(x)): 15 y[i]=(10-randomarray[i]) 16 17 18 return y 12 13 14 15 #def initial_area(x): 16 # y = zeros(len(x),'f') 17 # for i in range(len(x)): 18 # y[i]=(10-randomarray[i]) 19 # return y 20 21 19 22 20 23 def bed(x): 21 24 22 y = zeros(len(x), Float)25 y = zeros(len(x),'f') 23 26 for i in range(len(x)): 24 27 y[i]=randomarray[i] … … 26 29 return y 27 30 28 def stage(x): 29 30 y = zeros(len(x),Float) 31 for i in range(len(x)): 32 y[i]=3+x[i]/2000 31 32 33 34 35 def width(x): 36 y = zeros(len(x),'f') 37 return y+1 38 39 40 stage = 6.0 41 42 def initial_area(x): 43 44 a_bed = bed(x) 45 a_width = width(x) 46 47 a_height = 6.0 - a_bed 48 49 y = a_height*a_width 50 33 51 return y 34 35 36 37 def width(x): 38 return 1 39 52 53 40 54 import time 41 55 … … 53 67 print "Evaluating domain with %d cells" %N 54 68 cell_len = L/N # Origin = 0.0 55 points = zeros(N+1, Float)69 points = zeros(N+1,'f') 56 70 57 71 # Define the centroid points … … 64 78 # Define random array for width 65 79 66 randomarray=zeros(len(points), Float)80 randomarray=zeros(len(points),'f') 67 81 for j in range(N+1): 68 82 randomarray[j]=random.normalvariate(3,1) … … 70 84 71 85 # Set initial values of quantities - default to zero 72 domain.set_quantity('stage',stage) 73 #domain.set_quantity('area', initial_area) 74 domain.set_quantity('elevation',bed) 75 domain.set_quantity('width',width) 76 domain.setstageflag = True 86 #domain.set_quantity('stage',6.0) 87 domain.set_quantity('elevation',bed, location='centroids') 88 domain.set_quantity('width',width, location='centroids') 89 domain.set_quantity('area', initial_area, location='centroids') 90 91 92 #domain.setstageflag = True 77 93 # Set boundry type, order, timestepping method and limiter 78 94 domain.set_boundary({'exterior':Reflective_boundary(domain)}) … … 83 99 #domain.h0=0.0001 84 100 101 102 #domain.distribute_to_vertices_and_edges() 103 104 105 #AreaC = domain.quantities['area'].centroid_values 106 #BedC = domain.quantities['elevation'].centroid_values 107 #WidthC = domain.quantities['width'].centroid_values 108 # 109 #AreaC[:] = (6.0 - BedC)* WidthC 110 111 #domain.set_quantity('area', initial_area) 112 113 #domain.distribute_to_vertices_and_edges() 114 85 115 # Start timer 86 116 t0 = time.time() 87 i=1 88 117 i=0 118 119 print 'elevation vertex values' 89 120 print domain.quantities['elevation'].vertex_values 121 print 'stage vertex values' 122 print domain.quantities['stage'].vertex_values 123 print 'area vertex values' 124 print domain.quantities['area'].vertex_values 125 print 'width vertex values' 126 print domain.quantities['width'].vertex_values 127 128 HeightC = domain.quantities['height'].centroid_values 129 DischargeC = domain.quantities['discharge'].centroid_values 130 C = domain.centroids 131 print 'That took %.2f seconds' %(time.time()-t0) 132 X = (domain.vertices).flatten() 133 HeightQ = (domain.quantities['height'].vertex_values).flatten() 134 VelocityQ = (domain.quantities['velocity'].vertex_values).flatten() 135 Z = (domain.quantities['elevation'].vertex_values).flatten() 136 stage=HeightQ + Z 137 138 139 from pylab import plot,title,xlabel,ylabel,legend,savefig,show,hold,subplot,savefig 140 import pylab as p 141 #hold(False) 142 143 plot1=subplot(211) 144 145 plot(X,stage,X,Z) 146 147 plot1.set_ylim([-1,11]) 148 xlabel('Position') 149 ylabel('Stage') 150 legend(('Solution', 'Bed Height'), 151 'upper right', shadow=True) 152 153 plot2=subplot(212) 154 155 156 plot(X,VelocityQ) 157 plot2.set_ylim([-5,5]) 158 xlabel('Position') 159 ylabel('Velocity') 160 161 savefig(str(i)+'.png') 162 i+=1 163 plot1.clear() 164 plot2.clear() 165 90 166 91 167 for t in domain.evolve(yieldstep = yieldstep, finaltime = finaltime): 92 168 domain.write_time() 93 169 94 N = float(N)95 170 HeightC = domain.quantities['height'].centroid_values 96 171 DischargeC = domain.quantities['discharge'].centroid_values 97 172 C = domain.centroids 98 173 print 'That took %.2f seconds' %(time.time()-t0) 99 X = domain.vertices 100 HeightQ = domain.quantities['height'].vertex_values 101 VelocityQ = domain.quantities['velocity'].vertex_values 102 x = X.flat 103 z = domain.quantities['elevation'].vertex_values.flat 104 stage=HeightQ.flat+z 105 106 174 X = (domain.vertices).flatten() 175 HeightQ = (domain.quantities['height'].vertex_values).flatten() 176 VelocityQ = (domain.quantities['velocity'].vertex_values).flatten() 177 Z = (domain.quantities['elevation'].vertex_values).flatten() 178 stage=HeightQ + Z 107 179 108 180 … … 113 185 plot1=subplot(211) 114 186 115 plot( x,stage,x,z)187 plot(X,stage,X,Z) 116 188 117 189 plot1.set_ylim([-1,11]) … … 122 194 123 195 plot2=subplot(212) 124 plot(x,VelocityQ.flat) 196 197 198 plot(X,VelocityQ) 125 199 plot2.set_ylim([-5,5]) 126 200 xlabel('Position') … … 130 204 i+=1 131 205 plot1.clear() 206 plot2.clear() 132 207 133 208 print domain.quantities['elevation'].vertex_values … … 138 213 C = domain.centroids 139 214 print 'That took %.2f seconds' %(time.time()-t0) 140 X = domain.vertices 141 HeightQ = domain.quantities['height'].vertex_values 142 VelocityQ = domain.quantities['velocity'].vertex_values 143 x = X.flat 144 z = domain.quantities['elevation'].vertex_values.flat 145 stage=HeightQ.flat+z 146 147 215 X = (domain.vertices).flatten() 216 HeightQ = (domain.quantities['height'].vertex_values).flatten() 217 VelocityQ = (domain.quantities['velocity'].vertex_values).flatten() 218 Z = (domain.quantities['elevation'].vertex_values).flatten() 219 stage=HeightQ + Z 148 220 149 221 … … 154 226 plot1=subplot(211) 155 227 156 plot( x,stage,x,z)228 plot(X,stage,X,Z) 157 229 158 230 plot1.set_ylim([-1,11]) … … 163 235 164 236 plot2=subplot(212) 165 plot( x,VelocityQ.flat)237 plot(X,VelocityQ) 166 238 plot2.set_ylim([-5,5]) 167 239 xlabel('Position')
Note: See TracChangeset
for help on using the changeset viewer.