Index: anuga_core/source/anuga/culvert_flows/culvert_class.py
===================================================================
--- anuga_core/source/anuga/culvert_flows/culvert_class.py	(revision 6928)
+++ anuga_core/source/anuga/culvert_flows/culvert_class.py	(revision 6948)
@@ -15,5 +15,5 @@
 
 import Numeric as num
-
+from math import sqrt
 
 class Below_interval(Exception): pass 
@@ -137,5 +137,5 @@
                  manning=None,          # Mannings Roughness for Culvert 
                  sum_loss=None,
-                 use_velocity_head=False, 
+                 use_velocity_head=False, # FIXME(Ole): Get rid of - always True
                  use_momentum_jet=False, # FIXME(Ole): Not yet implemented
                  label=None,
@@ -464,5 +464,5 @@
             log_filename = self.log_filename
             
-        # Compute stage and energy at the 
+        # Compute stage, energy and velocity at the 
         # enquiry points at each end of the culvert
         openings = self.openings
@@ -474,18 +474,21 @@
             depth = h = stage-opening.elevation
                                                            
-                                               
+            
+            # Get velocity                                 
+            xmomentum = dq['xmomentum'].get_values(location='centroids',
+                                                   indices=[idx])[0]
+            ymomentum = dq['xmomentum'].get_values(location='centroids',
+                                                   indices=[idx])[0]
+            
+            if h > minimum_allowed_height:
+                u = xmomentum/(h + velocity_protection/h)
+                v = ymomentum/(h + velocity_protection/h)
+            else:
+                u = v = 0.0
+                
+            v_squared = u*u + v*v
+            
             if self.use_velocity_head is True:
-                xmomentum = dq['xmomentum'].get_values(location='centroids',
-                                                       indices=[idx])[0]
-                ymomentum = dq['xmomentum'].get_values(location='centroids',
-                                                       indices=[idx])[0]
-            
-                if h > minimum_allowed_height:
-                    u = xmomentum/(h + velocity_protection/h)
-                    v = ymomentum/(h + velocity_protection/h)
-                else:
-                    u = v = 0.0
-                
-                velocity_head = 0.5*(u*u + v*v)/g    
+                velocity_head = 0.5*v_squared/g    
             else:
                 velocity_head = 0.0
@@ -495,4 +498,5 @@
             opening.stage = stage
             opening.depth = depth
+            opening.velocity = sqrt(v_squared)
             
 
@@ -583,4 +587,6 @@
                     self.culvert_routine(inlet.depth,
                                          outlet.depth,
+                                         inlet.velocity,
+                                         outlet.velocity,
                                          inlet.specific_energy, 
                                          delta_total_energy, 
Index: anuga_core/source/anuga/culvert_flows/culvert_routines.py
===================================================================
--- anuga_core/source/anuga/culvert_flows/culvert_routines.py	(revision 6928)
+++ anuga_core/source/anuga/culvert_flows/culvert_routines.py	(revision 6948)
@@ -19,4 +19,6 @@
 def boyd_generalised_culvert_model(inlet_depth, 
                                    outlet_depth,
+                                   inlet_velocity,
+                                   outlet_velocity,
                                    inlet_specific_energy, 
                                    delta_total_energy, 
Index: anuga_core/source/anuga/culvert_flows/test_culvert_routines.py
===================================================================
--- anuga_core/source/anuga/culvert_flows/test_culvert_routines.py	(revision 6928)
+++ anuga_core/source/anuga/culvert_flows/test_culvert_routines.py	(revision 6948)
@@ -31,4 +31,7 @@
         inlet_depth=2.0
         outlet_depth=0.0
+        
+        inlet_velocity=0.0,
+        outlet_velocity=0.0,        
 
         culvert_length=4.0
@@ -49,4 +52,6 @@
         Q, v, d = boyd_generalised_culvert_model(inlet_depth, 
                                                  outlet_depth,
+                                                 inlet_velocity,
+                                                 outlet_velocity,
                                                  inlet_specific_energy, 
                                                  delta_total_energy, 
@@ -80,4 +85,7 @@
         outlet_depth=0.0
 
+        inlet_velocity=0.0,
+        outlet_velocity=0.0,                
+        
         culvert_length=4.0
         culvert_width=1.2
@@ -97,4 +105,6 @@
         Q, v, d = boyd_generalised_culvert_model(inlet_depth, 
                                                  outlet_depth,
+                                                 inlet_velocity,
+                                                 outlet_velocity,
                                                  inlet_specific_energy, 
                                                  delta_total_energy, 
@@ -189,4 +199,6 @@
         Q, v, d = boyd_generalised_culvert_model(inlet_depth, 
                                                  outlet_depth,
+                                                 inlet_velocity,
+                                                 outlet_velocity,
                                                  inlet_specific_energy, 
                                                  delta_total_energy, 
