Changeset 9056

Timestamp:

Mar 3, 2014, 9:58:02 AM (11 years ago)

Author:

davies

Message:

Minor updates to validation tests

Location:

trunk

Files:

• anuga_core/source/anuga/config.py (modified) (1 diff)
• anuga_core/source/anuga_validation_tests/analytical_exact/rundown_mild_slope_coarse/numerical_rundown_channel_coarse.py (modified) (2 diffs)
• anuga_core/source/anuga_validation_tests/analytical_exact/rundown_mild_slope_coarse/plot_results.py (modified) (3 diffs)
• anuga_core/source/anuga_validation_tests/case_studies/patong/project.py (modified) (2 diffs)
• anuga_core/source/anuga_validation_tests/case_studies/patong/run_model.py (modified) (1 diff)
• anuga_core/source/anuga_validation_tests/saved_parameters.tex (modified) (1 diff)
• anuga_work/development/gareth/experimental/bal_and/swb2_domain.py (modified) (1 diff)
• anuga_work/development/gareth/experimental/bal_and/swb2_domain_ext.c (modified) (7 diffs)
• anuga_work/development/gareth/tests/channel_floodplain/channel_floodplain1.py (modified) (1 diff)
• anuga_work/development/gareth/tests/channel_floodplain/plotme.py (modified) (3 diffs)
• anuga_work/development/gareth/tests/runup_sinusoid/runup_sinusoid.py (modified) (1 diff)
• anuga_work/development/gareth/tests/runup_sinusoid/runup_sinusoidplot.py (modified) (1 diff)
• anuga_work/development/gareth/tests/shallow_steep_slope/channel_SU_2plot.py (modified) (4 diffs)
• anuga_work/development/gareth/tests/shallow_steep_slope/channel_SU_sparse.py (modified) (1 diff)
• anuga_work/development/gareth/tests/wave_runup/plotme.py (modified) (2 diffs)
• anuga_work/development/gareth/tests/wave_runup/runup.py (modified) (2 diffs)

trunk/anuga_core/source/anuga/config.py

@@ -126 +126 @@
 #flow_algorithm = '2_0'    # 2nd order TVD scheme and more aggressive piecewise linear spatial reconstruction
 #flow_algorithm = '2.5'  # 3rd order TVD scheme and more aggressive piecewise linear spatial reconstruction
-
+#flow_algorithm = 'tsunami' # 2nd order space and time, well balanced inc at wet-dry fronts, porosity-type alg
+#flow_algorithm = 'DE1' # 2nd order space and time, discontinuous elevation, well balanced + better shallow flows than 'tsunami'
 
 

trunk/anuga_core/source/anuga_validation_tests/analytical_exact/rundown_mild_slope_coarse/numerical_rundown_channel_coarse.py

@@ -41 +41 @@
 
 def init_stage(x,y):
-    stg= -x/10.+0.01 # Constant depth: 1 cm.
+    stg= -x/10.+0.10 # Constant depth: 10 cm.
     return stg
 #line0=[ [0.,0.], [0., 100.] ]
@@ -79 +79 @@
 # Evolve system through time
 #------------------------------------------------------------------------------
-for t in domain.evolve(yieldstep=10.0, finaltime=800.0):
+for t in domain.evolve(yieldstep=10.0, finaltime=2000.0):
     print domain.timestepping_statistics()
     print (domain.areas*(domain.quantities['stage'].centroid_values - domain.quantities['elevation'].centroid_values)).sum()

trunk/anuga_core/source/anuga_validation_tests/analytical_exact/rundown_mild_slope_coarse/plot_results.py

@@ -31 +31 @@
 
 ll=p2.vel.shape[0]-1
+
 #-------------------------------
 # Define variables of case study
@@ -55 +56 @@
 line.set_label('numerical')
 pyplot.plot( (0,100),(dana,dana), 'r',label='analytical' )
+pyplot.yscale('log')
 pyplot.legend()
 #pyplot.plot(x[v],elev[v],'green')
@@ -84 +86 @@
 pyplot.title('Final Xvelocity around the line x=50.')
 pyplot.legend(loc='best')
+pyplot.ylim((0,1.5*uana))
 pyplot.savefig('x_velocity.png')
 

trunk/anuga_core/source/anuga_validation_tests/case_studies/patong/project.py

@@ -35 +35 @@
 finaltime = 15000       # final time for simulation
 
+nameFlag='xxNameFlagxx'
+
 use_buildings = False   # FIXME: 6/6/2013 -- Code needs to be edited to support 'True'
 setup = 'basic'         # This can be one of three values
@@ -53 +55 @@
 #------------------------------------------------------------------------------
 
-output_comment = [setup, tide, event_number, 'final']
+output_comment = [setup, tide, event_number, nameFlag]
 
 #------------------------------------------------------------------------------

trunk/anuga_core/source/anuga_validation_tests/case_studies/patong/run_model.py

@@ -118 +118 @@
     domain.set_datadir(project.output_run) 
     domain.set_flow_algorithm('tsunami')
-    domain.minimum_allowed_height=0.001
+    #domain.set_flow_algorithm('DE1')
 
     #-------------------------------------------------------------------------------

trunk/anuga_core/source/anuga_validation_tests/saved_parameters.tex

@@ -3 +3 @@
 \newcommand{\majorR}{\UScore{1.3.0-beta}}
 \newcommand{\minorR}{\UScore{9003}}
-\newcommand{\timeR}{{Thu Dec  5 15:39:13 2013}}
+\newcommand{\timeR}{{Sun Feb 16 20:16:29 2014}}

trunk/anuga_work/development/gareth/experimental/bal_and/swb2_domain.py

@@ -73 +73 @@
 
         self.beta_w=1.0
-        self.beta_w_dry=0.4
+        self.beta_w_dry=0.0
         self.beta_uh=1.0
         self.beta_uh_dry=0.0

trunk/anuga_work/development/gareth/experimental/bal_and/swb2_domain_ext.c

@@ -242 +242 @@
       if(i==0) edgeflux[i] += (s_max*s_min)*(max(q_right_rotated[i],ze) - max(q_left_rotated[i],ze));
       //if(i==0) edgeflux[i] += (s_max*s_min)*(hre - hle);
-      if(i==1) edgeflux[i] += (s_max*s_min)*(uh_right - uh_left);
-      if(i==2) edgeflux[i] += (s_max*s_min)*(vh_right - vh_left);
+      if((i==1)) edgeflux[i] += (s_max*s_min)*(uh_right - uh_left);
+      //if((i==1)) edgeflux[i] += (s_max*s_min)*(uh_right/(h_right+1.0e-08)*hre - uh_left/(h_left+1.0e-08)*hle);
+      if((i==2)) edgeflux[i] += (s_max*s_min)*(vh_right - vh_left);
+      //if((i==2)) edgeflux[i] += (s_max*s_min)*(vh_right/(h_right+1.0e-08)*hre - vh_left/(h_left+1.0e-08)*hle);
 
       edgeflux[i] *= inverse_denominator;
@@ -538 +540 @@
     for(k=0; k<number_of_elements; k++){
         hc = max(stage_centroid_values[k] - bed_centroid_values[k],0.);
-        stage_explicit_update[k]=0.;
-        xmom_explicit_update[k]=0.;
-        ymom_explicit_update[k]=0.;
+        //stage_explicit_update[k]=0.;
+        //xmom_explicit_update[k]=0.;
+        //ymom_explicit_update[k]=0.;
 
         for(i=0;i<3;i++){
@@ -751 +753 @@
   double dx1, dx2, dy1, dy2, dxv0, dxv1, dxv2, dyv0, dyv1, dyv2, dq0, dq1, dq2, area2, inv_area2, dpth,momnorm;
   double dqv[3], qmin, qmax, hmin, hmax, bedmax,bedmin, stagemin;
-  double hc, h0, h1, h2, beta_tmp, hfactor, xtmp, ytmp, weight, tmp;
+  double hc, h0, h1, h2, beta_tmp, hfactor, xtmp, ytmp, weight, tmp, v0, v1, v2, vmax,vmin, vc;
   double dk, dv0, dv1, dv2, de[3], demin, dcmax, r0scale, vel_norm, l1, l2, a_tmp, b_tmp, c_tmp,d_tmp;
   
@@ -845 +847 @@
   // condition) set its momentum to zero too. Aim is to prevent local 'pits'
   // with a non-zero water surface gradient building up momentum over time
-  for (k=0; k<number_of_elements;k++){
-
-      k3=k*3;
-      k0 = surrogate_neighbours[k3];
-      k1 = surrogate_neighbours[k3 + 1];
-      k2 = surrogate_neighbours[k3 + 2];
-
-      if((height_centroid_values[k0] < minimum_allowed_height | k0==k) &
-         (height_centroid_values[k1] < minimum_allowed_height | k1==k) &
-         (height_centroid_values[k2] < minimum_allowed_height | k2==k)){
-              xmom_centroid_store[k] = 0.;
-              xmom_centroid_values[k] = 0.;
-              ymom_centroid_store[k] = 0.;
-              ymom_centroid_values[k] = 0.;
-
-      }
-  }
+  //for (k=0; k<number_of_elements;k++){
+  //
+  //    k3=k*3;
+  //    k0 = surrogate_neighbours[k3];
+  //    k1 = surrogate_neighbours[k3 + 1];
+  //    k2 = surrogate_neighbours[k3 + 2];
+  //
+  //    if((height_centroid_values[k0] < minimum_allowed_height | k0==k) &
+  //       (height_centroid_values[k1] < minimum_allowed_height | k1==k) &
+  //       (height_centroid_values[k2] < minimum_allowed_height | k2==k)){
+  //            xmom_centroid_store[k] = 0.;
+  //            xmom_centroid_values[k] = 0.;
+  //            ymom_centroid_store[k] = 0.;
+  //            ymom_centroid_values[k] = 0.;
+
+  //    }
+  //}
 
   // Begin extrapolation routine
@@ -1052 +1054 @@
       hmin = min(min(h0, min(h1, h2)), hc);
       hmax = max(max(h0, max(h1, h2)), hc);
+
+      //vc = sqrt(xmom_centroid_values[k]*xmom_centroid_values[k]+
+      //          ymom_centroid_values[k]*ymom_centroid_values[k]);
+      //v0 = sqrt(xmom_centroid_values[k0]*xmom_centroid_values[k0]+
+      //          ymom_centroid_values[k0]*ymom_centroid_values[k0]);
+      //v1 = sqrt(xmom_centroid_values[k1]*xmom_centroid_values[k1]+
+      //          ymom_centroid_values[k1]*ymom_centroid_values[k1]);
+      //v2 = sqrt(xmom_centroid_values[k2]*xmom_centroid_values[k2]+
+      //          ymom_centroid_values[k2]*ymom_centroid_values[k2]);
+      //vmax=max(max(v0*h0,max(v1*h1,v2*h2)),vc*hc);
+      //vmin=min(min(v0*h0,min(v1*h1,v2*h2)),vc*hc);
 
       // Look for strong changes in cell depth as an indicator of near-wet-dry
@@ -1075 +1088 @@
       // Set hfactor to 0 smoothly as hmin --> minumum_allowed_height 
       hfactor=min( 1.2*max(hmin-minimum_allowed_height,0.)/(max(hmin,0.)+1.*minimum_allowed_height), hfactor);
+
+      //hfactor=min(1.5*max(vmin,0.)/(max(vmax,0.)+1.0e-06), hfactor);
       //hfactor=min( 1.2*max(hc-minimum_allowed_height,0.)/(max(hc,0.)+1.*minimum_allowed_height), hfactor);
 
@@ -1252 +1267 @@
 
       // Limit the gradient
-      limit_gradient(dqv, qmin, qmax, beta_tmp);
+       limit_gradient(dqv, qmin, qmax, beta_tmp);
       
       for (i=0;i<3;i++)

trunk/anuga_work/development/gareth/tests/channel_floodplain/channel_floodplain1.py

@@ -88 +88 @@
 domain.set_name('channel_floodplain1') # Output name
 #domain.set_store_vertices_uniquely(True)
+domain.set_flow_algorithm('DE1')
 #domain.use_edge_limiter=False
 #domain.extrapolate_velocity_second_order=False

trunk/anuga_work/development/gareth/tests/channel_floodplain/plotme.py

@@ -1 +1 @@
 #import util
-from bal_and import plot_utils as util
+#from bal_and import plot_utils as util
+from anuga.utilities import plot_utils as util
 from matplotlib import pyplot as pyplot
 #import pylab
@@ -7 +8 @@
 #index=40
 index=570
+#index=225
 
 #p2 = util.get_output('channel_floodplain1_bal_dev.sww', minimum_allowed_height=0.01)
@@ -13 +15 @@
 
 #p2 = util.get_output('channel_floodplain1_bal_dev_lowvisc.sww', 0.01)
-p2 = util.get_output('channel_floodplain1.sww', 0.01)
+p2 = util.get_output('channel_floodplain1.sww', 0.001)
 p=util.get_centroids(p2, velocity_extrapolation=True)
 

trunk/anuga_work/development/gareth/tests/runup_sinusoid/runup_sinusoid.py

@@ -28 +28 @@
 #domain.set_store_centroids(True)
 #domain.set_timestepping_method('euler')
-domain.set_flow_algorithm('DE1')
+#domain.set_flow_algorithm('DE1')
+domain.set_flow_algorithm('tsunami')
 domain.set_quantities_to_be_stored({'stage': 2, 'xmomentum': 2, 
          'ymomentum': 2, 'elevation': 2})
-domain.minimum_storable_height=1.0e-03
+#domain.minimum_storable_height=1.0e-03
 domain.set_store_vertices_uniquely()
 #------------------

trunk/anuga_work/development/gareth/tests/runup_sinusoid/runup_sinusoidplot.py

@@ -5 +5 @@
 #---------------
 #import util
-from bal_and import plot_utils as util
+#from bal_and import plot_utils as util
+from anuga.utilities import plot_utils as util
 from matplotlib import pyplot as pyplot
 import numpy

trunk/anuga_work/development/gareth/tests/shallow_steep_slope/channel_SU_2plot.py

@@ -9 +9 @@
 import scipy
 from matplotlib import pyplot as pyplot
-#from anuga.utilities import plot_utils as util
-from bal_and import plot_utils as util
+from anuga.utilities import plot_utils as util
+#from bal_and import plot_utils as util
 #--------------
 # Get variables
@@ -61 +61 @@
     pyplot.ylabel('depth (m)')
 
-pyplot.title('Flow depth at 400s -- there should be a central region with steady uniform flow ' )
+pyplot.title('Flow depth at end of simulation -- there should be a central region with steady uniform flow ' )
 pyplot.savefig('final_depth_v2.png')
 
@@ -83 +83 @@
 v=(p2.x==p2.x[tmp])
 
+ll=p2.vel.shape[0]-1
 #--------------------------------------------------
 # Check if y velocity is zero everywhere (it should be)
 #--------------------------------------------------
 pyplot.clf()
-pyplot.plot(p2.y[v],p2.yvel[100,v],'o', label='Numerical')
+pyplot.plot(p2.y[v],p2.yvel[ll,v],'o', label='Numerical')
 pyplot.plot((0,100),(0.0, 0.0),label='Analytical')
 pyplot.xlabel('Distance along the line x==50 (across the slope) (m)')
@@ -97 +98 @@
 
 pyplot.clf()
-pyplot.plot(p2.y[v],p2.xvel[100,v], 'o', label='Numerical')
+pyplot.plot(p2.y[v],p2.xvel[ll,v], 'o', label='Numerical')
 pyplot.plot((0,100),(uana,uana),label='Analytical')
 #pyplot.ylim([1.0,3.0])

trunk/anuga_work/development/gareth/tests/shallow_steep_slope/channel_SU_sparse.py

@@ -37 +37 @@
 #------------------------------------------------------------------------------
 def topography(x, y):
-        return -x/10.  + 1.*(numpy.sin(x/10.) +abs(y-50.)/10.) -0.*(x>80.) # linear bed slope
+        return -x/10.  #+ 1.*(numpy.sin(x/10.) +abs(y-50.)/10.) -0.*(x>80.) # linear bed slope
 
 def stagetopo(x,y):

trunk/anuga_work/development/gareth/tests/wave_runup/plotme.py

@@ -66 +66 @@
 for i in range(len(model_shore_x)):
     # Compute index of shoreline
-    vtmp = p2.stage[i,:]>p2.elev+shoreline_depth_thresh1
+    #vtmp = p2.stage[i,:]>p2.elev+shoreline_depth_thresh1
+    #model_shore_x[i]=p2.x[vtmp].min()
     
-    model_shore_x[i]=p2.x[vtmp].min()
+    vtmp = p2.stage[i,:]<p2.elev+shoreline_depth_thresh1
+    model_shore_x[i]=p2.x[vtmp].max()
+    vtmp2=(p2.x>model_shore_x[i])
+    model_shore_x[i]=p2.x[vtmp2].min()
     #v_ind=(p2.stage[i,vtmp]-p2.elev[vtmp]).argmin()
     # Store x and xvel
@@ -76 +80 @@
     # Extract shoreline velocity. It is tricky, to avoid getting
     # zero velocity points -- might be a better way to do it?
-    vtmp2= (p2.stage[i,:]>p2.elev+shoreline_depth_thresh1)*\
-           (p2.stage[i,:]<p2.elev+shoreline_depth_thresh2)
-    mloc=abs(p2.xvel[i,vtmp2]).argmax()
+    #vtmp2= (p2.stage[i,:]>p2.elev+shoreline_depth_thresh1)*\
+    #       (p2.stage[i,:]<p2.elev+shoreline_depth_thresh2)
+    #mloc=abs(p2.xvel[i,vtmp2]).argmax()
+
     #print abs(p2.xvel[i,vtmp2]).max(), mloc
-    model_shore_u[i]=p2.xvel[i,vtmp2][mloc]
+    #model_shore_u[i]=p2.xvel[i,vtmp2][mloc]
     #vtmp2 = (p2.stage[i,:]>p2.elev+0.01)
     #mloc= (vtmp2*(p2.stage[i,:]-p2.elev)).argmin()
     #model_shore_u[i] = p2.xvel[i,mloc]
+    
+    vtmp2=vtmp
+    mloc=model_shore_x[i]
+    model_shore_u[i] = p2.xvel[i,p2.x==model_shore_x[i]].max()
 
 pyplot.plot(p2.time, model_shore_x-200.,'-o')

trunk/anuga_work/development/gareth/tests/wave_runup/runup.py

@@ -10 +10 @@
 
 from math import sin, pi, exp
-from bal_and import *
-#from anuga import Domain as Domain
+#from bal_and import *
+from anuga import Domain as Domain
 
 #from anuga_tsunami import *
@@ -32 +32 @@
 domain.set_datadir('.')                          # Use current folder
 #domain.set_quantities_to_be_stored({'stage': 2, 'xmomentum': 2, 'ymomentum': 2, 'elevation': 1})
-#domain.set_flow_algorithm('DE1')
+domain.set_flow_algorithm('DE1')
 #------------------
 # Define topography