-
Notifications
You must be signed in to change notification settings - Fork 33
/
test_okada.py
108 lines (93 loc) · 3.55 KB
/
test_okada.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
from okada_wrapper import dc3d0wrapper, dc3dwrapper
from numpy import linspace, zeros, log
from matplotlib.pyplot import contourf, contour, \
xlabel, ylabel, colorbar, show, savefig
import matplotlib
import time
matplotlib.rcParams['mathtext.fontset'] = 'cm'
matplotlib.rcParams['font.size'] = 14
matplotlib.rcParams['xtick.direction'] = 'out'
matplotlib.rcParams['ytick.direction'] = 'out'
matplotlib.rcParams['lines.linewidth'] = 1
def get_params():
source_depth = 3.0
obs_depth = 3.0
poisson_ratio = 0.25
mu = 1.0
dip = 90
lmda = 2 * mu * poisson_ratio / (1 - 2 * poisson_ratio)
alpha = (lmda + mu) / (lmda + 2 * mu)
return source_depth, obs_depth, poisson_ratio, mu, dip, alpha
def test_dc3d0():
source_depth, obs_depth, poisson_ratio, mu, dip, alpha = get_params()
n = (100, 100)
x = linspace(-1, 1, n[0])
y = linspace(-1, 1, n[1])
ux = zeros((n[0], n[1]))
for i in range(100):
for j in range(100):
success, u, grad_u = dc3d0wrapper(alpha,
[x[i], y[j], -obs_depth],
source_depth, dip,
[1.0, 0.0, 0.0, 0.0])
assert (success == 0)
ux[i, j] = u[0]
cntrf = contourf(x, y, log(abs(ux.T)))
contour(x, y, log(abs(ux.T)), colors='k', linestyles='solid')
xlabel(r'$\mathrm{x}$')
ylabel(r'$\mathrm{y}$')
cbar = colorbar(cntrf)
cbar.set_label(r'$\log(u_{\mathrm{x}})$')
show()
def test_dc3d():
source_depth, obs_depth, poisson_ratio, mu, dip, alpha = get_params()
n = (100, 100)
x = linspace(-1, 1, n[0])
y = linspace(-1, 1, n[1])
ux = zeros((n[0], n[1]))
for i in range(n[0]):
for j in range(n[1]):
success, u, grad_u = dc3dwrapper(alpha,
[x[i], y[j], -obs_depth],
source_depth, dip,
[-0.6, 0.6], [-0.6, 0.6],
[1.0, 0.0, 0.0])
assert (success == 0)
ux[i, j] = u[0]
levels = linspace(-0.5, 0.5, 21)
cntrf = contourf(x, y, ux.T, levels=levels)
contour(x, y, ux.T, colors='k', levels=levels, linestyles='solid')
xlabel(r'$\mathrm{x}$')
ylabel(r'$\mathrm{y}$')
cbar = colorbar(cntrf)
tick_locator = matplotlib.ticker.MaxNLocator(nbins=5)
cbar.locator = tick_locator
cbar.update_ticks()
cbar.set_label(r'$u_{\mathrm{x}}$')
savefig("strike_slip.png")
show()
def test_success():
# should fail because z is positive
success, u, grad_u = dc3d0wrapper(1.0, [0.0, 0.0, 1.0],
1.0, 90,
[1.0, 0.0, 0.0, 0.0])
assert (success == 2)
success, u, grad_u = dc3dwrapper(1.0, [0.0, 0.0, 1.0],
0.0, 90, [-0.7, 0.7], [-0.7, 0.7],
[1.0, 0.0, 0.0])
assert (success == 2)
def benchmark():
n = 1000000
start = time.time()
for i in range(n):
success, u, grad_u = dc3dwrapper(1.0, [0.0, 0.0, 1.0],
0.0, 90, [-0.7, 0.7], [-0.7, 0.7],
[1.0, 0.0, 0.0])
end = time.time()
T = end - start
print(str(n) + " DC3D evaluations took: " + str(T) + " seconds")
print("Giving a time of " + str(T / n) + " seconds per evaluation.")
if __name__ == '__main__':
# test_dc3d0()
test_dc3d()
# benchmark()