# source:anuga_core/source/anuga/shallow_water/eqf.py@5611

Last change on this file since 5611 was 4320, checked in by sexton, 17 years ago

new earthquake source function (converted fortran to python), plus new figures for slide report

File size: 4.0 KB
Line
1#
2# earthquake_tsunami function
3#
4
5"""This function returns a callable object representing an initial water
6   displacement generated by a submarine earthqauke.
7
8Using input parameters:
9
10Required
11 length  along-stike length of rupture area
12 width   down-dip width of rupture area
13 strike  azimuth (degrees, measured from north) of fault axis
14 dip     angle of fault dip in degrees w.r.t. horizontal
15 depth   depth to base of rupture area
16
17Optional
18 x0      x origin (0)
19 y0      y origin (0)
20 slip    metres of fault slip (1)
21 rake    angle of slip (w.r.t. horizontal) in fault plane (90 degrees)
22
23The returned object is a callable okada function that represents
24the initial water displacement generated by a submarine earthuake.
25
26"""
27
28def earthquake_tsunami(length, width, strike, depth, \
29                       dip, x0=0.0, y0=0.0, slip=1.0, rake=90.,\
30                       domain=None, verbose=False):
31
32    from math import sin, radians
33
34    if domain is not None:
35        xllcorner = domain.geo_reference.get_xllcorner()
36        yllcorner = domain.geo_reference.get_yllcorner()
37        x0 = x0 - xllcorner  # fault origin (relative)
38        y0 = y0 - yllcorner
39
40    #a few temporary print statements
41    if verbose is True:
42        print '\nThe Earthquake ...'
43        print '\tLength: ', length
44        print '\tDepth: ', depth
45        print '\tStrike: ', strike
46        print '\tWidth: ', width
47        print '\tDip: ', dip
48        print '\tSlip: ', slip
49        print '\tx0: ', x0
50        print '\ty0: ', y0
51
52    # warning test
53    test = width*1000.0*sin(radians(dip)) - depth
54
55    if verbose is True:
56        if test > 0.0:
57            print 'Earthquake source not located below seafloor'
59
60    return Okada_func(length=length, width=width, dip=dip, \
61                      x0=x0, y0=y0, strike=strike, depth=depth, \
62                      slip=slip, rake=rake, test=test)
63
64
65#
67#
68
69"""This is a callable class representing the initial water displacment
70   generated by an earthquake.
71
72Using input parameters:
73
74Required
75 length  along-stike length of rupture area
76 width   down-dip width of rupture area
77 strike  azimuth (degrees, measured from north) of fault axis
78 dip     angle of fault dip in degrees w.r.t. horizontal
79 depth   depth to base of rupture area
80
81Optional
82 x0      x origin (0)
83 y0      y origin (0)
84 slip    metres of fault slip (1)
85 rake    angle of slip (w.r.t. horizontal) in fault plane (90 degrees)
86
87"""
88
90
91    def __init__(self, length, width, dip, x0, y0, strike, \
92                 depth, slip, rake, test):
93        self.dip = dip
94        self.length = length
95        self.width = width
96        self.x0 = x0
97        self.y0 = y0
98        self.strike = strike
99        self.depth = depth
100        self.slip = slip
101        self.rake = rake
102        self.test = test
103
104
105    def __call__(self, x, y):
106        """Make Okada_func a callable object.
107
108        If called as a function, this object returns z values representing
109        the initial 3D distribution of water heights at the points (x,y)
110        produced by a submarine mass failure.
111        """
112
113        from math import sin, cos, radians, exp, cosh
114        from Numeric import zeros, Float
116
117        #ensure vectors x and y have the same length
118        N = len(x)
119        assert N == len(y)
120
121        depth = self.depth
122        dip = self.dip
123        length = self.length
124        width = self.width
125        x0 = self.x0
126        y0 = self.y0
127        strike = self.strike
128        dip = self.dip
129        rake = self.rake
130        slip = self.slip
131
132        #double Gaussian calculation assumes water displacement is oriented
133        #E-W, so, for displacement at some angle alpha clockwise from the E-W
134        #direction, rotate (x,y) coordinates anti-clockwise by alpha
135