1 | """Example of shallow water wave equation. |
---|
2 | |
---|
3 | This version is for profiling |
---|
4 | (The original) |
---|
5 | """ |
---|
6 | |
---|
7 | ###################### |
---|
8 | # Module imports |
---|
9 | # |
---|
10 | from shallow_water import Domain, Reflective_boundary, Dirichlet_boundary,\ |
---|
11 | Transmissive_boundary, Time_boundary,\ |
---|
12 | Weir_simple as Weir, Constant_height |
---|
13 | |
---|
14 | from mesh_factory import rectangular |
---|
15 | from Numeric import array |
---|
16 | |
---|
17 | |
---|
18 | ###################### |
---|
19 | # Domain |
---|
20 | # |
---|
21 | |
---|
22 | N = 16 |
---|
23 | |
---|
24 | print 'Creating domain' |
---|
25 | #Create basic mesh |
---|
26 | points, vertices, boundary = rectangular(N, N) |
---|
27 | |
---|
28 | #Create shallow water domain |
---|
29 | domain = Domain(points, vertices, boundary) |
---|
30 | domain.default_order = 2 |
---|
31 | #domain.visualise = True |
---|
32 | |
---|
33 | |
---|
34 | print 'Field values' |
---|
35 | def slope(x, y): |
---|
36 | return -x |
---|
37 | |
---|
38 | domain.set_quantity('elevation', slope) |
---|
39 | domain.set_quantity('friction', 0.3) |
---|
40 | |
---|
41 | |
---|
42 | # Boundary conditions |
---|
43 | # |
---|
44 | Br = Reflective_boundary(domain) |
---|
45 | domain.set_boundary({'left': Br, 'right': Br, 'bottom': Br, 'top': Br}) |
---|
46 | |
---|
47 | |
---|
48 | ###################### |
---|
49 | #Initial condition |
---|
50 | domain.set_quantity('level', Constant_height(slope, 0.3)) |
---|
51 | |
---|
52 | domain.check_integrity() |
---|
53 | |
---|
54 | |
---|
55 | ###################### |
---|
56 | #Evolution |
---|
57 | |
---|
58 | import time |
---|
59 | t0 = time.time() |
---|
60 | |
---|
61 | s = 'for t in domain.evolve(yieldstep = 0.02, finaltime = 0.2): domain.write_time()' |
---|
62 | |
---|
63 | |
---|
64 | |
---|
65 | import profile, pstats |
---|
66 | FN = 'profile.dat' |
---|
67 | |
---|
68 | profile.run(s, FN) |
---|
69 | |
---|
70 | print 'That took %.2f seconds' %(time.time()-t0) |
---|
71 | |
---|
72 | S = pstats.Stats(FN) |
---|
73 | #S.sort_stats('time').print_stats(20) |
---|
74 | s = S.sort_stats('cumulative').print_stats(30) |
---|
75 | |
---|
76 | print s |
---|