source: anuga_validation/circular_island/simple_circular/circular.py @ 5317

Last change on this file since 5317 was 5317, checked in by kristy, 15 years ago

update of circular.py.

File size: 4.2 KB
Line 
1"""Simple water flow example using ANUGA
2
3Water flowing down a channel with more complex topography
4"""
5
6#------------------------------------------------------------------------------
7# Import necessary modules
8#------------------------------------------------------------------------------
9from anuga.abstract_2d_finite_volumes.mesh_factory import rectangular_cross
10from anuga.shallow_water import Domain
11from anuga.shallow_water import Reflective_boundary
12from anuga.shallow_water import Dirichlet_boundary
13from anuga.shallow_water import Time_boundary
14from anuga.pmesh.mesh_interface import create_mesh_from_regions
15from anuga.shallow_water import Transmissive_Momentum_Set_Stage_boundary
16from math import tan, sqrt, sin, pi
17
18#------------------------------------------------------------------------------
19# Setup computational domain
20#------------------------------------------------------------------------------
21length = 30.
22width = 25.
23Cx = 12.96                  # centre of island on the x axis
24Cy = 13.8                   # centre of island on the y axis
25dx = dy = 1                 # Resolution: Length of subdivisions on both axes
26water_depth = 0.32          # Can be 0.32 or 0.42
27
28#boundary
29poly_domain = [[0,0],[length,0],[length,width],[0,width]]
30meshname = 'test.msh'
31
32# exporting asc grid
33xmin = 0
34xmax = length   
35ymin = 0
36ymax = width
37
38#Create interior region
39Dome = [[(Cx)-4,(Cy)-4],[(Cx)+4,(Cy)-4,],
40        [(Cx)+4,(Cy)+4],[(Cx)-4,(Cy)+4]]
41
42remainder_res=1
43Dome_res = .01
44
45interior_dome = [[Dome, Dome_res]]
46
47#Create mesh
48create_mesh_from_regions(poly_domain,
49                         boundary_tags={'wavemaker': [0], 'right': [1],
50                                        'top': [2], 'left': [3]},
51                         maximum_triangle_area = remainder_res,
52                         filename=meshname,
53                         interior_regions = interior_dome,
54                         use_cache=False,
55                         verbose=True)
56
57domain = Domain(meshname, use_cache=False, verbose = True)
58domain.set_name('circular')                  # Output name
59print domain.statistics()
60
61#------------------------------------------------------------------------------
62# Setup initial conditions
63#------------------------------------------------------------------------------
64def topography(x,y):
65    """Complex topography defined by a function of vectors x and y
66    """
67   
68                                   
69    z= 0*x                      # defining z for all values other than the if statements
70    r= 3.6                      # radius, provided in document
71    angle = 14                  # angle, provided in document
72    h= r*tan(angle/57.2957795)  # finding height of cone if not truncated
73   
74
75    N = len(x)
76    for i in range(N):
77       
78        #truncated top
79        if (x[i]-Cx)**2 + (y[i]-Cy)**2 <1.1**2:
80            z[i] += 0.625
81
82        # cone
83        if (x[i]-Cx)**2 + (y[i]-Cy)**2 <r**2 and (x[i]-Cx)**2 + (y[i]-Cy)**2 >1.1**2:
84           z[i] = -(sqrt(((x[i]-Cx)**2+(y[i]-Cy)**2)/((r/h)**2))-h)
85    return z   
86
87domain.set_quantity('elevation', topography, verbose=True)  # Use function for elevation
88domain.set_quantity('friction', 0.01)         # Constant friction
89domain.set_quantity('stage',water_depth)   # Dry initial condition
90
91
92#------------------------------------------------------------------------------
93# Setup boundary conditions
94#------------------------------------------------------------------------------
95# Create boundary function from timeseries provided in file
96
97##boundary_filename="ts2cnew1_input_20_80sec_new"
98##prepare_timeboundary(boundary_filename+'.txt')
99##
100##function = file_function(boundary_filename+'.tms',
101##                         domain, verbose=True)
102def wave_form(t):
103    return 0.1*sin(2*pi*t/50.)
104
105# Create and assign boundary objects
106Bw = Dirichlet_boundary([water_depth, 0, 0])          #wall
107Bt = Transmissive_Momentum_Set_Stage_boundary(domain, wave_form) #wavemaker
108domain.set_boundary({'left': Bw, 'right': Bw, 'top': Bw, 'wavemaker': Bt})
109
110
111
112
113#------------------------------------------------------------------------------
114# Evolve system through time
115#------------------------------------------------------------------------------
116for t in domain.evolve(yieldstep = 0.2, finaltime = 100.0):
117    print domain.timestepping_statistics()
118
119
120   
Note: See TracBrowser for help on using the repository browser.