Changeset 5827 for anuga_work/development/anuga_1d/dry_dam_sudi.py

Timestamp:

Oct 9, 2008, 12:54:08 PM (16 years ago)

Author:

sudi

Message:

Adding new versions of shallow_water_domain

File:

• anuga_work/development/anuga_1d/dry_dam_sudi.py (modified) (9 diffs)

anuga_work/development/anuga_1d/dry_dam_sudi.py

@@ -1 +1 @@
 import os
 from math import sqrt, pi
-from shallow_water_domain import *
+from shallow_water_domain_suggestion1 import *
 from Numeric import allclose, array, zeros, ones, Float, take, sqrt
 from config import g, epsilon
@@ -8 +8 @@
     
     #t  = 0.0     # time (s)
-    g  = 9.81    # gravity (m/s^2)
+    # gravity (m/s^2)
     h1 = 10.0    # depth upstream (m)
     h0 = 0.0     # depth downstream (m)
@@ -43 +43 @@
             h[i] = h0 
 
-    return h , u*h
+    return h , u*h, u
 
 #def newLinePlot(title='Simple Plot'):
@@ -65 +65 @@
 print "TEST 1D-SOLUTION III -- DRY BED"
 
-L = 2000.0     # Length of channel (m)
-N = 800        # Number of computational cells
-cell_len = L/N # Origin = 0.0
-
-points = zeros(N+1,Float)
-for i in range(N+1):
-    points[i] = i*cell_len
-
-domain = Domain(points)
-
 def stage(x):
     h1 = 10.0
@@ -92 +82 @@
 
 finaltime = 20.0
-yieldstep = 1.0
+yieldstep = 20.0
 L = 2000.0     # Length of channel (m)
-number_of_cells = [200]#,200,500,1000,2000,5000,10000,20000]
+number_of_cells = [810]#,200,500,1000,2000,5000,10000,20000]
 h_error = zeros(len(number_of_cells),Float)
 uh_error = zeros(len(number_of_cells),Float)
@@ -110 +100 @@
     domain.set_quantity('stage', stage)
     domain.set_boundary({'exterior': Reflective_boundary(domain)})
-    domain.default_order = 2
-    domain.default_time_order = 2
+    domain.order = 2
+    domain.set_timestepping_method('rk2')
     domain.cfl = 1.0
-    domain.limiter = "minmod"
+    domain.limiter = "vanleer"
 
     t0 = time.time()
@@ -124 +114 @@
     XmomC = domain.quantities['xmomentum'].centroid_values
     C = domain.centroids
-    h, uh = analytical_sol(C,domain.time)
+    h, uh, u = analytical_sol(C,domain.time)
     h_error[k] = 1.0/(N)*sum(abs(h-StageC))
     uh_error[k] = 1.0/(N)*sum(abs(uh-XmomC))
@@ -134 +124 @@
     StageQ = domain.quantities['stage'].vertex_values
     XmomQ = domain.quantities['xmomentum'].vertex_values
-    h, uh = analytical_sol(X.flat,domain.time)
+    velQ = domain.quantities['velocity'].vertex_values
+    
+    h, uh, u = analytical_sol(X.flat,domain.time)
     x = X.flat
     
@@ -151 +143 @@
            'upper right', shadow=True)
     plot2 = subplot(212)
-    plot(x,uh,x,XmomQ.flat)
+    plot(x,u,x,velQ.flat)
     plot2.set_ylim([-35,35])
     
     xlabel('Position')
-    ylabel('Xmomentum')
+    ylabel('Velocity')
    
     file = "dry_bed_"