1 | """Example of shallow water wave equation |
---|
2 | consisting of an asymetrical converging channel. |
---|
3 | |
---|
4 | Copyright 2005 |
---|
5 | Christopher Zoppou, Stephen Roberts |
---|
6 | ANU |
---|
7 | |
---|
8 | Specific methods pertaining to the 2D shallow water equation |
---|
9 | are imported from shallow_water |
---|
10 | for use with the generic finite volume framework |
---|
11 | |
---|
12 | Conserved quantities are h, uh and vh stored as elements 0, 1 and 2 in the |
---|
13 | numerical vector named conserved_quantities. |
---|
14 | """ |
---|
15 | |
---|
16 | #------------------------------- |
---|
17 | # Set up path and import modules |
---|
18 | # import visualise2_chris as visualise |
---|
19 | # import Image, ImageGrab |
---|
20 | import sys |
---|
21 | from os import sep |
---|
22 | sys.path.append('..'+sep+'pyvolution') |
---|
23 | |
---|
24 | from shallow_water import Domain, Constant_height |
---|
25 | from shallow_water import Transmissive_boundary, Reflective_boundary,\ |
---|
26 | Dirichlet_boundary |
---|
27 | from math import sqrt, cos, sin, pi |
---|
28 | from mesh_factory import oblique |
---|
29 | |
---|
30 | #-------------- |
---|
31 | # Define domain |
---|
32 | n = 60 |
---|
33 | m = 80 |
---|
34 | leny = 30. |
---|
35 | lenx = 40. |
---|
36 | n = 50 |
---|
37 | m = 60 |
---|
38 | points, elements, boundary = oblique(m, n, lenx, leny) |
---|
39 | domain = Domain(points, elements, boundary) |
---|
40 | |
---|
41 | #---------------- |
---|
42 | # Order of scheme |
---|
43 | domain.default_order=1 |
---|
44 | |
---|
45 | #--------------------------------- |
---|
46 | # Store output format and location |
---|
47 | domain.store = True |
---|
48 | domain.format = "sww" #"sww" for NET.CDF binary format or "dat" for ASCII |
---|
49 | domain.filename = "oblique_first_order" |
---|
50 | |
---|
51 | #------------------------ |
---|
52 | # Visualization smoothing |
---|
53 | domain.smooth=True |
---|
54 | |
---|
55 | #-------------- |
---|
56 | # Set bed slope |
---|
57 | def x_slope(x, y): |
---|
58 | return 0*x |
---|
59 | domain.set_quantity('elevation', x_slope) |
---|
60 | |
---|
61 | #------------- |
---|
62 | # Set friction |
---|
63 | domain.set_quantity('friction', 0.0) |
---|
64 | |
---|
65 | #-------------------- |
---|
66 | # Boundary conditions |
---|
67 | R = Reflective_boundary(domain) |
---|
68 | T = Transmissive_boundary(domain) |
---|
69 | D = Dirichlet_boundary([1.0, 8.57, 0.0]) |
---|
70 | domain.set_boundary({'left': D, 'right': T, 'top': R, 'bottom': R}) |
---|
71 | |
---|
72 | #------------------ |
---|
73 | # Initial condition |
---|
74 | h = 0.5 |
---|
75 | domain.set_quantity('level', Constant_height(x_slope, h) ) |
---|
76 | |
---|
77 | #---------------------------------------------------------- |
---|
78 | # Decide which quantities are to be stored at each timestep |
---|
79 | domain.quantities_to_be_stored = ['level', 'xmomentum', 'ymomentum'] |
---|
80 | |
---|
81 | |
---|
82 | #---------- |
---|
83 | # Evolution |
---|
84 | import time |
---|
85 | t0 = time.time() |
---|
86 | for t in domain.evolve(yieldstep = 1.0, finaltime = 50): |
---|
87 | domain.write_time() |
---|
88 | |
---|
89 | print 'That took %.2f seconds' %(time.time()-t0) |
---|
90 | |
---|
91 | #----------------------------------- |
---|
92 | #Save the last frame as an EPS file |
---|
93 | #filename = 'ccube.eps' |
---|
94 | #print 'Saving last frame in EPS format in ', filemame |
---|
95 | #im = ImageGrab.grab() |
---|
96 | #im.save(filename) |
---|
97 | |
---|
98 | |
---|
99 | |
---|