Changeset 3913
- Timestamp:
- Nov 3, 2006, 11:35:32 AM (18 years ago)
- Location:
- anuga_validation/okushiri_2005
- Files:
-
- 4 edited
Legend:
- Unmodified
- Added
- Removed
-
anuga_validation/okushiri_2005/compare_timeseries.py
r3903 r3913 5 5 """ 6 6 7 from Numeric import allclose 7 from Numeric import allclose, argmin, argmax 8 8 from Scientific.IO.NetCDF import NetCDFFile 9 9 … … 67 67 68 68 69 expected_covariances = {'Boundary': 5.299487474489660856e-05, 70 'ch5': 1.169650160375980370e-04, 71 'ch7': 1.123947836450141360e-04, 72 'ch9': 1.248329330436513066e-04} 73 74 expected_differences = {'Boundary': 8.269932106586290379e-04, 75 'ch5': 3.411769502897898498e-03, 76 'ch7': 2.777024544282339583e-03, 77 'ch9': 3.183530359567784788e-03} 78 79 69 #expected_covariances = {'Boundary': 5.299487474489660856e-05, 70 # 'ch5': 1.169650160375980370e-04, 71 # 'ch7': 1.123947836450141360e-04, 72 # 'ch9': 1.248329330436513066e-04} 73 # 74 #expected_differences = {'Boundary': 8.269932106586290379e-04, 75 # 'ch5': 3.411769502897898498e-03, 76 # 'ch7': 2.777024544282339583e-03, 77 # 'ch9': 3.183530359567784788e-03} 78 # 79 80 81 expected_covariance = {'Boundary': 5.265850426352561254e-05, 82 'ch5': 1.171913883449135149e-04, 83 'ch7': 1.130834789656423632e-04, 84 'ch9': 1.257648236914514148e-04} 85 86 expected_difference = {'Boundary': 9.188561500629542633e-04, 87 'ch5': 3.439740513194261117e-03, 88 'ch7': 2.868914369295409543e-03, 89 'ch9': 3.287722754612562841e-03} 90 91 expected_maximum_diff = {'Boundary': 4.865498851439054029e-05, 92 'ch5': 1.808007760316858448e-03, 93 'ch7': 5.464055727564115506e-04, 94 'ch9': 1.856547605566957748e-03} 95 96 expected_minimum_diff = {'Boundary': 2.293007809784416290e-04, 97 'ch5': 1.494336068527772621e-04, 98 'ch7': 4.515664278503157131e-03, 99 'ch9': 3.406785970601119290e-03} 100 101 expected_argmax_timelag = {'Boundary': 3.000000000000007105e-01, 102 'ch5': 9.999999999999786837e-02, 103 'ch7': 0.000000000000000000e+00, 104 'ch9': 0.000000000000000000e+00} 105 106 expected_argmin_timelag = {'Boundary': 3.000000000000007105e-01, 107 'ch5': 8.499999999999996447e-01, 108 'ch7': 4.500000000000010658e-01, 109 'ch9': 9.499999999999992895e-01} 110 111 112 # Results from lwru2_variable_mesh.sww 113 # 114 #Validating Boundary 115 #Covariance = 5.265850426352561254e-05 116 #Accumulated difference = 9.188561500629542633e-04 117 #Difference in maxima = 4.865498851439054029e-05 118 #Difference in minima = 2.293007809784416290e-04 119 #Timelag between maxima = 3.000000000000007105e-01 120 #Timelag between minima = 3.000000000000007105e-01 121 #Next 122 # 123 #Validating ch5 124 #Covariance = 1.171913883449135149e-04 125 #Accumulated difference = 3.439740513194261117e-03 126 #Difference in maxima = 1.808007760316858448e-03 127 #Difference in minima = 1.494336068527772621e-04 128 #Timelag between maxima = 9.999999999999786837e-02 129 #Timelag between minima = 8.499999999999996447e-01 130 #Next# 131 # 132 #Validating ch7 133 #Covariance = 1.130834789656423632e-04 134 #Accumulated difference = 2.868914369295409543e-03 135 #Difference in maxima = 5.464055727564115506e-04 136 #Difference in minima = 4.515664278503157131e-03 137 #Timelag between maxima = 0.000000000000000000e+00 138 #Timelag between minima = 4.500000000000010658e-01 139 # 140 #Next 141 #Validating ch9 142 #Covariance = 1.257648236914514148e-04 143 #Accumulated difference = 3.287722754612562841e-03 144 #Difference in maxima = 1.856547605566957748e-03 145 #Difference in minima = 3.406785970601119290e-03 146 #Timelag between maxima = 0.000000000000000000e+00 147 #Timelag between minima = 9.499999999999992895e-01 80 148 81 149 … … 149 217 model.append(f(t, point_id=k)[0]) 150 218 219 151 220 # Covariance measures 152 221 res = cov(observed_timeseries, model) 153 222 print 'Covariance = %.18e' %res 154 223 155 if res < expected_covariances[name]-eps: 156 print 'Result is better than expected by: %.18e'\ 157 %(res-expected_covariances[name]) 158 print 'Expect = %.18e' %expected_covariances[name] 159 elif res > expected_covariances[name]+eps: 160 print 'FAIL: Result is worse than expected by: %.18e'\ 161 %(res-expected_covariances[name]) 162 print 'Expect = %.18e' %expected_covariances[name] 163 else: 164 pass 224 if res < expected_covariance[name]-eps: 225 print ' Result is better than expected by: %.18e'\ 226 %(res-expected_covariance[name]) 227 print ' Expect = %.18e' %expected_covariance[name] 228 elif res > expected_covariance[name]+eps: 229 print ' FAIL: Result is worse than expected by: %.18e'\ 230 %(res-expected_covariance[name]) 231 print ' Expect = %.18e' %expected_covariance[name] 232 else: 233 pass 234 165 235 166 236 # Difference measures … … 168 238 print 'Accumulated difference = %.18e' %res 169 239 170 if res < expected_differences[name]-eps: 171 print 'Result is better than expected by: %.18e'\ 172 %(res-expected_differences[name]) 173 print 'Expect = %.18e' %expected_differences[name] 174 elif res > expected_differences[name]+eps: 175 print 'FAIL: Result is worse than expected by: %.18e'\ 176 %(res-expected_differences[name]) 177 print 'Expect = %.18e' %expected_differences[name] 240 if res < expected_difference[name]-eps: 241 print ' Result is better than expected by: %.18e'\ 242 %(res-expected_difference[name]) 243 print ' Expect = %.18e' %expected_difference[name] 244 elif res > expected_difference[name]+eps: 245 print ' FAIL: Result is worse than expected by: %.18e'\ 246 %(res-expected_difference[name]) 247 print ' Expect = %.18e' %expected_difference[name] 248 else: 249 pass 250 251 252 # Extrema 253 res = abs(max(observed_timeseries)-max(model)) 254 print 'Difference in maxima = %.18e' %res 255 256 if res < expected_maximum_diff[name]-eps: 257 print ' Result is better than expected by: %.18e'\ 258 %(res-expected_maximum_diff[name]) 259 print ' Expect = %.18e' %expected_maximum_diff[name] 260 elif res > expected_maximum_diff[name]+eps: 261 print ' FAIL: Result is worse than expected by: %.18e'\ 262 %(res-expected_maximum_diff[name]) 263 print ' Expect = %.18e' %expected_maximum_diff[name] 264 else: 265 pass 266 267 268 269 res = abs(min(observed_timeseries)-min(model)) 270 print 'Difference in minima = %.18e' %res 271 272 if res < expected_minimum_diff[name]-eps: 273 print ' Result is better than expected by: %.18e'\ 274 %(res-expected_minimum_diff[name]) 275 print ' Expect = %.18e' %expected_minimum_diff[name] 276 elif res > expected_minimum_diff[name]+eps: 277 print ' FAIL: Result is worse than expected by: %.18e'\ 278 %(res-expected_minimum_diff[name]) 279 print ' Expect = %.18e' %expected_minimum_diff[name] 280 else: 281 pass 282 283 284 # Locations of extrema 285 i0 = argmax(observed_timeseries) 286 i1 = argmax(model) 287 res = abs(reference_time[i1] - reference_time[i0]) 288 print 'Timelag between maxima = %.18e' %res 289 290 if res < expected_argmax_timelag[name]-eps: 291 print ' Result is better than expected by: %.18e'\ 292 %(res-expected_argmax_timelag[name]) 293 print ' Expect = %.18e' %expected_argmax_timelag[name] 294 elif res > expected_argmax_timelag[name]+eps: 295 print ' FAIL: Result is worse than expected by: %.18e'\ 296 %(res-expected_argmax_timelag[name]) 297 print ' Expect = %.18e' %expected_argmax_timelag[name] 298 else: 299 pass 300 301 302 i0 = argmin(observed_timeseries) 303 i1 = argmin(model) 304 res = abs(reference_time[i1] - reference_time[i0]) 305 print 'Timelag between minima = %.18e' %res 306 if res < expected_argmin_timelag[name]-eps: 307 print ' Result is better than expected by: %.18e'\ 308 %(res-expected_argmin_timelag[name]) 309 print ' Expect = %.18e' %expected_argmin_timelag[name] 310 elif res > expected_argmin_timelag[name]+eps: 311 print ' FAIL: Result is worse than expected by: %.18e'\ 312 %(res-expected_argmin_timelag[name]) 313 print ' Expect = %.18e' %expected_argmin_timelag[name] 178 314 else: 179 315 pass … … 186 322 hold(False) 187 323 188 plot(reference_time, validation_data[name], 'r .-',324 plot(reference_time, validation_data[name], 'r-', 189 325 reference_time, model, 'k-') 190 326 title('Gauge %s' %name) -
anuga_validation/okushiri_2005/create_okushiri.py
r3901 r3913 152 152 153 153 154 # Island area and drawdown region 154 # Island area and drawdown region (original) 155 155 island_0 = [xleft + 2*(xright-xleft)/3+1.2, ytop-0.5] 156 156 island_1 = [xleft + 2*(xright-xleft)/3+0.5, ybottom + 2*(ytop-ybottom)/3] 157 island_2 = [xleft + (xright-xleft)/2+0.4, ybottom + 2*(ytop-ybottom)/3-0.3] 158 island_3 = [xleft + (xright-xleft)/2+0.4, ybottom + (ytop-ybottom)/3+0.3] 159 island_4 = [xleft + 2*(xright-xleft)/3+0.4, ybottom + (ytop-ybottom)/3-0.3] 160 island_5 = [xleft + 2*(xright-xleft)/3+1.2, ybottom+0.8] 161 island_6 = [xl-.01, yb] # Keep right edge just off the gulleys 162 island_7 = [xl-.01, yt] 157 island_2 = [xleft + (xright-xleft)/2+0.3, ybottom + 2*(ytop-ybottom)/3-0.3] 158 island_3 = [xleft + (xright-xleft)/2+0.3, ybottom + (ytop-ybottom)/3+0.3] 159 island_4 = [xleft + 2*(xright-xleft)/3+0.4, ybottom + (ytop-ybottom)/3-0.3] 160 island_5 = [xleft + 2*(xright-xleft)/3+1.2, ybottom+0.2] 161 island_6 = [xl-.01, yb] #OK 162 island_7 = [xl-.01, yt] #OK 163 164 165 # Island area and drawdown region 166 #island_0 = [xleft + 2*(xright-xleft)/3+1.2, ytop-0.5] 167 #island_1 = [xleft + 2*(xright-xleft)/3+0.5, ybottom + 2*(ytop-ybottom)/3] 168 #island_2 = [xleft + (xright-xleft)/2+0.4, ybottom + 2*(ytop-ybottom)/3-0.3] 169 #island_3 = [xleft + (xright-xleft)/2+0.4, ybottom + (ytop-ybottom)/3+0.3] 170 #island_4 = [xleft + 2*(xright-xleft)/3+0.4, ybottom + (ytop-ybottom)/3-0.3] 171 #island_5 = [xleft + 2*(xright-xleft)/3+1.2, ybottom+0.8] 172 #island_6 = [xl-.01, yb] # Keep right edge just off the gulleys 173 #island_7 = [xl-.01, yt] 163 174 164 175 island = [island_0, island_1, island_2, … … 167 178 168 179 180 # Region spanning half right hand side of domain just inside boundary (org) 181 rhs_nw = [xleft + (xright-xleft)/3+1, ytop-0.02] 182 rhs_sw = [xleft + (xright-xleft)/3+1, ybottom+0.02] 183 rhs_se = [xright-0.02, ybottom+0.02] 184 rhs_ne = [xright-0.02, ytop-0.02] 185 169 186 # Region spanning half right hand side of domain just inside boundary 170 rhs_nw = [xleft + (xright-xleft)/3+1, ytop-1.4]171 rhs_sw = [xleft + (xright-xleft)/3+1, ybottom+0.5]172 rhs_se = [xright-0.1, ybottom+0.2]173 rhs_ne = [xright-0.1, ytop-0.2]187 #rhs_nw = [xleft + (xright-xleft)/3+1, ytop-1.4] 188 #rhs_sw = [xleft + (xright-xleft)/3+1, ybottom+0.5] 189 #rhs_se = [xright-0.1, ybottom+0.2] 190 #rhs_ne = [xright-0.1, ytop-0.2] 174 191 175 192 rhs_region = [rhs_nw, rhs_ne, rhs_se, rhs_sw] … … 179 196 # Interior regions and creation of mesh 180 197 interior_regions = [[rhs_region, 0.0005], 181 [ gulleys, 0.00002*base_resolution],182 [ island, 0.00007*base_resolution]]198 [island, 0.0003*base_resolution], 199 [gulleys, 0.00003*base_resolution]] 183 200 184 201 meshname = project.mesh_filename + '.msh' -
anuga_validation/okushiri_2005/create_okushiri_original.py
r3854 r3913 196 196 197 197 198 #inner = m.addRegionEN(island_3[0]+.1, island_3[1]+.1) 199 #inner.setMaxArea(0.00007*base_resolution) 200 198 201 inner = m.addRegionEN(island_3[0]+.1, island_3[1]+.1) 199 inner.setMaxArea(0.000 07*base_resolution)202 inner.setMaxArea(0.0003*base_resolution) 200 203 201 204 202 205 gulleys = m.addRegionEN(p0[0]+0.1, p0[1]+0.1) 203 gulleys.setMaxArea(0.0000 2*base_resolution)206 gulleys.setMaxArea(0.00003*base_resolution) 204 207 205 208 -
anuga_validation/okushiri_2005/project.py
r3883 r3913 17 17 18 18 # Model output 19 output_filename = 'okushiri .sww'19 output_filename = 'okushiri_new.sww' 20 20 21 21
Note: See TracChangeset
for help on using the changeset viewer.