Changeset 5175 for anuga_core/source/anuga/shallow_water/shallow_water_ext.c

Timestamp:

Mar 19, 2008, 8:02:05 PM (16 years ago)

Author:

steve

Message:

Ole and I seem to have got edge limiteri working.
Set domain.use_edge_limiter = True
to test the method

File:

• anuga_core/source/anuga/shallow_water/shallow_water_ext.c (modified) (5 diffs)

anuga_core/source/anuga/shallow_water/shallow_water_ext.c

@@ -589 +589 @@
     // Calculate minimal depth across all three vertices
     hmin = min(hv[0], min(hv[1], hv[2]));
+
+    //if (hmin < 0.0 ) {
+    //  printf("hmin = %f\n",hmin);
+    //}
     
 
@@ -647 +651 @@
       }
 
-      
+      //printf("alpha=%f, tight_slope_limiters=%d\n", alpha, tight_slope_limiters);
+
       /*      
       // Experimental code for controlling velocities at vertices.
@@ -704 +709 @@
 
         // Update momentum at vertices
-        if (tight_slope_limiters == 1) {        
-          // FIXME(Ole): Here's what I think (as of 17 Nov 2007) 
+
+
+        // Update momentum at vertices
+        // Update momentum as a linear combination of
+        // xmomc and ymomc (shallow) and momentum
+        // from extrapolator xmomv and ymomv (deep).
+
+
+        //xmomv[k3+i] = (1-alpha)*xmomc[k] + alpha*xmomv[k3+i];
+        //ymomv[k3+i] = (1-alpha)*ymomc[k] + alpha*ymomv[k3+i];
+
+        if (tight_slope_limiters == 1) {
+          // FIXME(Ole): Here's what I think (as of 17 Nov 2007)
           // we need to do. Simple and efficient:
         
@@ -721 +737 @@
           hv[i] = wv[k3+i] - zv[k3+i]; // Recompute (balanced) vertex depth
           xmomv[k3+i] = uc*hv[i];
-          ymomv[k3+i] = vc*hv[i];       
+          ymomv[k3+i] = vc*hv[i];
         } else {
           // Update momentum as a linear combination of
           // xmomc and ymomc (shallow) and momentum
           // from extrapolator xmomv and ymomv (deep).
-          // FIXME (Ole): Is this really needed? Could we use the above 
+          // FIXME (Ole): Is this really needed? Could we use the above
           // instead?
           
@@ -1189 +1205 @@
 
       hfactor = 0.0;
-      if (hmin > 0.1 ) {
-          hfactor = (hmin-0.1)/(hmin+0.4);
+      if (hmin > 0.001 ) {
+          hfactor = (hmin-0.001)/(hmin+0.004);
       }