Index: /inundation/ga/storm_surge/pyvolution/HUSK.txt
===================================================================
--- /inundation/ga/storm_surge/pyvolution/HUSK.txt	(revision 268)
+++ /inundation/ga/storm_surge/pyvolution/HUSK.txt	(revision 269)
@@ -18,4 +18,5 @@
 
 Don't use default conserved qs - hardwire into shallow_water.py
+Even hardwire pertinent arrays from cons quan into shallow_water and use directly with flux etc
 
 
Index: /inundation/ga/storm_surge/pyvolution/domain.py
===================================================================
--- /inundation/ga/storm_surge/pyvolution/domain.py	(revision 268)
+++ /inundation/ga/storm_surge/pyvolution/domain.py	(revision 269)
@@ -91,17 +91,4 @@
 	self.filename = 'domain'
 	self.checkpoint = False
-
-        #Realtime visualisation
-        self.visualise = False
-        
-        #Stored output
-        self.store=False
-	self.format = 'dat'    
-        self.smooth = True
-
-        #Reduction operation for get_vertex_values              
-        #from pytools.stats import mean
-        #self.reduction = mean 
-        self.reduction = min  #Looks better near steep slopes
         
 
@@ -277,16 +264,6 @@
 	if self.checkpoint is True:
             self.goto_latest_checkpoint()
-            
-
-        #Store model data, e.g. for visualisation    
-        if self.store is True and self.time == 0.0:
-            self.store_bathymetry()        
-            self.store_conserved_quantities()	
-
-        if self.visualise is True and self.time == 0.0:
-            import realtime_visualisation as visualise
-            visualise.create_surface(self)
-        
-	
+
+            
         yield(self.time)  #Yield initial values
 	
@@ -294,23 +271,7 @@
             #Update boundary values
             self.update_boundary()            
-	    
-	    #print
-	    #for name in self.conserved_quantities:
-            #    Q = self.quantities[name] 
-	    #    #print 'Vertices (%s):' %name, Q.vertex_values[:]
-            #    print 'B_val (%s):' %name, Q.boundary_values[:]
-
-	    
-	    #print
-	    #for name in self.conserved_quantities:
-            #    Q = self.quantities[name] 
-	    #    print 'Edges (%s):' %name, Q.edge_values[:4]	    	    
-
-            #Compute all fluxes and timestep suitable for all volumes
+
+	    #Compute fluxes across each element edge 
             self.compute_fluxes()
-	    ##print
-	    ##for name in self.conserved_quantities:
-            ##    Q = self.quantities[name] 
-	    ##    print 'EU:', Q.explicit_update[:4]	    
 
             #Update timestep to fit yieldstep and finaltime
@@ -320,25 +281,7 @@
             self.update_conserved_quantities()            
 
-            #print
-	    #print 'Centroids'
-            #print self.quantities['level'].centroid_values[1:4],\
-            #      self.quantities['level'].centroid_values[13]
-            #print self.quantities['xmomentum'].centroid_values[1:4],\
-            #      self.quantities['xmomentum'].centroid_values[13]	    
-            #print self.quantities['ymomentum'].centroid_values[1:4],\
-            #      self.quantities['ymomentum'].centroid_values[13]
-		  					    
-	    
             #Update vertex and edge values
             self.distribute_to_vertices_and_edges()
 	    
-	    #print
-	    #for name in self.conserved_quantities:
-            #    Q = self.quantities[name] 
-	    #    print 'Vertices (%s):' %name, Q.vertex_values[1:4], Q.vertex_values[13]
-		 	    	    	    
-	    
-	    
-                
             #Update time    
             self.time += self.timestep
@@ -348,16 +291,10 @@
                 self.number_of_first_order_steps += 1        
 
-		
-            #print 'Time=', self.time, self.timestep	
-            #print self.quantities['level'].centroid_values[:4]
-            #print self.quantities['xmomentum'].centroid_values[:4]
-            #print self.quantities['ymomentum'].centroid_values[:4]				
-	    #print
-		
             #Yield results
             if finaltime is not None and abs(self.time - finaltime) < epsilon:
-                # Yield final time and allow inspection of domain
+                # Yield final time and stop
                 yield(self.time)
                 break
+
                 
             if abs(self.yieldtime - yieldstep) < epsilon: 
@@ -368,12 +305,5 @@
                     self.delete_old_checkpoints()
                     
-                #Store model data, e.g. for subsequent visualisation    
-                if self.store is True:
-                    self.store_conserved_quantities()
-
-                #Real time viz
-                if self.visualise is True:
-                    visualise.update(self)
-                    
+                #Pass control on to outer loop for more specific actions    
                 yield(self.time)  
 
Index: /inundation/ga/storm_surge/pyvolution/shallow_water.py
===================================================================
--- /inundation/ga/storm_surge/pyvolution/shallow_water.py	(revision 268)
+++ /inundation/ga/storm_surge/pyvolution/shallow_water.py	(revision 269)
@@ -36,4 +36,19 @@
         self.forcing_terms.append(gravity)
         self.forcing_terms.append(manning_friction)
+
+        #Realtime visualisation
+        self.visualise = False
+
+
+        #Stored output
+        self.store=False
+	self.format = 'dat'    
+        self.smooth = True
+
+        #Reduction operation for get_vertex_values              
+        #from pytools.stats import mean
+        #self.reduction = mean 
+        self.reduction = min  #Looks better near steep slopes
+        
 
         #Establish shortcuts to relevant quantities (for efficiency)
@@ -84,4 +99,30 @@
         #(either from this module or C-extension)        
         distribute_to_vertices_and_edges(self)
+
+
+    def evolve(self, yieldstep = None, finaltime = None):
+        #Call basic machinery from parent class
+
+
+        if self.visualise is True and self.time == 0.0:
+            import realtime_visualisation as visualise
+            visualise.create_surface(self)
+
+        #Store model data, e.g. for visualisation    
+        if self.store is True and self.time == 0.0:
+            self.store_bathymetry()        
+            ###self.store_conserved_quantities()	
+                    
+        for t in Generic_domain.evolve(self, yieldstep, finaltime):
+            #Real time viz
+            if self.visualise is True:
+                visualise.update(self)
+
+            #Store model data, e.g. for subsequent visualisation    
+            if self.store is True:
+                self.store_conserved_quantities()
+
+            #Pass control on to outer loop for more specific actions    
+            yield(t)
         
 
Index: /inundation/ga/storm_surge/pyvolution/show_balanced_limiters.py
===================================================================
--- /inundation/ga/storm_surge/pyvolution/show_balanced_limiters.py	(revision 268)
+++ /inundation/ga/storm_surge/pyvolution/show_balanced_limiters.py	(revision 269)
@@ -31,4 +31,5 @@
                                          origin=(-0.07, 0))
 
+print 'Number of elements', len(vertices)
 #Create shallow water domain
 domain = Domain(points, vertices, boundary)
@@ -80,7 +81,8 @@
 
 #Evolve
-for t in domain.evolve(yieldstep = 0.1, finaltime = 50):
+for t in domain.evolve(yieldstep = 0.1, finaltime = 30):
     domain.write_time()
     #print domain.quantities['level'].centroid_values[:6]
 print 'Done'    
     
+raw_input()
