Skip to content

Commit 334eeb9

Browse files
authored
Merge pull request #226 from IntelPython/feature/npbench_dpnp_implementation
Add dpnp implementation for npbench
2 parents d858ecd + bb04b74 commit 334eeb9

File tree

53 files changed

+1757
-0
lines changed

Some content is hidden

Large Commits have some content hidden by default. Use the searchbox below for content that may be hidden.

53 files changed

+1757
-0
lines changed
Lines changed: 19 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,19 @@
1+
# SPDX-FileCopyrightText: 2014 Jérôme Kieffer et al.
2+
# SPDX-FileCopyrightText: 2021 ETH Zurich and the NPBench authors
3+
# SPDX-FileCopyrightText: 2022 - 2023 Intel Corporation
4+
#
5+
# SPDX-License-Identifier: BSD-3-Clause
6+
7+
"""
8+
Jérôme Kieffer and Giannis Ashiotis. Pyfai: a python library for
9+
high performance azimuthal integration on gpu, 2014. In Proceedings of the
10+
7th European Conference on Python in Science (EuroSciPy 2014).
11+
"""
12+
13+
import dpnp as np
14+
15+
16+
def azimint_hist(data, radius, npt):
17+
histu = np.histogram(radius, npt)[0]
18+
histw = np.histogram(radius, npt, weights=data)[0]
19+
return histw / histu
Lines changed: 25 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,25 @@
1+
# SPDX-FileCopyrightText: 2014 Jérôme Kieffer et al.
2+
# SPDX-FileCopyrightText: 2021 ETH Zurich and the NPBench authors
3+
# SPDX-FileCopyrightText: 2022 - 2023 Intel Corporation
4+
#
5+
# SPDX-License-Identifier: BSD-3-Clause
6+
7+
"""
8+
Jérôme Kieffer and Giannis Ashiotis. Pyfai: a python library for
9+
high performance azimuthal integration on gpu, 2014. In Proceedings of the
10+
7th European Conference on Python in Science (EuroSciPy 2014).
11+
"""
12+
13+
import dpnp as np
14+
15+
16+
def azimint_naive(data, radius, npt):
17+
rmax = radius.max()
18+
res = np.zeros(npt, dtype=np.float64)
19+
for i in range(npt):
20+
r1 = rmax * i / npt
21+
r2 = rmax * (i + 1) / npt
22+
mask_r12 = np.logical_and((r1 <= radius), (radius < r2))
23+
values_r12 = data[mask_r12]
24+
res[i] = values_r12.mean()
25+
return res
Lines changed: 109 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,109 @@
1+
# SPDX-FileCopyrightText: 2017 Lorena A. Barba, Gilbert F. Forsyth.
2+
# SPDX-FileCopyrightText: 2018 Barba, Lorena A., and Forsyth, Gilbert F.
3+
# SPDX-FileCopyrightText: 2021 ETH Zurich and the NPBench authors
4+
# SPDX-FileCopyrightText: 2022 - 2023 Intel Corporation
5+
#
6+
# SPDX-License-Identifier: BSD-3-Clause
7+
8+
"""
9+
CFD Python: the 12 steps to Navier-Stokes equations.
10+
Journal of Open Source Education, 1(9), 21,
11+
https://doi.org/10.21105/jose.00021
12+
"""
13+
14+
import dpnp as np
15+
16+
17+
def build_up_b(b, rho, dt, u, v, dx, dy):
18+
b[1:-1, 1:-1] = rho * (
19+
1
20+
/ dt
21+
* (
22+
(u[1:-1, 2:] - u[1:-1, 0:-2]) / (2 * dx)
23+
+ (v[2:, 1:-1] - v[0:-2, 1:-1]) / (2 * dy)
24+
)
25+
- ((u[1:-1, 2:] - u[1:-1, 0:-2]) / (2 * dx)) ** 2
26+
- 2
27+
* (
28+
(u[2:, 1:-1] - u[0:-2, 1:-1])
29+
/ (2 * dy)
30+
* (v[1:-1, 2:] - v[1:-1, 0:-2])
31+
/ (2 * dx)
32+
)
33+
- ((v[2:, 1:-1] - v[0:-2, 1:-1]) / (2 * dy)) ** 2
34+
)
35+
36+
37+
def pressure_poisson(nit, p, dx, dy, b):
38+
pn = np.empty_like(p)
39+
pn = p.copy()
40+
41+
for q in range(nit):
42+
pn = p.copy()
43+
p[1:-1, 1:-1] = (
44+
(pn[1:-1, 2:] + pn[1:-1, 0:-2]) * dy**2
45+
+ (pn[2:, 1:-1] + pn[0:-2, 1:-1]) * dx**2
46+
) / (2 * (dx**2 + dy**2)) - dx**2 * dy**2 / (
47+
2 * (dx**2 + dy**2)
48+
) * b[
49+
1:-1, 1:-1
50+
]
51+
52+
p[:, -1] = p[:, -2] # dp/dx = 0 at x = 2
53+
p[0, :] = p[1, :] # dp/dy = 0 at y = 0
54+
p[:, 0] = p[:, 1] # dp/dx = 0 at x = 0
55+
p[-1, :] = 0 # p = 0 at y = 2
56+
57+
58+
def cavity_flow(nx, ny, nt, nit, u, v, dt, dx, dy, p, rho, nu):
59+
un = np.empty_like(u)
60+
vn = np.empty_like(v)
61+
b = np.zeros((ny, nx))
62+
63+
for n in range(nt):
64+
un = u.copy()
65+
vn = v.copy()
66+
67+
build_up_b(b, rho, dt, u, v, dx, dy)
68+
pressure_poisson(nit, p, dx, dy, b)
69+
70+
u[1:-1, 1:-1] = (
71+
un[1:-1, 1:-1]
72+
- un[1:-1, 1:-1] * dt / dx * (un[1:-1, 1:-1] - un[1:-1, 0:-2])
73+
- vn[1:-1, 1:-1] * dt / dy * (un[1:-1, 1:-1] - un[0:-2, 1:-1])
74+
- dt / (2 * rho * dx) * (p[1:-1, 2:] - p[1:-1, 0:-2])
75+
+ nu
76+
* (
77+
dt
78+
/ dx**2
79+
* (un[1:-1, 2:] - 2 * un[1:-1, 1:-1] + un[1:-1, 0:-2])
80+
+ dt
81+
/ dy**2
82+
* (un[2:, 1:-1] - 2 * un[1:-1, 1:-1] + un[0:-2, 1:-1])
83+
)
84+
)
85+
86+
v[1:-1, 1:-1] = (
87+
vn[1:-1, 1:-1]
88+
- un[1:-1, 1:-1] * dt / dx * (vn[1:-1, 1:-1] - vn[1:-1, 0:-2])
89+
- vn[1:-1, 1:-1] * dt / dy * (vn[1:-1, 1:-1] - vn[0:-2, 1:-1])
90+
- dt / (2 * rho * dy) * (p[2:, 1:-1] - p[0:-2, 1:-1])
91+
+ nu
92+
* (
93+
dt
94+
/ dx**2
95+
* (vn[1:-1, 2:] - 2 * vn[1:-1, 1:-1] + vn[1:-1, 0:-2])
96+
+ dt
97+
/ dy**2
98+
* (vn[2:, 1:-1] - 2 * vn[1:-1, 1:-1] + vn[0:-2, 1:-1])
99+
)
100+
)
101+
102+
u[0, :] = 0
103+
u[:, 0] = 0
104+
u[:, -1] = 0
105+
u[-1, :] = 1 # set velocity on cavity lid equal to 1
106+
v[0, :] = 0
107+
v[-1, :] = 0
108+
v[:, 0] = 0
109+
v[:, -1] = 0
Lines changed: 224 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,224 @@
1+
# SPDX-FileCopyrightText: 2017 Lorena A. Barba, Gilbert F. Forsyth.
2+
# SPDX-FileCopyrightText: 2018 Barba, Lorena A., and Forsyth, Gilbert F.
3+
# SPDX-FileCopyrightText: 2021 ETH Zurich and the NPBench authors
4+
# SPDX-FileCopyrightText: 2022 - 2023 Intel Corporation
5+
#
6+
# SPDX-License-Identifier: BSD-3-Clause
7+
8+
"""
9+
CFD Python: the 12 steps to Navier-Stokes equations.
10+
Journal of Open Source Education, 1(9), 21,
11+
https://doi.org/10.21105/jose.00021
12+
"""
13+
14+
import dpnp as np
15+
16+
17+
def build_up_b(rho, dt, dx, dy, u, v):
18+
b = np.zeros_like(u)
19+
b[1:-1, 1:-1] = rho * (
20+
1
21+
/ dt
22+
* (
23+
(u[1:-1, 2:] - u[1:-1, 0:-2]) / (2 * dx)
24+
+ (v[2:, 1:-1] - v[0:-2, 1:-1]) / (2 * dy)
25+
)
26+
- ((u[1:-1, 2:] - u[1:-1, 0:-2]) / (2 * dx)) ** 2
27+
- 2
28+
* (
29+
(u[2:, 1:-1] - u[0:-2, 1:-1])
30+
/ (2 * dy)
31+
* (v[1:-1, 2:] - v[1:-1, 0:-2])
32+
/ (2 * dx)
33+
)
34+
- ((v[2:, 1:-1] - v[0:-2, 1:-1]) / (2 * dy)) ** 2
35+
)
36+
37+
# Periodic BC Pressure @ x = 2
38+
b[1:-1, -1] = rho * (
39+
1
40+
/ dt
41+
* (
42+
(u[1:-1, 0] - u[1:-1, -2]) / (2 * dx)
43+
+ (v[2:, -1] - v[0:-2, -1]) / (2 * dy)
44+
)
45+
- ((u[1:-1, 0] - u[1:-1, -2]) / (2 * dx)) ** 2
46+
- 2
47+
* (
48+
(u[2:, -1] - u[0:-2, -1])
49+
/ (2 * dy)
50+
* (v[1:-1, 0] - v[1:-1, -2])
51+
/ (2 * dx)
52+
)
53+
- ((v[2:, -1] - v[0:-2, -1]) / (2 * dy)) ** 2
54+
)
55+
56+
# Periodic BC Pressure @ x = 0
57+
b[1:-1, 0] = rho * (
58+
1
59+
/ dt
60+
* (
61+
(u[1:-1, 1] - u[1:-1, -1]) / (2 * dx)
62+
+ (v[2:, 0] - v[0:-2, 0]) / (2 * dy)
63+
)
64+
- ((u[1:-1, 1] - u[1:-1, -1]) / (2 * dx)) ** 2
65+
- 2
66+
* (
67+
(u[2:, 0] - u[0:-2, 0])
68+
/ (2 * dy)
69+
* (v[1:-1, 1] - v[1:-1, -1])
70+
/ (2 * dx)
71+
)
72+
- ((v[2:, 0] - v[0:-2, 0]) / (2 * dy)) ** 2
73+
)
74+
75+
return b
76+
77+
78+
def pressure_poisson_periodic(nit, p, dx, dy, b):
79+
pn = np.empty_like(p)
80+
81+
for q in range(nit):
82+
pn = p.copy()
83+
p[1:-1, 1:-1] = (
84+
(pn[1:-1, 2:] + pn[1:-1, 0:-2]) * dy**2
85+
+ (pn[2:, 1:-1] + pn[0:-2, 1:-1]) * dx**2
86+
) / (2 * (dx**2 + dy**2)) - dx**2 * dy**2 / (
87+
2 * (dx**2 + dy**2)
88+
) * b[
89+
1:-1, 1:-1
90+
]
91+
92+
# Periodic BC Pressure @ x = 2
93+
p[1:-1, -1] = (
94+
(pn[1:-1, 0] + pn[1:-1, -2]) * dy**2
95+
+ (pn[2:, -1] + pn[0:-2, -1]) * dx**2
96+
) / (2 * (dx**2 + dy**2)) - dx**2 * dy**2 / (
97+
2 * (dx**2 + dy**2)
98+
) * b[
99+
1:-1, -1
100+
]
101+
102+
# Periodic BC Pressure @ x = 0
103+
p[1:-1, 0] = (
104+
(pn[1:-1, 1] + pn[1:-1, -1]) * dy**2
105+
+ (pn[2:, 0] + pn[0:-2, 0]) * dx**2
106+
) / (2 * (dx**2 + dy**2)) - dx**2 * dy**2 / (
107+
2 * (dx**2 + dy**2)
108+
) * b[
109+
1:-1, 0
110+
]
111+
112+
# Wall boundary conditions, pressure
113+
p[-1, :] = p[-2, :] # dp/dy = 0 at y = 2
114+
p[0, :] = p[1, :] # dp/dy = 0 at y = 0
115+
116+
117+
def channel_flow(nit, u, v, dt, dx, dy, p, rho, nu, F):
118+
udiff = 1
119+
stepcount = 0
120+
121+
while udiff > 0.001:
122+
un = u.copy()
123+
vn = v.copy()
124+
125+
b = build_up_b(rho, dt, dx, dy, u, v)
126+
pressure_poisson_periodic(nit, p, dx, dy, b)
127+
128+
u[1:-1, 1:-1] = (
129+
un[1:-1, 1:-1]
130+
- un[1:-1, 1:-1] * dt / dx * (un[1:-1, 1:-1] - un[1:-1, 0:-2])
131+
- vn[1:-1, 1:-1] * dt / dy * (un[1:-1, 1:-1] - un[0:-2, 1:-1])
132+
- dt / (2 * rho * dx) * (p[1:-1, 2:] - p[1:-1, 0:-2])
133+
+ nu
134+
* (
135+
dt
136+
/ dx**2
137+
* (un[1:-1, 2:] - 2 * un[1:-1, 1:-1] + un[1:-1, 0:-2])
138+
+ dt
139+
/ dy**2
140+
* (un[2:, 1:-1] - 2 * un[1:-1, 1:-1] + un[0:-2, 1:-1])
141+
)
142+
+ F * dt
143+
)
144+
145+
v[1:-1, 1:-1] = (
146+
vn[1:-1, 1:-1]
147+
- un[1:-1, 1:-1] * dt / dx * (vn[1:-1, 1:-1] - vn[1:-1, 0:-2])
148+
- vn[1:-1, 1:-1] * dt / dy * (vn[1:-1, 1:-1] - vn[0:-2, 1:-1])
149+
- dt / (2 * rho * dy) * (p[2:, 1:-1] - p[0:-2, 1:-1])
150+
+ nu
151+
* (
152+
dt
153+
/ dx**2
154+
* (vn[1:-1, 2:] - 2 * vn[1:-1, 1:-1] + vn[1:-1, 0:-2])
155+
+ dt
156+
/ dy**2
157+
* (vn[2:, 1:-1] - 2 * vn[1:-1, 1:-1] + vn[0:-2, 1:-1])
158+
)
159+
)
160+
161+
# Periodic BC u @ x = 2
162+
u[1:-1, -1] = (
163+
un[1:-1, -1]
164+
- un[1:-1, -1] * dt / dx * (un[1:-1, -1] - un[1:-1, -2])
165+
- vn[1:-1, -1] * dt / dy * (un[1:-1, -1] - un[0:-2, -1])
166+
- dt / (2 * rho * dx) * (p[1:-1, 0] - p[1:-1, -2])
167+
+ nu
168+
* (
169+
dt / dx**2 * (un[1:-1, 0] - 2 * un[1:-1, -1] + un[1:-1, -2])
170+
+ dt / dy**2 * (un[2:, -1] - 2 * un[1:-1, -1] + un[0:-2, -1])
171+
)
172+
+ F * dt
173+
)
174+
175+
# Periodic BC u @ x = 0
176+
u[1:-1, 0] = (
177+
un[1:-1, 0]
178+
- un[1:-1, 0] * dt / dx * (un[1:-1, 0] - un[1:-1, -1])
179+
- vn[1:-1, 0] * dt / dy * (un[1:-1, 0] - un[0:-2, 0])
180+
- dt / (2 * rho * dx) * (p[1:-1, 1] - p[1:-1, -1])
181+
+ nu
182+
* (
183+
dt / dx**2 * (un[1:-1, 1] - 2 * un[1:-1, 0] + un[1:-1, -1])
184+
+ dt / dy**2 * (un[2:, 0] - 2 * un[1:-1, 0] + un[0:-2, 0])
185+
)
186+
+ F * dt
187+
)
188+
189+
# Periodic BC v @ x = 2
190+
v[1:-1, -1] = (
191+
vn[1:-1, -1]
192+
- un[1:-1, -1] * dt / dx * (vn[1:-1, -1] - vn[1:-1, -2])
193+
- vn[1:-1, -1] * dt / dy * (vn[1:-1, -1] - vn[0:-2, -1])
194+
- dt / (2 * rho * dy) * (p[2:, -1] - p[0:-2, -1])
195+
+ nu
196+
* (
197+
dt / dx**2 * (vn[1:-1, 0] - 2 * vn[1:-1, -1] + vn[1:-1, -2])
198+
+ dt / dy**2 * (vn[2:, -1] - 2 * vn[1:-1, -1] + vn[0:-2, -1])
199+
)
200+
)
201+
202+
# Periodic BC v @ x = 0
203+
v[1:-1, 0] = (
204+
vn[1:-1, 0]
205+
- un[1:-1, 0] * dt / dx * (vn[1:-1, 0] - vn[1:-1, -1])
206+
- vn[1:-1, 0] * dt / dy * (vn[1:-1, 0] - vn[0:-2, 0])
207+
- dt / (2 * rho * dy) * (p[2:, 0] - p[0:-2, 0])
208+
+ nu
209+
* (
210+
dt / dx**2 * (vn[1:-1, 1] - 2 * vn[1:-1, 0] + vn[1:-1, -1])
211+
+ dt / dy**2 * (vn[2:, 0] - 2 * vn[1:-1, 0] + vn[0:-2, 0])
212+
)
213+
)
214+
215+
# Wall BC: u,v = 0 @ y = 0,2
216+
u[0, :] = 0
217+
u[-1, :] = 0
218+
v[0, :] = 0
219+
v[-1, :] = 0
220+
221+
udiff = (np.sum(u) - np.sum(un)) / np.sum(u)
222+
stepcount += 1
223+
224+
return stepcount
Lines changed: 14 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,14 @@
1+
# SPDX-FileCopyrightText: 2023 Stefan Behnel, Robert Bradshaw,
2+
# Dag Sverre Seljebotn, Greg Ewing, William Stein, Gabriel Gellner, et al.
3+
# SPDX-FileCopyrightText: 2021 ETH Zurich and the NPBench authors
4+
# SPDX-FileCopyrightText: 2022 - 2023 Intel Corporation
5+
#
6+
# SPDX-License-Identifier: BSD-3-Clause
7+
8+
# https://cython.readthedocs.io/en/latest/src/userguide/numpy_tutorial.html
9+
10+
import dpnp as np
11+
12+
13+
def compute(array_1, array_2, a, b, c):
14+
return np.clip(array_1, 2, 10) * a + array_2 * b + c

0 commit comments

Comments
 (0)