Index: /trunk/anuga_core/validation_tests/Tests/Experimental_data/dam_break_Yeh_Petroff/for_exp.txt
===================================================================
--- /trunk/anuga_core/validation_tests/Tests/Experimental_data/dam_break_Yeh_Petroff/for_exp.txt	(revision 8673)
+++ /trunk/anuga_core/validation_tests/Tests/Experimental_data/dam_break_Yeh_Petroff/for_exp.txt	(revision 8673)
@@ -0,0 +1,416 @@
+4.3E-19		-4.18E0
+1.14E-2		-2.40E0
+2.66E-2		-2.12E0
+4.34E-19	-5.82E-2
+2.28E-2		-1.18E0
+4.56E-2		-1.18E0
+1.90E-2		-7.14E-1
+3.80E-2		-7.16E-1
+4.56E-2		-1.55E-1
+6.08E-2		-6.20E-2
+7.98E-2		3.03E-2
+8.75E-2		3.11E-1
+9.89E-2		5.91E-1
+1.10E-1		1.22E-1
+1.10E-1		-3.46E-1
+1.33E-1		2.14E-1
+1.29E-1		-6.63E-2
+1.44E-1		-4.42E-1
+1.60E-1		1.19E-1
+1.79E-1		-6.94E-2
+1.98E-1		-2.58E-1
+2.13E-1		-3.52E-1
+2.17E-1		-6.33E-1
+2.21E-1		-9.14E-1
+2.28E-1		-1.29E0
+2.43E-1		-7.28E-1
+2.40E-1		-4.47E-1
+2.47E-1		-1.67E-1
+2.59E-1		-7.44E-2
+2.81E-1		-7.58E-2
+2.97E-1		-3.57E-1
+3.16E-1		-3.59E-1
+3.08E-1		1.61E-2
+2.97E-1		2.98E-1
+3.23E-1		-7.48E-2
+3.35E-1		2.02E-1
+3.42E-1		8.56E-1
+3.42E-1		1.60E0
+3.46E-1		2.54E0
+3.42E-1		4.32E0
+3.57E-1		4.60E0
+3.54E-1		7.97E0
+3.61E-1		7.97E0
+3.54E-1		1.13E1
+3.69E-1		1.43E1
+3.65E-1		1.67E1
+3.73E-1		1.75E1
+3.73E-1		2.13E1
+3.69E-1		2.51E1
+3.84E-1		2.62E1
+3.80E-1		2.89E1
+3.92E-1		3.09E1
+3.95E-1		3.32E1
+4.14E-1		3.33E1
+4.07E-1		3.20E1
+4.37E-1		3.24E1
+4.30E-1		3.19E1
+4.26E-1		3.09E1
+4.26E-1		3.04E1
+4.30E-1		3.00E1
+4.52E-1		2.88E1
+4.37E-1		2.86E1
+4.37E-1		2.80E1
+4.45E-1		2.51E1
+4.56E-1		2.47E1
+4.45E-1		2.40E1
+4.60E-1		2.40E1
+4.79E-1		2.33E1
+4.56E-1		2.26E1
+4.56E-1		2.21E1
+4.68E-1		2.22E1
+4.79E-1		2.15E1
+4.71E-1		2.08E1
+4.90E-1		2.10E1
+4.83E-1		1.89E1
+4.94E-1		1.97E1
+5.06E-1		1.97E1
+5.17E-1		1.95E1
+5.25E-1		1.94E1
+5.40E-1		1.95E1
+5.59E-1		1.96E1
+5.48E-1		1.91E1
+5.44E-1		1.87E1
+5.55E-1		1.87E1
+4.98E-1		1.81E1
+5.13E-1		1.78E1
+5.32E-1		1.77E1
+5.70E-1		1.82E1
+6.01E-1		1.84E1
+5.93E-1		1.82E1
+5.70E-1		1.75E1
+6.31E-1		1.82E1
+6.35E-1		1.77E1
+6.20E-1		1.74E1
+6.08E-1		1.72E1
+6.62E-1		1.80E1
+6.54E-1		1.77E1
+6.50E-1		1.74E1
+6.88E-1		1.80E1
+6.81E-1		1.76E1
+6.81E-1		1.72E1
+6.73E-1		1.69E1
+7.11E-1		1.70E1
+7.38E-1		1.66E1
+7.57E-1		1.63E1
+7.79E-1		1.61E1
+7.98E-1		1.55E1
+8.21E-1		1.53E1
+8.40E-1		1.60E1
+8.67E-1		1.55E1
+7.87E-1		1.52E1
+8.52E-1		1.51E1
+8.44E-1		1.48E1
+8.63E-1		1.47E1
+8.52E-1		1.43E1
+8.82E-1		1.43E1
+8.67E-1		1.36E1
+8.59E-1		1.34E1
+8.94E-1		1.34E1
+8.82E-1		1.33E1
+8.97E-1		1.30E1
+9.09E-1		1.28E1
+9.16E-1		1.24E1
+9.09E-1		1.21E1
+9.28E-1		1.19E1
+9.51E-1		1.13E1
+9.51E-1		1.09E1
+9.70E-1		1.03E1
+9.85E-1		1.02E1
+9.66E-1		9.61E0
+1.00E0		9.52E0
+9.96E-1		8.96E0
+1.01E0		8.49E0
+1.02E0		8.39E0
+1.03E0		7.83E0
+1.04E0		7.46E0
+1.05E0		8.30E0
+1.06E0		8.39E0
+1.07E0		8.48E0
+1.05E0		6.61E0
+1.06E0		6.89E0
+1.08E0		7.45E0
+1.09E0		6.99E0
+1.10E0		6.70E0
+1.11E0		6.33E0
+1.13E0		6.05E0
+1.16E0		5.58E0
+1.15E0		5.58E0
+1.14E0		5.11E0
+1.16E0		4.45E0
+1.17E0		4.83E0
+1.18E0		5.29E0
+1.19E0		5.58E0
+1.19E0		4.73E0
+1.20E0		4.73E0
+1.22E0		4.36E0
+1.23E0		4.08E0
+1.25E0		3.70E0
+1.28E0		3.79E0
+1.29E0		3.60E0
+1.33E0		3.32E0
+1.32E0		3.60E0
+1.32E0		2.95E0
+1.36E0		3.13E0
+1.40E0		2.38E0
+1.38E0		2.01E0
+1.32E0		1.73E0
+1.41E0		1.07E0
+1.43E0		1.35E0
+1.42E0		3.21E-1
+1.43E0		-2.41E-1
+1.45E0		-3.36E-1
+1.34E0		-7.03E-1
+1.35E0		-1.83E0
+1.37E0		-2.11E0
+1.38E0		-3.23E0
+1.39E0		-3.79E0
+1.40E0		-3.61E0
+1.44E0		-1.18E0
+1.46E0		-1.46E0
+1.46E0		-2.21E0
+1.47E0		-3.33E0
+1.48E0		-3.52E0
+1.48E0		-3.89E0
+1.42E0		-5.11E0
+1.43E0		-5.95E0
+1.44E0		-6.88E0
+1.48E0		-4.83E0
+1.49E0		-5.02E0
+1.49E0		-6.70E0
+1.51E0		-6.80E0
+1.49E0		-7.54E0
+1.46E0		-8.48E0
+1.48E0		-8.95E0
+1.52E0		-8.48E0
+1.49E0		-9.41E0
+1.50E0		-1.02E1
+1.52E0		-1.02E1
+1.53E0		-9.88E0
+1.54E0		-9.60E0
+1.56E0		-9.42E0
+1.58E0		-8.95E0
+1.59E0		-9.42E0
+1.58E0		-9.89E0
+1.57E0		-1.02E1
+1.56E0		-1.04E1
+1.55E0		-1.06E1
+1.51E0		-1.12E1
+1.53E0		-1.13E1
+1.55E0		-1.12E1
+1.57E0		-1.16E1
+1.59E0		-1.10E1
+1.60E0		-1.07E1
+1.61E0		-1.04E1
+1.62E0		-1.01E1
+1.63E0		-9.70E0
+1.65E0		-9.52E0
+1.67E0		-9.14E0
+1.68E0		-8.86E0
+1.59E0		-7.74E0
+1.61E0		-8.39E0
+1.62E0		-8.02E0
+1.63E0		-8.02E0
+1.66E0		-8.02E0
+1.67E0		-7.74E0
+1.69E0		-7.56E0
+1.70E0		-7.37E0
+1.71E0		-7.09E0
+1.72E0		-6.81E0
+1.61E0		-6.90E0
+1.62E0		-6.71E0
+1.62E0		-6.52E0
+1.64E0		-6.52E0
+1.66E0		-6.80E0
+1.68E0		-6.62E0
+1.70E0		-5.97E0
+1.71E0		-5.69E0
+1.72E0		-5.59E0
+1.73E0		-5.31E0
+1.74E0		-5.31E0
+1.71E0		-8.68E0
+1.72E0		-8.40E0
+1.73E0		-8.12E0
+1.75E0		-7.75E0
+1.75E0		-7.56E0
+1.77E0		-7.47E0
+1.79E0		-7.47E0
+1.80E0		-7.19E0
+1.75E0		-6.62E0
+1.76E0		-6.34E0
+1.77E0		-6.16E0
+1.81E0		-6.25E0
+1.82E0		-5.97E0
+1.83E0		-5.79E0
+1.84E0		-5.32E0
+1.85E0		-5.04E0
+1.81E0		-5.04E0
+1.81E0		-5.41E0
+1.79E0		-5.78E0
+1.83E0		-4.38E0
+1.85E0		-3.82E0
+1.87E0		-4.67E0
+1.87E0		-4.29E0
+1.90E0		-4.01E0
+1.90E0		-3.73E0
+1.92E0		-3.36E0
+1.94E0		-3.17E0
+1.96E0		-2.89E0
+1.98E0		-2.52E0
+1.98E0		-2.24E0
+1.96E0		-2.52E0
+1.87E0		-3.17E0
+1.88E0		-2.89E0
+1.90E0		-2.61E0
+1.93E0		-2.14E0
+1.95E0		-1.86E0
+1.97E0		-1.49E0
+1.98E0		-1.30E0
+2.00E0		-9.32E-1
+2.02E0		-6.52E-1
+2.04E0		-3.73E-1
+2.05E0		-2.79E-1
+2.09E0		-9.51E-2
+2.11E0		-9.63E-2
+2.12E0		1.84E-1
+2.17E0		2.74E-1
+2.20E0		8.54E-2
+1.88E0		-5.60E0
+1.90E0		-5.04E0
+1.91E0		-4.86E0
+1.93E0		-4.39E0
+1.95E0		-4.11E0
+1.97E0		-3.74E0
+1.98E0		-3.36E0
+2.00E0		-2.99E0
+2.03E0		-2.62E0
+2.04E0		-2.24E0
+2.06E0		-1.96E0
+2.09E0		-1.03E0
+2.10E0		-1.41E0
+2.08E0		-1.50E0
+2.05E0		-1.59E0
+2.04E0		-1.40E0
+2.03E0		-1.49E0
+2.02E0		-1.77E0
+2.02E0		-2.06E0
+2.11E0		-1.13E0
+2.13E0		-1.03E0
+2.13E0		-9.40E-1
+2.16E0		-7.55E-1
+2.17E0		-7.55E-1
+2.17E0		-7.55E-1
+2.20E0		-4.76E-1
+2.21E0		-3.83E-1
+2.23E0		-1.97E-1
+2.24E0		-1.05E-1
+2.26E0		-1.06E-1
+2.26E0		-1.06E-1
+2.25E0		-3.86E-1
+2.26E0		-2.93E-1
+2.26E0		-2.93E-1
+2.29E0		-2.95E-1
+2.29E0		-4.82E-1
+2.30E0		-3.89E-1
+2.31E0		-4.83E-1
+2.31E0		-5.77E-1
+2.33E0		-5.78E-1
+2.33E0		-9.52E-1
+2.33E0		-7.65E-1
+2.33E0		-7.65E-1
+2.37E0		-7.67E-1
+2.37E0		-9.55E-1
+2.37E0		-9.55E-1
+2.37E0		-9.55E-1
+2.40E0		-1.14E0
+2.40E0		-1.14E0
+2.43E0		-1.24E0
+2.42E0		-1.43E0
+2.42E0		-1.43E0
+2.45E0		-1.71E0
+2.47E0		-1.80E0
+2.49E0		-1.90E0
+2.49E0		-1.90E0
+2.51E0		-2.18E0
+2.51E0		-2.18E0
+2.54E0		-2.37E0
+2.55E0		-2.46E0
+2.55E0		-2.46E0
+2.56E0		-2.46E0
+2.59E0		-2.93E0
+2.59E0		-2.93E0
+2.60E0		-2.84E0
+2.62E0		-2.65E0
+2.57E0		-2.93E0
+2.55E0		-1.62E0
+2.57E0		-1.81E0
+2.60E0		-1.81E0
+2.61E0		-1.81E0
+2.62E0		-1.72E0
+2.64E0		-1.63E0
+2.67E0		-1.35E0
+2.67E0		-1.35E0
+2.68E0		-1.44E0
+2.71E0		-1.35E0
+2.75E0		-1.45E0
+2.75E0		-1.45E0
+2.75E0		-1.35E0
+2.65E0		-3.12E0
+2.68E0		-3.31E0
+2.68E0		-3.03E0
+2.69E0		-3.03E0
+2.74E0		-3.32E0
+2.73E0		-3.04E0
+2.70E0		-2.94E0
+2.70E0		-2.66E0
+2.72E0		-1.54E0
+2.68E0		-2.56E0
+2.67E0		-2.56E0
+2.66E0		-2.47E0
+2.65E0		-2.09E0
+2.76E0		-2.76E0
+2.76E0		-2.76E0
+2.79E0		-2.67E0
+2.81E0		-2.48E0
+2.82E0		-2.39E0
+2.84E0		-2.20E0
+2.85E0		-2.20E0
+2.87E0		-2.01E0
+2.88E0		-1.83E0
+2.89E0		-1.55E0
+2.89E0		-1.55E0
+2.90E0		-1.36E0
+2.79E0		-1.45E0
+2.82E0		-1.45E0
+2.82E0		-1.45E0
+2.83E0		-1.45E0
+2.86E0		-1.28E0
+2.86E0		-1.08E0
+2.87E0		-9.86E-1
+2.88E0		-7.06E-1
+2.91E0		-7.08E-1
+2.99E0		-1.51E-1
+2.98E0		-1.51E-1
+2.97E0		-1.50E-1
+2.97E0		-1.50E-1
+2.95E0		-2.42E-1
+2.94E0		-5.23E-1
+2.93E0		-7.09E-1
+2.92E0		-7.08E-1
+2.99E0		-9.00E-1
+2.98E0		-8.06E-1
+2.96E0		-8.98E-1
+2.96E0		-1.27E0
+2.95E0		-1.37E0
+2.95E0		-1.46E0
+2.94E0		-1.65E0
Index: /trunk/anuga_core/validation_tests/Tests/Experimental_data/dam_break_Yeh_Petroff/numerical_Yeh_Petroff.py
===================================================================
--- /trunk/anuga_core/validation_tests/Tests/Experimental_data/dam_break_Yeh_Petroff/numerical_Yeh_Petroff.py	(revision 8673)
+++ /trunk/anuga_core/validation_tests/Tests/Experimental_data/dam_break_Yeh_Petroff/numerical_Yeh_Petroff.py	(revision 8673)
@@ -0,0 +1,166 @@
+"""Simple water flow example using ANUGA
+
+Water driven up a linear slope and time varying boundary,
+similar to a beach environment
+"""
+
+#------------------------------------------------------------------------------
+# Import necessary modules
+#------------------------------------------------------------------------------
+import sys
+import anuga
+from anuga import Domain as Domain
+from math import cos
+from numpy import zeros, float
+from time import localtime, strftime, gmtime
+#from balanced_dev import *
+from anuga.geometry.polygon import inside_polygon, is_inside_triangle
+
+
+#-------------------------------------------------------------------------------
+# Copy scripts to time stamped output directory and capture screen
+# output to file
+#-------------------------------------------------------------------------------
+time = strftime('%Y%m%d_%H%M%S',localtime())
+
+#output_dir = 'dam_break_'+time
+output_dir = '.'
+output_file = 'dam_break'
+
+#anuga.copy_code_files(output_dir,__file__)
+#start_screen_catcher(output_dir+'_')
+
+
+#------------------------------------------------------------------------------
+# Setup domain
+#------------------------------------------------------------------------------
+dx = 0.01
+dy = dx
+L = 1.6
+W = 0.61
+
+# structured mesh
+points, vertices, boundary = anuga.rectangular_cross(int(L/dx), int(W/dy), L, W, (0.0, 0.0))
+
+#domain = anuga.Domain(points, vertices, boundary) 
+domain = Domain(points, vertices, boundary) 
+
+domain.set_name(output_file)                
+domain.set_datadir(output_dir) 
+
+#------------------------------------------------------------------------------
+# Setup Algorithm, either using command line arguments
+# or override manually yourself
+#------------------------------------------------------------------------------
+from anuga.utilities.argparsing import parse_standard_args
+alg, cfl = parse_standard_args()
+domain.set_flow_algorithm(alg)
+domain.set_CFL(cfl)
+domain.set_minimum_allowed_height(0.01)
+
+#------------------------------------------------------------------------------
+# Setup initial conditions
+#------------------------------------------------------------------------------
+def stage(x,y):
+    h = zeros(len(x), float)
+    for i in range(len(x)):
+        if x[i] <= 0.4:
+            h[i] = 0.3
+        elif x[i] <= 0.4+0.5:
+            h[i] = 0.01
+        elif x[i] <= 0.4+0.5+0.12:
+            if 0.25 <= y[i] <= 0.25+0.12:
+                h[i] = 0.75
+            else:
+                h[i] = 0.01
+        else:
+            h[i] = 0.01
+    return h
+
+def elevation(x,y):
+    z = zeros(len(x), float)
+    for i in range(len(x)):
+        if 0.4+0.5 <= x[i] <= 0.4+0.5+0.12:
+            if 0.25 <= y[i] <= 0.25+0.12:
+                z[i] = 0.75
+    return z
+domain.set_quantity('stage', stage)
+domain.set_quantity('elevation',elevation)
+domain.set_quantity('friction', 0.0)
+
+#-----------------------------------------------------------------------------
+# Setup boundary conditions
+#------------------------------------------------------------------------------
+from math import sin, pi, exp
+Br = anuga.Reflective_boundary(domain)      # Solid reflective wall
+Bt = anuga.Transmissive_boundary(domain)    # Continue all values on boundary 
+Bd = anuga.Dirichlet_boundary([1,0.,0.]) # Constant boundary values
+
+# Associate boundary tags with boundary objects
+domain.set_boundary({'left': Br, 'right': Br, 'top': Br, 'bottom': Br})
+
+
+#===============================================================================
+##from anuga.visualiser import RealtimeVisualiser
+##vis = RealtimeVisualiser(domain)
+##vis.render_quantity_height("stage", zScale =h0*500, dynamic=True)
+##vis.colour_height_quantity('stage', (0.0, 0.5, 1.0))
+##vis.start()
+#===============================================================================
+
+
+#---------------------------------------------
+# Find triangle that contains the point Point
+# and print to file
+#---------------------------------------------
+##PointL = []
+##print "0.12/dx=",int(round(0.12/dx))
+##STOP
+##PointL=[(0.868, 0.255),
+##          (0.868, 0.265),
+##          (0.868, 0.275),
+##          (0.868, 0.285),
+##          (0.868, 0.295),
+##          (0.868, 0.305),
+##          (0.868, 0.315),
+##          (0.868, 0.325),
+##          (0.868, 0.335),
+##          (0.868, 0.345),
+##          (0.868, 0.355),
+##          (0.868, 0.365)]
+###PointR = []
+##PointR=[(1.032, 0.255),
+##          (1.032, 0.265),
+##          (1.032, 0.275),
+##          (1.032, 0.285),
+##          (1.032, 0.295),
+##          (1.032, 0.305),
+##          (1.032, 0.315),
+##          (1.032, 0.325),
+##          (1.032, 0.335),
+##          (1.032, 0.345),
+##          (1.032, 0.355),
+##          (1.032, 0.365)]
+##for i in range(len(PointR)):
+##    Point = PointR[i]
+##    for n in range(len(domain.triangles)):
+##        tri = domain.get_vertex_coordinates(n)
+##        if is_inside_triangle(Point,tri):
+##            #print 'Point is within triangle with vertices '+'%s'%tri
+##            n_point = n
+##            break
+##    print 'n_point = ',n_point
+
+#------------------------------------------------------------------------------
+# Evolve system through time
+#------------------------------------------------------------------------------
+for t in domain.evolve(yieldstep = 0.05, finaltime = 3.0):
+    #print domain.timestepping_statistics(track_speeds=True)
+    print domain.timestepping_statistics()
+    #vis.update()
+
+
+#test against know data
+    
+#vis.evolveFinished()
+
Index: /trunk/anuga_core/validation_tests/Tests/Experimental_data/dam_break_Yeh_Petroff/plot_results_force.py
===================================================================
--- /trunk/anuga_core/validation_tests/Tests/Experimental_data/dam_break_Yeh_Petroff/plot_results_force.py	(revision 8673)
+++ /trunk/anuga_core/validation_tests/Tests/Experimental_data/dam_break_Yeh_Petroff/plot_results_force.py	(revision 8673)
@@ -0,0 +1,131 @@
+"""
+    Quick plot of the Yeh-Petroff dam break outputs
+"""
+#---------------
+# Import Modules
+#---------------
+import numpy
+from matplotlib import pyplot as pyplot
+from anuga.utilities import plot_utils as util
+from numpy import sqrt, sin, cos, pi
+import csv, pprint
+from numpy import zeros
+
+
+#################Opening the reference solution################
+infile = open('for_exp.txt', 'r')
+table = []
+for row in csv.reader(infile, delimiter="\t"):
+    table.append(row)
+infile.close()
+
+tF = zeros((len(table),2)) #t for time, F for Force
+for r in range(len(table)):
+    if len(table[r])==2:
+        for c in range(len(table[r])):        
+            tF[r][c] = float(table[r][c])
+    elif len(table[r])==3:
+        for c in range(len(table[r])):        
+            if c==0:
+                tF[r][0] = float(table[r][c])
+            elif c==2:
+                tF[r][1] = float(table[r][c])
+            #else:
+            #    print "table[r][c]=",table[r][c]
+    else:
+        STOP
+        
+t = zeros(len(tF))        
+F = zeros(len(tF))
+for r in range(len(tF)):
+    for c in range(len(tF[0])):
+        if c==0:
+            t[r] = tF[r][c]
+        else:
+            F[r] = tF[r][c]
+
+
+indicesL = [21330,
+21334,
+21338,
+21342,
+21346,
+21350,
+21354,
+21358,
+21362,
+21366,
+21370,
+21374]
+##
+##indicesR = [24988,
+##24992,
+##24996,
+##25000,
+##25004,
+##25008,
+##25012,
+##25016,
+##25020,
+##25024,
+##25028,
+##25032]
+
+
+##indicesL = [
+##21086,
+##21090,
+##21094,
+##21098,
+##21102,
+##21106,
+##21110,
+##21114,
+##21118,
+##21122,
+##21126,
+##21130
+##    ]
+
+indicesR = [
+25232,
+25236,
+25240,
+25244,
+25248,
+25252,
+25256,
+25260,
+25264,
+25268,
+25272,
+25276
+    ]
+
+
+# Get the sww file
+p = util.get_output('dam_break.sww')
+p2=util.get_centroids(p)
+
+
+N = len(p2.time)
+Sh = zeros(N)
+for i in range(N):
+    for k in range(len(indicesL)):
+        Sh[i] += p2.stage[i,indicesL[k]]
+        Sh[i] -= p2.stage[i,indicesR[k]]
+force = 9.8*Sh
+    
+
+#Plot stages
+pyplot.clf()
+pyplot.plot(t,F,'bo', label='experimental')
+pyplot.plot(p2.time+0.1,force,'k-', label='numerical')
+pyplot.legend(loc='best')
+pyplot.title('Force acting to the column over time')
+pyplot.xlabel('Time')
+pyplot.ylabel('Force')
+pyplot.savefig('Force.png')
+pyplot.show()
+
+
