Changeset 5827 for anuga_work/development/anuga_1d/dry_dam_sudi.py
- Timestamp:
- Oct 9, 2008, 12:54:08 PM (15 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
anuga_work/development/anuga_1d/dry_dam_sudi.py
r5587 r5827 1 1 import os 2 2 from math import sqrt, pi 3 from shallow_water_domain import *3 from shallow_water_domain_suggestion1 import * 4 4 from Numeric import allclose, array, zeros, ones, Float, take, sqrt 5 5 from config import g, epsilon … … 8 8 9 9 #t = 0.0 # time (s) 10 g = 9.81# gravity (m/s^2)10 # gravity (m/s^2) 11 11 h1 = 10.0 # depth upstream (m) 12 12 h0 = 0.0 # depth downstream (m) … … 43 43 h[i] = h0 44 44 45 return h , u*h 45 return h , u*h, u 46 46 47 47 #def newLinePlot(title='Simple Plot'): … … 65 65 print "TEST 1D-SOLUTION III -- DRY BED" 66 66 67 L = 2000.0 # Length of channel (m)68 N = 800 # Number of computational cells69 cell_len = L/N # Origin = 0.070 71 points = zeros(N+1,Float)72 for i in range(N+1):73 points[i] = i*cell_len74 75 domain = Domain(points)76 77 67 def stage(x): 78 68 h1 = 10.0 … … 92 82 93 83 finaltime = 20.0 94 yieldstep = 1.084 yieldstep = 20.0 95 85 L = 2000.0 # Length of channel (m) 96 number_of_cells = [ 200]#,200,500,1000,2000,5000,10000,20000]86 number_of_cells = [810]#,200,500,1000,2000,5000,10000,20000] 97 87 h_error = zeros(len(number_of_cells),Float) 98 88 uh_error = zeros(len(number_of_cells),Float) … … 110 100 domain.set_quantity('stage', stage) 111 101 domain.set_boundary({'exterior': Reflective_boundary(domain)}) 112 domain. default_order = 2113 domain. default_time_order = 2102 domain.order = 2 103 domain.set_timestepping_method('rk2') 114 104 domain.cfl = 1.0 115 domain.limiter = " minmod"105 domain.limiter = "vanleer" 116 106 117 107 t0 = time.time() … … 124 114 XmomC = domain.quantities['xmomentum'].centroid_values 125 115 C = domain.centroids 126 h, uh = analytical_sol(C,domain.time)116 h, uh, u = analytical_sol(C,domain.time) 127 117 h_error[k] = 1.0/(N)*sum(abs(h-StageC)) 128 118 uh_error[k] = 1.0/(N)*sum(abs(uh-XmomC)) … … 134 124 StageQ = domain.quantities['stage'].vertex_values 135 125 XmomQ = domain.quantities['xmomentum'].vertex_values 136 h, uh = analytical_sol(X.flat,domain.time) 126 velQ = domain.quantities['velocity'].vertex_values 127 128 h, uh, u = analytical_sol(X.flat,domain.time) 137 129 x = X.flat 138 130 … … 151 143 'upper right', shadow=True) 152 144 plot2 = subplot(212) 153 plot(x,u h,x,XmomQ.flat)145 plot(x,u,x,velQ.flat) 154 146 plot2.set_ylim([-35,35]) 155 147 156 148 xlabel('Position') 157 ylabel(' Xmomentum')149 ylabel('Velocity') 158 150 159 151 file = "dry_bed_"
Note: See TracChangeset
for help on using the changeset viewer.