"""Simple water flow example using ANUGA Benchmark problem from teh third International Workshop on Long-Wave Runup Models http://www.cee.cornell.edu/longwave/index.cfm?page=benchmark&problem=3 """ #------------------------------------------------------------------------------ # Import necessary modules #------------------------------------------------------------------------------ from anuga.abstract_2d_finite_volumes.mesh_factory import rectangular from anuga.shallow_water import Domain from anuga.shallow_water import Reflective_boundary from anuga.shallow_water import Dirichlet_boundary from anuga.shallow_water import Time_boundary from anuga.shallow_water import Transmissive_boundary #------------------------------------------------------------------------------ # Setup computational domain #------------------------------------------------------------------------------ from anuga.pmesh.mesh_interface import create_mesh_from_regions bounding_polygon = [[0,10],[10,10],[10,0],[0,0]] remainder_res = 0.25 meshname = 'test.msh' create_mesh_from_regions(bounding_polygon, boundary_tags={'top': [0], 'right': [1], 'bottom': [2], 'left': [3]}, maximum_triangle_area=remainder_res, filename=meshname, #interior_regions=interior_regions, use_cache=False, verbose=True) from anuga.abstract_2d_finite_volumes.pmesh2domain import pmesh_to_domain_instance from anuga.caching import cache domain = cache(pmesh_to_domain_instance, (meshname, Domain), dependencies = [meshname]) # Create domain domain.set_name('slide') # Output to file runup.sww domain.set_datadir('.') # Use current directory for output domain.set_quantities_to_be_stored(['stage',# Store all conserved quantities 'xmomentum', 'ymomentum']) #------------------------------------------------------------------------------ # Setup initial conditions #------------------------------------------------------------------------------ from math import tan, radians beta = radians(5.7) delta = 1. mu = 0.01 g = 9.81 def topography(x,y): return x*tan(beta) # H(x) = x tan(beta) domain.set_quantity('elevation', topography) # Use function for elevation domain.set_quantity('friction', 0.1) # Constant friction domain.set_quantity('stage', 0.0) # Constant negative initial stage from Numeric import sqrt, exp t = 0 def ho(x,y): return delta*exp(x) #return delta*( -(2*sqrt(x*mu*mu/(delta*tan(beta))) - sqrt(g/delta)*mu*t)**2) #------------------------------------------------------------------------------ # Setup boundary conditions #------------------------------------------------------------------------------ from math import sin, pi, exp Br = Reflective_boundary(domain) # Solid reflective wall Bt = Transmissive_boundary(domain) # Continue all values on boundary Bd = Dirichlet_boundary([-0.2,0.,0.]) # Constant boundary values Bw = Time_boundary(domain=domain, # Time dependent boundary f=lambda t: [(.1*sin(t*2*pi)-0.3) * exp(-t), 0.0, 0.0]) # Associate boundary tags with boundary objects domain.set_boundary({'left': Br, 'right': Bw, 'top': Br, 'bottom': Br}) #------------------------------------------------------------------------------ # Evolve system through time #------------------------------------------------------------------------------ for t in domain.evolve(yieldstep = 0.1, finaltime = 10.0): t1 = domain.get_time() t = t1 domain.set_quantity('elevation', ho) domain.write_time()