# Changeset 7818 for anuga_work/development/anuga_1d/test_quantity.py

Ignore:
Timestamp:
Jun 11, 2010, 10:27:11 AM (12 years ago)
Message:

working test_quantity

File:
1 edited

Unmodified
Added
Removed
• ## anuga_work/development/anuga_1d/test_quantity.py

 r7815 quantity2.set_values(2*quantity1) print quantity2.vertex_values assert allclose(quantity2.vertex_values, [[0,2], [2,4], [4,6], [6,8]]) #Assert that quantities are conserved from numpy import sum from Numeric import sum for k in range(quantity.centroid_values.shape[0]): assert allclose (quantity.centroid_values[k], quantity.update(timestep) sem = array([1.,1.,1.,1.])/array([1, 2, 3, 4]) sem = array([1.,1.,1.,1.]) denom = ones(4, Float)-timestep*sem assert allclose( quantity.centroid_values, x) #Test smoothing def test_smoothing(self): from shallow_water import Domain, Transmissive_boundary from Numeric import zeros, Float from anuga.utilities.numerical_tools import mean #Create shallow water domain domain = Domain(points10) domain.default_order=2 domain.reduction = mean #Set some field values domain.set_quantity('elevation', lambda x: x) domain.set_quantity('friction', 0.03) ###################### # Boundary conditions B = Transmissive_boundary(domain) domain.set_boundary( {'left': B, 'right': B, 'top': B, 'bottom': B}) ###################### #Initial condition - with jumps bed = domain.quantities['elevation'].vertex_values stage = zeros(bed.shape, Float) h = 0.03 for i in range(stage.shape[0]): if i % 2 == 0: stage[i,:] = bed[i,:] + h else: stage[i,:] = bed[i,:] domain.set_quantity('stage', stage) stage = domain.quantities['stage'] #Get smoothed stage A, V = stage.get_vertex_values(xy=False, smooth=True) Q = stage.vertex_values assert A.shape[0] == 9 assert V.shape[0] == 8 assert V.shape[1] == 3 #First four points assert allclose(A[0], (Q[0,2] + Q[1,1])/2) assert allclose(A[1], (Q[1,0] + Q[3,1] + Q[2,2])/3) assert allclose(A[2], Q[3,0]) assert allclose(A[3], (Q[0,0] + Q[5,1] + Q[4,2])/3) #Center point assert allclose(A[4], (Q[0,1] + Q[1,2] + Q[2,0] +\ Q[5,0] + Q[6,2] + Q[7,1])/6) #Check V assert allclose(V[0,:], [3,4,0]) assert allclose(V[1,:], [1,0,4]) assert allclose(V[2,:], [4,5,1]) assert allclose(V[3,:], [2,1,5]) assert allclose(V[4,:], [6,7,3]) assert allclose(V[5,:], [4,3,7]) assert allclose(V[6,:], [7,8,4]) assert allclose(V[7,:], [5,4,8]) #Get smoothed stage with XY X, Y, A1, V1 = stage.get_vertex_values(xy=True, smooth=True) assert allclose(A, A1) assert allclose(V, V1) #Check XY assert allclose(X[4], 0.5) assert allclose(Y[4], 0.5) assert allclose(X[7], 1.0) assert allclose(Y[7], 0.5) def test_vertex_values_no_smoothing(self): from mesh_factory import rectangular from shallow_water import Domain, Transmissive_boundary from Numeric import zeros, Float from anuga.utilities.numerical_tools import mean #Create basic mesh points, vertices, boundary = rectangular(2, 2) #Create shallow water domain domain = Domain(points, vertices, boundary) domain.default_order=2 domain.reduction = mean #Set some field values domain.set_quantity('elevation', lambda x,y: x) domain.set_quantity('friction', 0.03) ###################### #Initial condition - with jumps bed = domain.quantities['elevation'].vertex_values stage = zeros(bed.shape, Float) h = 0.03 for i in range(stage.shape[0]): if i % 2 == 0: stage[i,:] = bed[i,:] + h else: stage[i,:] = bed[i,:] domain.set_quantity('stage', stage) #Get stage stage = domain.quantities['stage'] A, V = stage.get_vertex_values(xy=False, smooth=False) Q = stage.vertex_values.flat for k in range(8): assert allclose(A[k], Q[k]) for k in range(8): assert V[k, 0] == 3*k assert V[k, 1] == 3*k+1 assert V[k, 2] == 3*k+2 X, Y, A1, V1 = stage.get_vertex_values(xy=True, smooth=False) assert allclose(A, A1) assert allclose(V, V1) #Check XY assert allclose(X[1], 0.5) assert allclose(Y[1], 0.5) assert allclose(X[4], 0.0) assert allclose(Y[4], 0.0) assert allclose(X[12], 1.0) assert allclose(Y[12], 0.0) #------------------------------------------------------------- if __name__ == "__main__": suite = unittest.makeSuite(Test_Quantity, 'test_set') suite = unittest.makeSuite(Test_Quantity, 'test') runner = unittest.TextTestRunner() runner.run(suite)
Note: See TracChangeset for help on using the changeset viewer.