Changeset 8893

Timestamp:

Jun 10, 2013, 11:11:02 AM (12 years ago)

Author:

davies

Message:

Updates to bal_dev

Location:

trunk/anuga_work/development/gareth

Files:

• experimental/balanced_dev/swb2_domain_ext.c (modified) (8 diffs)
• tests/merimbula_steve/run_sw_merimbula.py (modified) (1 diff)
• tests/shallow_steep_slope/channel_SU_sparse.py (modified) (4 diffs)

trunk/anuga_work/development/gareth/experimental/balanced_dev/swb2_domain_ext.c ¶

@@ -948 +948 @@
   ymom_centroid_store = malloc(number_of_elements*sizeof(double));
   stage_centroid_store = malloc(number_of_elements*sizeof(double));
-  //min_elevation_edgevalue = malloc(number_of_elements*sizeof(double));
+  min_elevation_edgevalue = malloc(number_of_elements*sizeof(double));
   max_elevation_edgevalue = malloc(number_of_elements*sizeof(double));
   cell_wetness_scale = malloc(number_of_elements*sizeof(double));
@@ -965 +965 @@
           ymom_centroid_values[k] = ymom_centroid_values[k]/dk;
 
-          //min_elevation_edgevalue[k] = min(elevation_edge_values[3*k], 
-          //                                 min(elevation_edge_values[3*k+1],
-          //                                     elevation_edge_values[3*k+2]));
+          min_elevation_edgevalue[k] = min(elevation_edge_values[3*k], 
+                                           min(elevation_edge_values[3*k+1],
+                                               elevation_edge_values[3*k+2]));
           max_elevation_edgevalue[k] = max(elevation_edge_values[3*k], 
                                            max(elevation_edge_values[3*k+1],
@@ -1005 +1005 @@
   }
 
-
   // Alternative 'PROTECT' step
   for(k=0; k<number_of_elements;k++){ 
-    if(count_wet_neighbours[k]<=3){
-        bedmax = max(elevation_centroid_values[k], 
-                     max(elevation_vertex_values[3*k],
+    //if(count_wet_neighbours[k]<=3){
+    if(1==1){
+        bedmax = max(elevation_vertex_values[3*k],
                          max(elevation_vertex_values[3*k+1], 
-                             elevation_vertex_values[3*k+2])));
-        if(cell_wetness_scale[k]==0. ||
-            stage_centroid_values[k] - elevation_centroid_values[k] < 0.05*(bedmax-elevation_centroid_values[k])){
+                             elevation_vertex_values[3*k+2]));
+        if((cell_wetness_scale[k]==0. ) ||
+            ((stage_centroid_values[k] - elevation_centroid_values[k] < 0.05*(bedmax-elevation_centroid_values[k])) &
+             count_wet_neighbours[k]<=3) ){
             xmom_centroid_store[k]=0.;
             xmom_centroid_values[k]=0.;
@@ -1033 +1033 @@
         // For partially wet cells, we now know that the edges of neighbouring
         // fully wet cells have been defined
-        if( cell_wetness_scale[k]==(1.0*k_wetdry) ){
+        if( cell_wetness_scale[k]==1.0*(1-k_wetdry) ){
           continue;
         }
@@ -1144 +1144 @@
                            max(elevation_centroid_values[k1], 
                                elevation_centroid_values[k2])));
-          //bedmin = min(elevation_centroid_values[k], 
-          //             min(elevation_centroid_values[k0],
-          //                 min(elevation_centroid_values[k1], 
-          //                     elevation_centroid_values[k2])));
+          bedmin = min(elevation_centroid_values[k], 
+                       min(elevation_centroid_values[k0],
+                           min(elevation_centroid_values[k1], 
+                               elevation_centroid_values[k2])));
           stagemin = min(max(stage_centroid_values[k], elevation_centroid_values[k]), 
                          min(max(stage_centroid_values[k0], elevation_centroid_values[k0]),
@@ -1280 +1280 @@
               //    //    //}else{
               //    //        // Constant depth extrapolation
-              //        //if(stage_centroid_values[k]<stage_centroid_values[ktmp]){
-              //          //  stage_edge_values[k3]= stage_centroid_values[k]-elevation_centroid_values[k]+elevation_edge_values[k3];
-              //          //  stage_edge_values[k3+1]= stage_centroid_values[k] -elevation_centroid_values[k]+elevation_edge_values[k3+1];
-              //          //  stage_edge_values[k3+2]= stage_centroid_values[k] -elevation_centroid_values[k]+elevation_edge_values[k3+2];
-              //        stage_edge_values[k3] = stage_centroid_values[k];
-              //        stage_edge_values[k3+1] = stage_centroid_values[k];
-              //        stage_edge_values[k3+2] = stage_centroid_values[k];
-              //        //}else{
+                    if(stage_centroid_values[k]<stage_centroid_values[ktmp]){
+                      // Constant stage extrapolation reduces the issue of having the bed-slope term playing up on nearly dry cells
+                      // stage_edge_values[k3]= stage_centroid_values[k]-elevation_centroid_values[k]+elevation_edge_values[k3];
+                      // stage_edge_values[k3+1]= stage_centroid_values[k] -elevation_centroid_values[k]+elevation_edge_values[k3+1];
+                      // stage_edge_values[k3+2]= stage_centroid_values[k] -elevation_centroid_values[k]+elevation_edge_values[k3+2];
+                      stage_edge_values[k3] = stage_centroid_values[k];
+                      stage_edge_values[k3+1] = stage_centroid_values[k];
+                      stage_edge_values[k3+2] = stage_centroid_values[k];
+                    }else{
               //           // stage_edge_values[k3] = stage_centroid_values[ktmp];
               //           // stage_edge_values[k3+1] = stage_centroid_values[ktmp];
               //           // stage_edge_values[k3+2] = stage_centroid_values[ktmp];
-                     stage_edge_values[k3]= max(stage_centroid_values[k]-elevation_centroid_values[k]+elevation_edge_values[k3], stage_centroid_values[ktmp]);
-                     stage_edge_values[k3+1]= max(stage_centroid_values[k] -elevation_centroid_values[k]+elevation_edge_values[k3+1],stage_centroid_values[ktmp]);
-                     stage_edge_values[k3+2]= max(stage_centroid_values[k] -elevation_centroid_values[k]+elevation_edge_values[k3+2], stage_centroid_values[ktmp]);
-              //        //}
+                     //stage_edge_values[k3]= max(stage_centroid_values[k]-elevation_centroid_values[k]+elevation_edge_values[k3], stage_centroid_values[ktmp]);
+                     //stage_edge_values[k3+1]= max(stage_centroid_values[k] -elevation_centroid_values[k]+elevation_edge_values[k3+1],stage_centroid_values[ktmp]);
+                     //stage_edge_values[k3+2]= max(stage_centroid_values[k] -elevation_centroid_values[k]+elevation_edge_values[k3+2], stage_centroid_values[ktmp]);
+                     ktmp = 3*ktmp+neighbour_edges[k3+ii];
+                     stage_edge_values[k3]= max(stage_centroid_values[k]-elevation_centroid_values[k]+elevation_edge_values[k3], stage_edge_values[ktmp]);
+                     stage_edge_values[k3+1]= max(stage_centroid_values[k] -elevation_centroid_values[k]+elevation_edge_values[k3+1],stage_edge_values[ktmp]);
+                     stage_edge_values[k3+2]= max(stage_centroid_values[k] -elevation_centroid_values[k]+elevation_edge_values[k3+2], stage_edge_values[ktmp]);
+                    }
 
               //    }else{
@@ -1414 +1419 @@
           //if( (area2>0) ){ 
           //if(count_wet_neighbours[k]>0){
-          //if(bedmax<=stagemin){
+          //if(min_elevation_edgevalue[k]<=stagemin ){
+          //if(bedmax<=stagemin ){
              limit_gradient(dqv, qmin, qmax, beta_tmp);
           //}else{
@@ -1823 +1829 @@
   free(ymom_centroid_store);
   free(stage_centroid_store);
-  //free(min_elevation_edgevalue);
+  free(min_elevation_edgevalue);
   free(max_elevation_edgevalue);
   free(cell_wetness_scale);

trunk/anuga_work/development/gareth/tests/merimbula_steve/run_sw_merimbula.py ¶

@@ -98 +98 @@
 #domain.set_beta(1.5)
 domain.set_name('merimbula')
-domain.minimum_allowed_height=0.01
+#domain.minimum_allowed_height=0.01
 
 #------------------------------------------------------------------------------

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

@@ -22 +22 @@
 from anuga.shallow_water.shallow_water_domain import Domain as Domain
 from balanced_dev import *
-from balanced_dev import Domain as Domain
+#from balanced_dev import Domain as Domain
 #from anuga_tsunami import *
 #from anuga_tsunami import Domain as Domain
@@ -36 +36 @@
 #------------------------------------------------------------------------------
 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):
@@ -44 +44 @@
 
 line1=[ [9.,0.], [9., 100.] ]
-Qin=2.
+Qin=0.5
 Inlet_operator(domain, line1,Qin)
 
@@ -73 +73 @@
 # Evolve system through time
 #------------------------------------------------------------------------------
-for t in domain.evolve(yieldstep=4.0, finaltime=400.0):
+for t in domain.evolve(yieldstep=4.0, finaltime=800.0):
     print domain.timestepping_statistics()
-    print domain.boundary_flux_integral
+    #print domain.boundary_flux_integral
     print (domain.areas*(domain.quantities['stage'].centroid_values - domain.quantities['elevation'].centroid_values)).sum()
     #print domain.quantities['stage'].centroid_values[myindex] - domain.quantities['elevation'].centroid_values[myindex]