# source:branches/anuga_1_1/anuga_work/production/west_tas_2008/smf_approx.py@7799

Last change on this file since 7799 was 4968, checked in by sexton, 16 years ago

SMF modelling for West Tas

File size: 3.9 KB
Line
1"""Script for examing potential tsunami amplitudes for a submarine mass failure for West Tasmania, Australia.
2
3Ole Nielsen and Duncan Gray, GA - 2005 and Adrian Hitchman and Jane Sexton, GA - 2006
4"""
5
6#-------------------------------------------------------------------------------# Import necessary modules
7#-------------------------------------------------------------------------------
8
9# Standard modules
10import os
11import time
12
13from pylab import *
14ion()
15hold(True)
16
17def stats_slide(depth,length,width,thickness,slope,gamma):
18    from math import sin, tan, radians, pi, sqrt, exp
19
20    gravity = 9.8
21    massco = 1.0
22    dragco = 1.0
23    psi = 0.0
27
28    a0 = gravity * sint * ((gamma-1)/(gamma+massco)) * (1-(tanp/tant))
29    ut = sqrt((gravity*depth) * (length*sint/depth) \
30                    * (pi*(gamma-1)/(2*dragco)) * (1-(tanp/tant)))
31    s0 = ut**2 / a0
32    t0 = ut / a0
33
34    #calculate some parameters of the water displacement produced by the slide
35
36    w = t0 * sqrt(gravity*depth)
37    a2D = s0 * (0.0574 - (0.0431*sint)) \
38             * thickness/length \
39             * ((length*sint/depth)**1.25) \
40             * (1 - exp(-2.2*(gamma-1)))
41    a3D = a2D / (1 + (15.5*sqrt(depth/(length*sint))))
42
43    #scaled_amp = a3D/a2D
44
45    return a3D, w
46
47def stats_slump(depth,length,width,thickness,slope,gamma):
48    from math import sin, radians, sqrt
49
50    radius = length**2 / (8.0 * thickness)
51    massco = 1.0
52    dphi = 0.48
53    gravity = 9.8
55
56    s0 = radius * dphi / 2
57    t0 = sqrt((radius*(gamma+massco)) / (gravity*(gamma-1)))
58    a0 = s0 / t0**2
59    um = s0 / t0
60
61    #calculate some parameters of the water displacement produced by the slump
62
63    w = t0 * sqrt(gravity*depth)
64    a2D = s0 * (0.131/sint) \
65             * thickness/length \
66             * (length*sint/depth)**1.25 \
67             * (length/radius)**0.63 * dphi**0.39 \
68             * (1.47 - (0.35*(gamma-1))) * (gamma-1)
69    a3D = a2D / (1. + (2.06*sqrt(depth/length)))
70
71    #scaled_amp = a3D/a2D
72
73    return a3D, w
74
75which_one = 'slide'
76#which_one = 'slump'
77par = 'depth'
78
79# as supplied by Cameron Mitchell end of 2007
80depths = range(750,2500,50)
81lengths = [25000]
82widths = [50000]
83thicknesses = [400]
84slopes = [4]
85gammas = [1.46]
86smftype = ['West Tas']
87smfdepths = [1800]
88
89if par == 'depth':
90
91    for i in range(len(lengths)):
92        amp = []
93        wavelength = []
94        for depth in depths:
95
96            if which_one == 'slide':
97                a, w = stats_slide(depth,lengths[i],widths[i],thicknesses[i],slopes[i],gammas[i])
98            else:
99                a, w = stats_slump(depth,lengths[i],widths[i],thicknesses[i],slopes[i],gammas[i])
100
101            amp.append(a)
102            wavelength.append(w/1000.)
103
104        figure(1)
105        plot(depths,amp)
106        xlabel('Depth (m)')
107        ylabel('amplitude (m)')
108
109        figure(2)
110        plot(depths,wavelength)
111        xlabel('Depth (m)')
112        ylabel('wavelength (km)')
113
114    figure(1)
115    legend(smftype)
116    savefig('depth_vs_amp')
117
118    figure(2)
119    legend(smftype)
120    savefig('depth_vs_wavelength')
121
123    wdata = []
124    for i in range(len(smfdepths)):
125        if which_one == 'slide':
126            a, w = stats_slide(smfdepths[i],lengths[i],widths[i],thicknesses[i],slopes[i],gammas[i])
127        else:
128            a, w = stats_slump(smfdepths[i],lengths[i],widths[i],thicknesses[i],slopes[i],gammas[i])