Changeset 9118

Timestamp:

May 19, 2014, 9:21:40 AM (10 years ago)

Author:

davies

Message:

Small change to riverwall numerics [motivated by problem in a coarse channel example]

File:

• trunk/anuga_core/source/anuga/shallow_water/swDE1_domain_ext.c (modified) (3 diffs)

trunk/anuga_core/source/anuga/shallow_water/swDE1_domain_ext.c

@@ -280 +280 @@
 
     minhd=min(h_left, h_right);
-    maxhd=max(h_left,h_right);
+    maxhd=max(h_left, h_right);
     // 'Raw' weir discharge = weirScale*2/3*H*(2/3*g*H)**0.5
     rw=weirScale*2./3.*maxhd*sqrt(2./3.*g*maxhd);
@@ -480 +480 @@
                 // riverwall_elevation + riverwall_rowIndex
                 RiverWall_count+=1;
-                // Since there is a wall, use first order extrapolation for this edge
-                // This also makes more sense from the viewpoint of weir relations
-                ql[0]=stage_centroid_values[k];
-                ql[1]=xmom_centroid_values[k];
-                ql[2]=ymom_centroid_values[k];
-                hle=hc;
-                zl=zc;
-                if(n>=0){
-                    qr[0]=stage_centroid_values[n];
-                    qr[1]=xmom_centroid_values[n];
-                    qr[2]=ymom_centroid_values[n];
-                    hre=hc_n;
-                    zr = zc_n;
-                }
                 
                 // Set central bed to riverwall elevation
                 z_half=max(riverwall_elevation[RiverWall_count-1], max(zl, zr)) ;
+
+                if(min(ql[0], qr[0]) < z_half){
+                    // Since there is a wall blocking the flow connection, use first order extrapolation for this edge
+                    ql[0]=stage_centroid_values[k];
+                    ql[1]=xmom_centroid_values[k];
+                    ql[2]=ymom_centroid_values[k];
+                    hle=hc;
+                    zl=zc;
+
+                    if(n>=0){
+                      qr[0]=stage_centroid_values[n];
+                      qr[1]=xmom_centroid_values[n];
+                      qr[2]=ymom_centroid_values[n];
+                      hre=hc_n;
+                      zr = zc_n;
+                    }else{
+                      hre=hc;
+                      zr = zc;
+                    }
+                    // Re-set central bed to riverwall elevation
+                    z_half=max(riverwall_elevation[RiverWall_count-1], max(zl, zr)) ;
+                }
+                
 
             }
@@ -516 +525 @@
             if(edge_flux_type[ki]==1){
                 //printf("%e \n", z_half);
-                weir_height=riverwall_elevation[RiverWall_count-1]-min(zl,zr); // Reference weir height  
+                weir_height=max(riverwall_elevation[RiverWall_count-1]-min(zl,zr), 0.); // Reference weir height  
+                //weir_height=max(z_half-max(zl,zr), 0.); // Reference weir height  
 
                 // Get Qfactor index - multiply the idealised weir discharge by this constant factor