Changeset 273

Timestamp:

Sep 6, 2004, 11:43:36 AM (21 years ago)

Author:

ole

Message:

Implemented protect in c and noted improvement.
Pyvolution is now faster than the previous version.

Location:

inundation/ga/storm_surge/pyvolution

Files:

• config.py (modified) (1 diff)
• shallow_water.py (modified) (2 diffs)
• shallow_water_ext.c (modified) (4 diffs)

inundation/ga/storm_surge/pyvolution/config.py

@@ -59 +59 @@
 
 use_psyco = True  #Use psyco optimisations
-use_psyco = False  #Do not use psyco optimisations
+#use_psyco = False  #Do not use psyco optimisations
 
 

inundation/ga/storm_surge/pyvolution/shallow_water.py

@@ -426 +426 @@
             #Control momentum
             xmomc[k] = ymomc[k] = 0.0
+            
+
+
+def protect_against_infintesimal_and_negative_heights_c(domain):
+    """Protect against infinitesimal heights and associated high velocities
+    """
+    
+    #Shortcuts
+    wc = domain.quantities['level'].centroid_values
+    zc = domain.quantities['elevation'].centroid_values
+    xmomc = domain.quantities['xmomentum'].centroid_values
+    ymomc = domain.quantities['ymomentum'].centroid_values            
+
+    from shallow_water_ext import protect
+
+    protect(domain.minimum_allowed_height, wc, zc, xmomc, ymomc)
+    
+
 
 def extrapolate_first_order(domain):
@@ -954 +972 @@
     gravity = gravity_c
     manning_friction = manning_friction_c
-    balance_deep_and_shallow = balance_deep_and_shallow_c 
+    balance_deep_and_shallow = balance_deep_and_shallow_c
+    protect_against_infintesimal_and_negative_heights = protect_against_infintesimal_and_negative_heights_c
     
     #distribute_to_vertices_and_edges = distribute_to_vertices_and_edges_c

inundation/ga/storm_surge/pyvolution/shallow_water_ext.c

@@ -253 +253 @@
   return 0;
 }
+
+
+
+int _protect(int N, 
+             double minimum_allowed_height,        
+             double* wc,
+             double* zc, 
+             double* xmomc, 
+             double* ymomc) {
+  
+  int k; 
+  double hc;
+  
+  //Compute linear combination between constant levels and and
+  //levels parallel to the bed elevation.     
+  
+  for (k=0; k<N; k++) {
+    hc = wc[k] - zc[k];
+    if (hc < minimum_allowed_height) {
+      if (hc < 0.0) {
+        //Control level and height
+        wc[k] = zc[k];
+      }
+
+      //Control momentum
+      xmomc[k] = 0.0;
+      ymomc[k] = 0.0;
+    }
+  }
+  return 0;
+}
+
+
 
 ///////////////////////////////////////////////////////////////////
@@ -559 +592 @@
 
 
+PyObject *protect(PyObject *self, PyObject *args) {
+  //
+  //    protect(minimum_allowed_height, wc, zc, xmomc, ymomc)
+  
+
+  PyArrayObject 
+  *wc,            //Level at centroids
+  *zc,            //Elevation at centroids    
+  *xmomc,         //Momentums at centroids
+  *ymomc; 
+
+    
+  int N; 
+  double minimum_allowed_height;
+  
+  // Convert Python arguments to C  
+  if (!PyArg_ParseTuple(args, "dOOOO",  
+                        &minimum_allowed_height,
+                        &wc, &zc, &xmomc, &ymomc))
+    return NULL;
+
+  N = wc -> dimensions[0];
+    
+  _protect(N,
+           minimum_allowed_height,
+           (double*) wc -> data,
+           (double*) zc -> data, 
+           (double*) xmomc -> data, 
+           (double*) ymomc -> data);
+ 
+  return Py_BuildValue("");  
+}
+
+
+
 PyObject *balance_deep_and_shallow(PyObject *self, PyObject *args) {
   //
-  //    balance_deep_and_shallow(domain.number_of_elements,
-  //                             wc, zc, hc, wv, zv, hv,
+  //    balance_deep_and_shallow(wc, zc, hc, wv, zv, hv,
   //                             xmomc, ymomc, xmomv, ymomv)
   
@@ -607 +674 @@
 
 
-
-
 //////////////////////////////////////////      
 // Method table for python module
@@ -623 +688 @@
   {"balance_deep_and_shallow", balance_deep_and_shallow, 
    METH_VARARGS, "Print out"},          
+  {"protect", protect, METH_VARARGS | METH_KEYWORDS, "Print out"},   
   //{"distribute_to_vertices_and_edges", 
   // distribute_to_vertices_and_edges, METH_VARARGS},