forked from python-hydro/pyro2
-
Notifications
You must be signed in to change notification settings - Fork 0
Expand file tree
/
Copy pathtest_mask.py
More file actions
109 lines (70 loc) · 2.38 KB
/
Copy pathtest_mask.py
File metadata and controls
109 lines (70 loc) · 2.38 KB
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
109
#!/usr/bin/env python3
import time
import numpy as np
import pyro.mesh.boundary as bnd
from pyro.mesh import patch
class Mask(object):
def __init__(self, nx, ny, ng):
self.nx = nx
self.ny = ny
self.ng = ng
ilo = ng
ihi = ng+nx-1
jlo = ng
jhi = ng+ny-1
# just the interior cells
self.valid = self._mask_array(nx, ny, ng)
# shifts in x
self.ip1 = self._mask_array(nx, ny, ng)
self.im1 = self._mask_array(nx, ny, ng)
self.ip2 = self._mask_array(nx, ny, ng)
self.im2 = self._mask_array(nx, ny, ng)
arrays = [self.valid, self.ip1, self.im1, self.ip2, self.im2]
shifts = [0, 1, -1, 2, -2]
for a, s in zip(arrays, shifts):
a[ilo+s:ihi+1+s, jlo:jhi+1] = True
# shifts in y
self.jp1 = self._mask_array(nx, ny, ng)
self.jm1 = self._mask_array(nx, ny, ng)
self.jp2 = self._mask_array(nx, ny, ng)
self.jm2 = self._mask_array(nx, ny, ng)
arrays = [self.jp1, self.jm1, self.jp2, self.jm2]
shifts = [1, -1, 2, -2]
for a, s in zip(arrays, shifts):
a[ilo:ihi+1, jlo+s:jhi+1+s] = True
def _mask_array(self, nx, ny, ng):
return np.zeros((nx+2*ng, ny+2*ng), dtype=bool)
n = 1024
myg = patch.Grid2d(n, 2*n, xmax=1.0, ymax=2.0)
myd = patch.CellCenterData2d(myg)
bc = bnd.BC()
myd.register_var("a", bc)
myd.create()
a = myd.get_var("a")
a[:, :] = np.random.rand(myg.qx, myg.qy)
# slicing method
start = time.time()
da = myg.scratch_array()
da[myg.ilo:myg.ihi+1, myg.jlo:myg.jhi+1] = \
a[myg.ilo+1:myg.ihi+2, myg.jlo:myg.jhi+1] - \
a[myg.ilo-1:myg.ihi, myg.jlo:myg.jhi+1]
print("slice method: ", time.time() - start)
# mask method
m = Mask(myg.nx, myg.ny, myg.ng)
start = time.time()
da2 = myg.scratch_array()
da2[m.valid] = a[m.ip1] - a[m.im1]
print("mask method: ", time.time() - start)
print(np.max(np.abs(da2 - da)))
# roll -- note, we roll in the opposite direction of the shift
start = time.time()
da3 = myg.scratch_array()
da3[:] = np.roll(a, -1, axis=0) - np.roll(a, 1, axis=0)
print("roll method: ", time.time() - start)
print(np.max(np.abs(da3[m.valid] - da[m.valid])))
# ArrayIndex
start = time.time()
da4 = myg.scratch_array()
da4.v()[:, :] = a.ip(1) - a.ip(-1)
print("ArrayIndex method: ", time.time() - start)
print(np.max(np.abs(da4[m.valid] - da[m.valid])))