Changeset 5162 for anuga_validation/convergence_study/convergence_structured.py

Timestamp:

Mar 11, 2008, 8:43:22 PM (17 years ago)

Author:

steve

Message:

Updated some methods for quantity. Looks like we can use old
limiting system with larger values of beta.

File:

• anuga_validation/convergence_study/convergence_structured.py (modified) (5 diffs)

anuga_validation/convergence_study/convergence_structured.py

@@ -17 +17 @@
 from anuga.shallow_water import Transmissive_boundary
 from anuga.shallow_water import Transmissive_Momentum_Set_Stage_boundary
-from anuga.geospatial_data.geospatial_data import *
+#from anuga.geospatial_data.geospatial_data import *
 from math import cos
+from Numeric import zeros, Float
 
 #------------------------------------------------------------------------------
 # Setup computational domain
 #------------------------------------------------------------------------------
-dx = 1000.
+dx = 200.
 dy = dx
 L = 100000.
@@ -33 +34 @@
 domain = Domain(points, vertices, boundary) 
 
-domain.set_timestepping_method('euler')
+domain.set_timestepping_method('rk2')
 domain.set_default_order(2)
 domain.set_name('myexample9')                
 domain.set_datadir('.')                     # Use current directory for output
 
-domain.beta_w      = 1.0
-domain.beta_w_dry  = 0.2
-domain.beta_uh     = 1.0
-domain.beta_uh_dry = 0.2
-domain.beta_vh     = 1.0
-domain.beta_vh_dry = 0.2
-domain.beta_h      = 1.0
+domain.set_all_limiters(0.9)
+
+print domain.beta_w
+domain.use_old_limiter = False
+domain.CFL = 1.0
 
 #------------------------------------------------------------------------------
@@ -50 +49 @@
 #------------------------------------------------------------------------------
 #domain.set_quantity('elevation', topography) # Use function for elevation
-domain.set_quantity('elevation',-100)
+domain.set_quantity('elevation',0.0)
 domain.set_quantity('friction', 0.00)
-domain.set_quantity('stage', 0.0)            
+
+h0 = 10.0
+h1 = 1.0
+
+def height(x,y):
+    z = zeros(len(x),Float)
+    for i in range(len(x)):
+        if x[i]<=50000.0:
+            z[i] = h0
+        else:
+            z[i] = h1
+    return z
+
+
+domain.set_quantity('stage', height)
+#domain.set_quantity('stage', 0.0)            
 
 #-----------------------------------------------------------------------------
@@ -65 +79 @@
 Bw = Time_boundary(domain=domain,     # Time dependent boundary  
 ## Sine wave
-                   f=lambda t: [(-amplitude*sin((1./300.)*t*2*pi)), 0.0, 0.0])
+#                   f=lambda t: [(-amplitude*sin((1./300.)*t*2*pi)), 0.0, 0.0])
+## Single wave
+                   f=lambda t: [h0, 0.0, 0.0])
 ## Sawtooth?
 #                   f=lambda t: [(-8.0*(sin((1./180.)*t*2*pi))+(1./2.)*sin((2./180.)*t*2*pi)+(1./3.)*sin((3./180.)*t*2*pi)), 0.0, 0.0])
@@ -92 +108 @@
 #------------------------------------------------------------------------------
 
-for t in domain.evolve(yieldstep = 20.0, finaltime = 10*40*60.):
+for t in domain.evolve(yieldstep = 50.0, finaltime = 10*40*60.):
     domain.write_time()
     vis.update()