Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
rhs_maker_device.F90
Go to the documentation of this file.
1! Copyright (c) 2022, The Neko Authors
2! All rights reserved.
3!
4! Redistribution and use in source and binary forms, with or without
5! modification, are permitted provided that the following conditions
6! are met:
7!
8! * Redistributions of source code must retain the above copyright
9! notice, this list of conditions and the following disclaimer.
10!
11! * Redistributions in binary form must reproduce the above
12! copyright notice, this list of conditions and the following
13! disclaimer in the documentation and/or other materials provided
14! with the distribution.
15!
16! * Neither the name of the authors nor the names of its
17! contributors may be used to endorse or promote products derived
18! from this software without specific prior written permission.
19!
20! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
21! "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
22! LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
23! FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
24! COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
25! INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
26! BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
27! LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
28! CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
29! LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
30! ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
31! POSSIBILITY OF SUCH DAMAGE.
32!
34 use rhs_maker
35 use device
36 use utils
38 use field, only : field_t
39 use num_types, only : rp, c_rp
40 use, intrinsic :: iso_c_binding, only : c_ptr
41 implicit none
42 private
43
45 contains
46 procedure, nopass :: compute_fluid => rhs_maker_sumab_device
48
49 type, public, extends(rhs_maker_ext_t) :: rhs_maker_ext_device_t
50 contains
51 procedure, nopass :: compute_fluid => rhs_maker_ext_device
52 procedure, nopass :: compute_scalar => scalar_rhs_maker_ext_device
54
55 type, public, extends(rhs_maker_bdf_t) :: rhs_maker_bdf_device_t
56 contains
57 procedure, nopass :: compute_fluid => rhs_maker_bdf_device
58 procedure, nopass :: compute_scalar => scalar_rhs_maker_bdf_device
60
61 type, public, extends(rhs_maker_oifs_t) :: rhs_maker_oifs_device_t
62 contains
63 procedure, nopass :: compute_fluid => rhs_maker_oifs_device
64 procedure, nopass :: compute_scalar => scalar_rhs_maker_oifs_device
66
67#ifdef HAVE_HIP
68 interface
69 subroutine rhs_maker_sumab_hip(u_d, v_d, w_d, uu_d, vv_d, ww_d, &
70 uulag1, uulag2, vvlag1, vvlag2, wwlag1, wwlag2, ab1, ab2, ab3, nab, n)&
71 bind(c, name = 'rhs_maker_sumab_hip')
72 use, intrinsic :: iso_c_binding
73 import c_rp
74 type(c_ptr), value :: u_d, v_d, w_d, uu_d, vv_d, ww_d
75 type(c_ptr), value :: uulag1, uulag2, vvlag1, vvlag2, wwlag1, wwlag2
76 real(c_rp) :: ab1, ab2, ab3
77 integer(c_int) :: nab, n
78 end subroutine rhs_maker_sumab_hip
79 end interface
80
81 interface
82 subroutine rhs_maker_ext_hip(abx1_d, aby1_d, abz1_d, &
83 abx2_d, aby2_d, abz2_d, &
84 bfx_d, bfy_d, bfz_d, &
85 rho, ab1, ab2, ab3, n) &
86 bind(c, name = 'rhs_maker_ext_hip')
87 use, intrinsic :: iso_c_binding
88 import c_rp
89 type(c_ptr), value :: abx1_d, aby1_d, abz1_d
90 type(c_ptr), value :: abx2_d, aby2_d, abz2_d
91 type(c_ptr), value :: bfx_d, bfy_d, bfz_d
92 real(c_rp) :: rho, ab1, ab2, ab3
93 integer(c_int) :: n
94 end subroutine rhs_maker_ext_hip
95 end interface
96
97 interface
98 subroutine scalar_rhs_maker_ext_hip(fs_lag_d, fs_laglag_d, fs_d, rho, &
99 ext1, ext2, ext3, n) &
100 bind(c, name = 'scalar_rhs_maker_ext_hip')
101 use, intrinsic :: iso_c_binding
102 import c_rp
103 type(c_ptr), value :: fs_lag_d, fs_laglag_d, fs_d
104 real(c_rp) :: rho, ext1, ext2, ext3
105 integer(c_int) :: n
106 end subroutine scalar_rhs_maker_ext_hip
107 end interface
108
109 interface
110 subroutine rhs_maker_bdf_hip(ulag1_d, ulag2_d, vlag1_d, vlag2_d, &
111 wlag1_d, wlag2_d, bfx_d, bfy_d, bfz_d, u_d, v_d, w_d, B_d, &
112 rho, dt, bd2, bd3, bd4, nbd, n) bind(c, name = 'rhs_maker_bdf_hip')
113 use, intrinsic :: iso_c_binding
114 import c_rp
115 type(c_ptr), value :: ulag1_d, ulag2_d, vlag1_d
116 type(c_ptr), value :: vlag2_d, wlag1_d, wlag2_d
117 type(c_ptr), value :: bfx_d, bfy_d, bfz_d, u_d, v_d, w_d, B_d
118 reaL(c_rp) :: rho, dt, bd2, bd3, bd4
119 integer(c_int) :: nbd, n
120 end subroutine rhs_maker_bdf_hip
121 end interface
122
123 interface
124 subroutine scalar_rhs_maker_bdf_hip(s_lag_d, s_laglag_d, fs_d, s_d, B_d, &
125 rho, dt, bd2, bd3, bd4, nbd, n) &
126 bind(c, name = 'scalar_rhs_maker_bdf_hip')
127 use, intrinsic :: iso_c_binding
128 import c_rp
129 type(c_ptr), value :: s_lag_d, s_laglag_d
130 type(c_ptr), value :: fs_d, s_d, B_d
131 reaL(c_rp) :: rho, dt, bd2, bd3, bd4
132 integer(c_int) :: nbd, n
133 end subroutine scalar_rhs_maker_bdf_hip
134 end interface
135
136 interface
137 subroutine rhs_maker_oifs_hip(phi_x_d, phi_y_d, phi_z_d, bf_x_d, &
138 bf_y_d, bf_z_d, rho, dt, n) bind(c, name = 'rhs_maker_oifs_hip')
139 use, intrinsic :: iso_c_binding
140 import c_rp
141 type(c_ptr), value :: bf_x_d, bf_y_d, bf_z_d
142 type(c_ptr), value :: phi_x_d, phi_y_d, phi_z_d
143 reaL(c_rp) :: rho, dt
144 integer(c_int) :: n
145 end subroutine rhs_maker_oifs_hip
146 end interface
147
148 interface
149 subroutine scalar_rhs_maker_oifs_hip(phi_s_d, bf_s_d, rho, dt, n) &
150 bind(c, name = 'scalar_rhs_maker_oifs_hip')
151 use, intrinsic :: iso_c_binding
152 import c_rp
153 type(c_ptr), value :: bf_s_d
154 type(c_ptr), value :: phi_s_d
155 reaL(c_rp) :: rho, dt
156 integer(c_int) :: n
157 end subroutine scalar_rhs_maker_oifs_hip
158 end interface
159#elif HAVE_CUDA
160 interface
161 subroutine rhs_maker_sumab_cuda(u_d, v_d, w_d, uu_d, vv_d, ww_d, &
162 uulag1, uulag2, vvlag1, vvlag2, wwlag1, wwlag2, ab1, ab2, ab3, nab, n)&
163 bind(c, name = 'rhs_maker_sumab_cuda')
164 use, intrinsic :: iso_c_binding
165 import c_rp
166 type(c_ptr), value :: u_d, v_d, w_d, uu_d, vv_d, ww_d
167 type(c_ptr), value :: uulag1, uulag2, vvlag1, vvlag2, wwlag1, wwlag2
168 real(c_rp) :: ab1, ab2, ab3
169 integer(c_int) :: nab, n
170 end subroutine rhs_maker_sumab_cuda
171 end interface
172
173 interface
174 subroutine rhs_maker_ext_cuda(abx1_d, aby1_d, abz1_d, &
175 abx2_d, aby2_d, abz2_d, &
176 bfx_d, bfy_d, bfz_d, &
177 rho, ab1, ab2, ab3, n) &
178 bind(c, name = 'rhs_maker_ext_cuda')
179 use, intrinsic :: iso_c_binding
180 import c_rp
181 type(c_ptr), value :: abx1_d, aby1_d, abz1_d
182 type(c_ptr), value :: abx2_d, aby2_d, abz2_d
183 type(c_ptr), value :: bfx_d, bfy_d, bfz_d
184 real(c_rp) :: rho, ab1, ab2, ab3
185 integer(c_int) :: n
186 end subroutine rhs_maker_ext_cuda
187 end interface
188
189 interface
190 subroutine scalar_rhs_maker_ext_cuda(fs_lag_d, fs_laglag_d, fs_d, rho, &
191 ext1, ext2, ext3, n) &
192 bind(c, name = 'scalar_rhs_maker_ext_cuda')
193 use, intrinsic :: iso_c_binding
194 import c_rp
195 type(c_ptr), value :: fs_lag_d, fs_laglag_d, fs_d
196 real(c_rp) :: rho, ext1, ext2, ext3
197 integer(c_int) :: n
198 end subroutine scalar_rhs_maker_ext_cuda
199 end interface
200
201 interface
202 subroutine rhs_maker_bdf_cuda(ulag1_d, ulag2_d, vlag1_d, vlag2_d, &
203 wlag1_d, wlag2_d, bfx_d, bfy_d, bfz_d, u_d, v_d, w_d, B_d, &
204 rho, dt, bd2, bd3, bd4, nbd, n) bind(c, name = 'rhs_maker_bdf_cuda')
205 use, intrinsic :: iso_c_binding
206 import c_rp
207 type(c_ptr), value :: ulag1_d, ulag2_d, vlag1_d
208 type(c_ptr), value :: vlag2_d, wlag1_d, wlag2_d
209 type(c_ptr), value :: bfx_d, bfy_d, bfz_d, u_d, v_d, w_d, B_d
210 reaL(c_rp) :: rho, dt, bd2, bd3, bd4
211 integer(c_int) :: nbd, n
212 end subroutine rhs_maker_bdf_cuda
213 end interface
214
215 interface
216 subroutine scalar_rhs_maker_bdf_cuda(s_lag_d, s_laglag_d, fs_d, s_d, B_d, &
217 rho, dt, bd2, bd3, bd4, nbd, n) &
218 bind(c, name = 'scalar_rhs_maker_bdf_cuda')
219 use, intrinsic :: iso_c_binding
220 import c_rp
221 type(c_ptr), value :: s_lag_d, s_laglag_d
222 type(c_ptr), value :: fs_d, s_d, B_d
223 reaL(c_rp) :: rho, dt, bd2, bd3, bd4
224 integer(c_int) :: nbd, n
225 end subroutine scalar_rhs_maker_bdf_cuda
226 end interface
227
228 interface
229 subroutine rhs_maker_oifs_cuda(phi_x_d, phi_y_d, phi_z_d, bf_x_d, &
230 bf_y_d, bf_z_d, rho, dt, n) bind(c, name = 'rhs_maker_oifs_cuda')
231 use, intrinsic :: iso_c_binding
232 import c_rp
233 type(c_ptr), value :: bf_x_d, bf_y_d, bf_z_d
234 type(c_ptr), value :: phi_x_d, phi_y_d, phi_z_d
235 reaL(c_rp) :: rho, dt
236 integer(c_int) :: n
237 end subroutine rhs_maker_oifs_cuda
238 end interface
239
240 interface
241 subroutine scalar_rhs_maker_oifs_cuda(phi_s_d, bf_s_d, rho, dt, n) &
242 bind(c, name = 'scalar_rhs_maker_oifs_cuda')
243 use, intrinsic :: iso_c_binding
244 import c_rp
245 type(c_ptr), value :: bf_s_d
246 type(c_ptr), value :: phi_s_d
247 reaL(c_rp) :: rho, dt
248 integer(c_int) :: n
249 end subroutine scalar_rhs_maker_oifs_cuda
250 end interface
251#elif HAVE_OPENCL
252 interface
253 subroutine rhs_maker_sumab_opencl(u_d, v_d, w_d, uu_d, vv_d, ww_d, &
254 uulag1, uulag2, vvlag1, vvlag2, wwlag1, wwlag2, ab1, ab2, ab3, nab, n)&
255 bind(c, name = 'rhs_maker_sumab_opencl')
256 use, intrinsic :: iso_c_binding
257 import c_rp
258 type(c_ptr), value :: u_d, v_d, w_d, uu_d, vv_d, ww_d
259 type(c_ptr), value :: uulag1, uulag2, vvlag1, vvlag2, wwlag1, wwlag2
260 real(c_rp) :: ab1, ab2, ab3
261 integer(c_int) :: nab, n
262 end subroutine rhs_maker_sumab_opencl
263 end interface
264
265 interface
266 subroutine rhs_maker_ext_opencl(abx1_d, aby1_d, abz1_d, &
267 abx2_d, aby2_d, abz2_d, &
268 bfx_d, bfy_d, bfz_d, &
269 rho, ab1, ab2, ab3, n) &
270 bind(c, name = 'rhs_maker_ext_opencl')
271 use, intrinsic :: iso_c_binding
272 import c_rp
273 type(c_ptr), value :: abx1_d, aby1_d, abz1_d
274 type(c_ptr), value :: abx2_d, aby2_d, abz2_d
275 type(c_ptr), value :: bfx_d, bfy_d, bfz_d
276 real(c_rp) :: rho, ab1, ab2, ab3
277 integer(c_int) :: n
278 end subroutine rhs_maker_ext_opencl
279 end interface
280
281 interface
282 subroutine scalar_rhs_maker_ext_opencl(fs_lag_d, fs_laglag_d, fs_d, rho, &
283 ext1, ext2, ext3, n) &
284 bind(c, name = 'scalar_rhs_maker_ext_opencl')
285 use, intrinsic :: iso_c_binding
286 import c_rp
287 type(c_ptr), value :: fs_lag_d, fs_laglag_d, fs_d
288 real(c_rp) :: rho, ext1, ext2, ext3
289 integer(c_int) :: n
290 end subroutine scalar_rhs_maker_ext_opencl
291 end interface
292
293 interface
294 subroutine rhs_maker_bdf_opencl(ulag1_d, ulag2_d, vlag1_d, vlag2_d, &
295 wlag1_d, wlag2_d, bfx_d, bfy_d, bfz_d, u_d, v_d, w_d, B_d, &
296 rho, dt, bd2, bd3, bd4, nbd, n) bind(c, name = 'rhs_maker_bdf_opencl')
297 use, intrinsic :: iso_c_binding
298 import c_rp
299 type(c_ptr), value :: ulag1_d, ulag2_d, vlag1_d
300 type(c_ptr), value :: vlag2_d, wlag1_d, wlag2_d
301 type(c_ptr), value :: bfx_d, bfy_d, bfz_d, u_d, v_d, w_d, B_d
302 reaL(c_rp) :: rho, dt, bd2, bd3, bd4
303 integer(c_int) :: nbd, n
304 end subroutine rhs_maker_bdf_opencl
305 end interface
306
307 interface
308 subroutine scalar_rhs_maker_bdf_opencl(s_lag_d, s_laglag_d, fs_d, s_d, B_d, &
309 rho, dt, bd2, bd3, bd4, nbd, n) &
310 bind(c, name = 'scalar_rhs_maker_bdf_opencl')
311 use, intrinsic :: iso_c_binding
312 import c_rp
313 type(c_ptr), value :: s_lag_d, s_laglag_d
314 type(c_ptr), value :: fs_d, s_d, B_d
315 reaL(c_rp) :: rho, dt, bd2, bd3, bd4
316 integer(c_int) :: nbd, n
317 end subroutine scalar_rhs_maker_bdf_opencl
318 end interface
319
320 interface
321 subroutine rhs_maker_oifs_opencl(phi_x_d, phi_y_d, phi_z_d, bf_x_d, &
322 bf_y_d, bf_z_d, rho, dt, n) bind(c, name = 'rhs_maker_oifs_opencl')
323 use, intrinsic :: iso_c_binding
324 import c_rp
325 type(c_ptr), value :: bf_x_d, bf_y_d, bf_z_d
326 type(c_ptr), value :: phi_x_d, phi_y_d, phi_z_d
327 reaL(c_rp) :: rho, dt
328 integer(c_int) :: n
329 end subroutine rhs_maker_oifs_opencl
330 end interface
331
332 interface
333 subroutine scalar_rhs_maker_oifs_opencl(phi_s_d, bf_s_d, rho, dt, n) &
334 bind(c, name = 'scalar_rhs_maker_oifs_opencl')
335 use, intrinsic :: iso_c_binding
336 import c_rp
337 type(c_ptr), value :: bf_s_d
338 type(c_ptr), value :: phi_s_d
339 reaL(c_rp) :: rho, dt
340 integer(c_int) :: n
341 end subroutine scalar_rhs_maker_oifs_opencl
342 end interface
343#endif
344
345contains
346
347 subroutine rhs_maker_sumab_device(u, v, w, uu, vv, ww, uulag, vvlag, wwlag, ab, nab)
348 type(field_t), intent(inout) :: u,v, w
349 type(field_t), intent(inout) :: uu, vv, ww
350 type(field_series_t), intent(inout) :: uulag, vvlag, wwlag
351 real(kind=rp), dimension(3), intent(in) :: ab
352 integer, intent(in) :: nab
353
354#ifdef HAVE_HIP
355 call rhs_maker_sumab_hip(u%x_d, v%x_d, w%x_d, uu%x_d, vv%x_d, ww%x_d, &
356 uulag%lf(1)%x_d, uulag%lf(2)%x_d, vvlag%lf(1)%x_d, vvlag%lf(2)%x_d, &
357 wwlag%lf(1)%x_d, wwlag%lf(2)%x_d, ab(1), ab(2), ab(3), nab, &
358 uu%dof%size())
359#elif HAVE_CUDA
360 call rhs_maker_sumab_cuda(u%x_d, v%x_d, w%x_d, uu%x_d, vv%x_d, ww%x_d, &
361 uulag%lf(1)%x_d, uulag%lf(2)%x_d, vvlag%lf(1)%x_d, vvlag%lf(2)%x_d, &
362 wwlag%lf(1)%x_d, wwlag%lf(2)%x_d, ab(1), ab(2), ab(3), nab, &
363 uu%dof%size())
364#elif HAVE_OPENCL
365 call rhs_maker_sumab_opencl(u%x_d, v%x_d, w%x_d, uu%x_d, vv%x_d, ww%x_d, &
366 uulag%lf(1)%x_d, uulag%lf(2)%x_d, vvlag%lf(1)%x_d, vvlag%lf(2)%x_d, &
367 wwlag%lf(1)%x_d, wwlag%lf(2)%x_d, ab(1), ab(2), ab(3), nab, &
368 uu%dof%size())
369#endif
370
371 end subroutine rhs_maker_sumab_device
372
373 subroutine rhs_maker_ext_device(fx_lag, fy_lag, fz_lag, &
374 fx_laglag, fy_laglag, fz_laglag, fx, fy, fz, &
375 rho, ext_coeffs, n)
376 type(field_t), intent(inout) :: fx_lag, fy_lag, fz_lag
377 type(field_t), intent(inout) :: fx_laglag, fy_laglag, fz_laglag
378 real(kind=rp), intent(in) :: rho, ext_coeffs(4)
379 integer, intent(in) :: n
380 real(kind=rp), intent(inout) :: fx(n), fy(n), fz(n)
381 type(c_ptr) :: fx_d, fy_d, fz_d
382
383 fx_d = device_get_ptr(fx)
384 fy_d = device_get_ptr(fy)
385 fz_d = device_get_ptr(fz)
386
387#ifdef HAVE_HIP
388 call rhs_maker_ext_hip(fx_lag%x_d, fy_lag%x_d, fz_lag%x_d, &
389 fx_laglag%x_d, fy_laglag%x_d, fz_laglag%x_d, &
390 fx_d, fy_d, fz_d, rho, &
391 ext_coeffs(1), ext_coeffs(2), ext_coeffs(3), n)
392#elif HAVE_CUDA
393 call rhs_maker_ext_cuda(fx_lag%x_d, fy_lag%x_d, fz_lag%x_d, &
394 fx_laglag%x_d, fy_laglag%x_d, fz_laglag%x_d, &
395 fx_d, fy_d, fz_d, rho, &
396 ext_coeffs(1), ext_coeffs(2), ext_coeffs(3), n)
397#elif HAVE_OPENCL
398 call rhs_maker_ext_opencl(fx_lag%x_d, fy_lag%x_d, fz_lag%x_d, &
399 fx_laglag%x_d, fy_laglag%x_d, fz_laglag%x_d, &
400 fx_d, fy_d, fz_d, rho, &
401 ext_coeffs(1), ext_coeffs(2), ext_coeffs(3), n)
402#endif
403
404 end subroutine rhs_maker_ext_device
405
406 subroutine scalar_rhs_maker_ext_device(fs_lag, fs_laglag, fs, &
407 rho, ext_coeffs, n)
408 type(field_t), intent(inout) :: fs_lag
409 type(field_t), intent(inout) :: fs_laglag
410 real(kind=rp), intent(in) :: rho, ext_coeffs(4)
411 integer, intent(in) :: n
412 real(kind=rp), intent(inout) :: fs(n)
413 type(c_ptr) :: fs_d
414
415 fs_d = device_get_ptr(fs)
416
417#ifdef HAVE_HIP
418 call scalar_rhs_maker_ext_hip(fs_lag%x_d, fs_laglag%x_d, fs_d, rho, &
419 ext_coeffs(1), ext_coeffs(2), ext_coeffs(3), n)
420#elif HAVE_CUDA
421 call scalar_rhs_maker_ext_cuda(fs_lag%x_d, fs_laglag%x_d, fs_d, rho, &
422 ext_coeffs(1), ext_coeffs(2), ext_coeffs(3), n)
423#elif HAVE_OPENCL
424 call scalar_rhs_maker_ext_opencl(fs_lag%x_d, fs_laglag%x_d, fs_d, rho, &
425 ext_coeffs(1), ext_coeffs(2), ext_coeffs(3), n)
426#endif
427
428 end subroutine scalar_rhs_maker_ext_device
429
430 subroutine rhs_maker_bdf_device(ulag, vlag, wlag, bfx, bfy, bfz, &
431 u, v, w, B, rho, dt, bd, nbd, n)
432 integer, intent(in) :: n, nbd
433 type(field_t), intent(in) :: u, v, w
434 type(field_series_t), intent(in) :: ulag, vlag, wlag
435 real(kind=rp), intent(inout) :: bfx(n), bfy(n), bfz(n)
436 real(kind=rp), intent(in) :: b(n)
437 real(kind=rp), intent(in) :: dt, rho, bd(4)
438 type(c_ptr) :: bfx_d, bfy_d, bfz_d, b_d
439
440 bfx_d = device_get_ptr(bfx)
441 bfy_d = device_get_ptr(bfy)
442 bfz_d = device_get_ptr(bfz)
443 b_d = device_get_ptr(b)
444
445#ifdef HAVE_HIP
446 call rhs_maker_bdf_hip(ulag%lf(1)%x_d, ulag%lf(2)%x_d, &
447 vlag%lf(1)%x_d, vlag%lf(2)%x_d, &
448 wlag%lf(1)%x_d, wlag%lf(2)%x_d, &
449 bfx_d, bfy_d, bfz_d, u%x_d, v%x_d, w%x_d, &
450 b_d, rho, dt, bd(2), bd(3), bd(4), nbd, n)
451#elif HAVE_CUDA
452 call rhs_maker_bdf_cuda(ulag%lf(1)%x_d, ulag%lf(2)%x_d, &
453 vlag%lf(1)%x_d, vlag%lf(2)%x_d, &
454 wlag%lf(1)%x_d, wlag%lf(2)%x_d, &
455 bfx_d, bfy_d, bfz_d, u%x_d, v%x_d, w%x_d, &
456 b_d, rho, dt, bd(2), bd(3), bd(4), nbd, n)
457#elif HAVE_OPENCL
458 call rhs_maker_bdf_opencl(ulag%lf(1)%x_d, ulag%lf(2)%x_d, &
459 vlag%lf(1)%x_d, vlag%lf(2)%x_d, &
460 wlag%lf(1)%x_d, wlag%lf(2)%x_d, &
461 bfx_d, bfy_d, bfz_d, u%x_d, v%x_d, w%x_d, &
462 b_d, rho, dt, bd(2), bd(3), bd(4), nbd, n)
463#endif
464
465 end subroutine rhs_maker_bdf_device
466
467 subroutine scalar_rhs_maker_bdf_device(s_lag, fs, s, B, rho, dt, &
468 bd, nbd, n)
469 integer, intent(in) :: n, nbd
470 type(field_t), intent(in) :: s
471 type(field_series_t), intent(in) :: s_lag
472 real(kind=rp), intent(inout) :: fs(n)
473 real(kind=rp), intent(in) :: b(n)
474 real(kind=rp), intent(in) :: dt, rho, bd(4)
475 type(c_ptr) :: fs_d, b_d
476
477 fs_d = device_get_ptr(fs)
478 b_d = device_get_ptr(b)
479
480#ifdef HAVE_HIP
481 call scalar_rhs_maker_bdf_hip(s_lag%lf(1)%x_d, s_lag%lf(2)%x_d, &
482 fs_d, s%x_d, b_d, rho, dt, bd(2), bd(3), bd(4), &
483 nbd, n)
484#elif HAVE_CUDA
485 call scalar_rhs_maker_bdf_cuda(s_lag%lf(1)%x_d, s_lag%lf(2)%x_d, &
486 fs_d, s%x_d, b_d, rho, dt, bd(2), bd(3), bd(4), &
487 nbd, n)
488#elif HAVE_OPENCL
489 call scalar_rhs_maker_bdf_opencl(s_lag%lf(1)%x_d, s_lag%lf(2)%x_d, &
490 fs_d, s%x_d, b_d, rho, dt, bd(2), bd(3), bd(4), &
491 nbd, n)
492#endif
493
494 end subroutine scalar_rhs_maker_bdf_device
495
496 subroutine rhs_maker_oifs_device(phi_x, phi_y, phi_z, bf_x, bf_y, bf_z, &
497 rho, dt, n)
498 real(kind=rp), intent(in) :: rho, dt
499 integer, intent(in) :: n
500 real(kind=rp), intent(inout) :: bf_x(n), bf_y(n), bf_z(n)
501 real(kind=rp), intent(inout) :: phi_x(n), phi_y(n), phi_z(n)
502
503 type(c_ptr) :: phi_x_d, phi_y_d, phi_z_d, bf_x_d, bf_y_d, bf_z_d
504
505 phi_x_d = device_get_ptr(phi_x)
506 phi_y_d = device_get_ptr(phi_y)
507 phi_z_d = device_get_ptr(phi_z)
508 bf_x_d = device_get_ptr(bf_x)
509 bf_y_d = device_get_ptr(bf_y)
510 bf_z_d = device_get_ptr(bf_z)
511
512#ifdef HAVE_HIP
513 call rhs_maker_oifs_hip(phi_x_d, phi_y_d, phi_z_d, &
514 bf_x_d, bf_y_d, bf_z_d, rho, dt, n)
515#elif HAVE_CUDA
516 call rhs_maker_oifs_cuda(phi_x_d, phi_y_d, phi_z_d, &
517 bf_x_d, bf_y_d, bf_z_d, rho, dt, n)
518#elif HAVE_OPENCL
519 call rhs_maker_oifs_opencl(phi_x_d, phi_y_d, phi_z_d, &
520 bf_x_d, bf_y_d, bf_z_d, rho, dt, n)
521#endif
522
523 end subroutine rhs_maker_oifs_device
524
525 subroutine scalar_rhs_maker_oifs_device(phi_s,bf_s, rho, dt, n)
526 real(kind=rp), intent(in) :: rho, dt
527 integer, intent(in) :: n
528 real(kind=rp), intent(inout) :: bf_s(n)
529 real(kind=rp), intent(inout) :: phi_s(n)
530 type(c_ptr) :: phi_s_d, bf_s_d
531
532 phi_s_d = device_get_ptr(phi_s)
533 bf_s_d = device_get_ptr(bf_s)
534
535#ifdef HAVE_HIP
536 call scalar_rhs_maker_oifs_hip(phi_s_d, bf_s_d, rho, dt, n)
537#elif HAVE_CUDA
538 call scalar_rhs_maker_oifs_cuda(phi_s_d, bf_s_d, rho, dt, n)
539#elif HAVE_OPENCL
540 call scalar_rhs_maker_oifs_opencl(phi_s_d, bf_s_d, rho, dt, n)
541#endif
542
543 end subroutine scalar_rhs_maker_oifs_device
544
545end module rhs_maker_device
Return the device pointer for an associated Fortran array.
Definition device.F90:96
Device abstraction, common interface for various accelerators.
Definition device.F90:34
Contains the field_serties_t type.
Defines a field.
Definition field.f90:34
integer, parameter, public c_rp
Definition num_types.f90:13
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
subroutine scalar_rhs_maker_bdf_device(s_lag, fs, s, b, rho, dt, bd, nbd, n)
subroutine rhs_maker_sumab_device(u, v, w, uu, vv, ww, uulag, vvlag, wwlag, ab, nab)
subroutine rhs_maker_oifs_device(phi_x, phi_y, phi_z, bf_x, bf_y, bf_z, rho, dt, n)
subroutine scalar_rhs_maker_ext_device(fs_lag, fs_laglag, fs, rho, ext_coeffs, n)
subroutine rhs_maker_ext_device(fx_lag, fy_lag, fz_lag, fx_laglag, fy_laglag, fz_laglag, fx, fy, fz, rho, ext_coeffs, n)
subroutine scalar_rhs_maker_oifs_device(phi_s, bf_s, rho, dt, n)
subroutine rhs_maker_bdf_device(ulag, vlag, wlag, bfx, bfy, bfz, u, v, w, b, rho, dt, bd, nbd, n)
Routines to generate the right-hand sides for the convection-diffusion equation. Employs the EXT/BDF ...
Definition rhs_maker.f90:38
Utilities.
Definition utils.f90:35
void scalar_rhs_maker_bdf_opencl(void *s_lag, void *s_laglag, void *fs, void *s, void *B, real *rho, real *dt, real *bd2, real *bd3, real *bd4, int *nbd, int *n)
Definition rhs_maker.c:201
void rhs_maker_sumab_opencl(void *u, void *v, void *w, void *uu, void *vv, void *ww, void *ulag1, void *ulag2, void *vlag1, void *vlag2, void *wlag1, void *wlag2, real *ext1, real *ext2, real *ext3, int *nab, int *n)
Definition rhs_maker.c:49
void scalar_rhs_maker_oifs_opencl(void *phi_s, void *bf_s, real *rho, real *dt, int *n)
Definition rhs_maker.c:267
void rhs_maker_ext_opencl(void *abx1, void *aby1, void *abz1, void *abx2, void *aby2, void *abz2, void *bfx, void *bfy, void *bfz, real *rho, real *ext1, real *ext2, real *ext3, int *n)
Definition rhs_maker.c:89
void rhs_maker_oifs_opencl(void *phi_x, void *phi_y, void *phi_z, void *bf_x, void *bf_y, void *bf_z, real *rho, real *dt, int *n)
Definition rhs_maker.c:236
void scalar_rhs_maker_ext_opencl(void *fs_lag, void *fs_laglag, void *fs, real *rho, real *ext1, real *ext2, real *ext3, int *n)
Definition rhs_maker.c:126
void rhs_maker_bdf_opencl(void *ulag1, void *ulag2, void *vlag1, void *vlag2, void *wlag1, void *wlag2, void *bfx, void *bfy, void *bfz, void *u, void *v, void *w, void *B, real *rho, real *dt, real *bd2, real *bd3, real *bd4, int *nbd, int *n)
Definition rhs_maker.c:156
void scalar_rhs_maker_oifs_cuda(void *phi_s, void *bf_s, real *rho, real *dt, int *n)
Definition rhs_maker.cu:157
void scalar_rhs_maker_bdf_cuda(void *s_lag, void *s_laglag, void *fs, void *s, void *B, real *rho, real *dt, real *bd2, real *bd3, real *bd4, int *nbd, int *n)
Definition rhs_maker.cu:120
void rhs_maker_oifs_cuda(void *phi_x, void *phi_y, void *phi_z, void *bf_x, void *bf_y, void *bf_z, real *rho, real *dt, int *n)
Definition rhs_maker.cu:140
void rhs_maker_ext_cuda(void *abx1, void *aby1, void *abz1, void *abx2, void *aby2, void *abz2, void *bfx, void *bfy, void *bfz, real *rho, real *ab1, real *ab2, real *ab3, int *n)
Definition rhs_maker.cu:63
void rhs_maker_bdf_cuda(void *ulag1, void *ulag2, void *vlag1, void *vlag2, void *wlag1, void *wlag2, void *bfx, void *bfy, void *bfz, void *u, void *v, void *w, void *B, real *rho, real *dt, real *bd2, real *bd3, real *bd4, int *nbd, int *n)
Definition rhs_maker.cu:98
void scalar_rhs_maker_ext_cuda(void *fs_lag, void *fs_laglag, void *fs, real *rho, real *ext1, real *ext2, real *ext3, int *n)
Definition rhs_maker.cu:82
void rhs_maker_sumab_cuda(void *u, void *v, void *w, void *uu, void *vv, void *ww, void *ulag1, void *ulag2, void *vlag1, void *vlag2, void *wlag1, void *wlag2, real *ab1, real *ab2, real *ab3, int *nab, int *n)
Definition rhs_maker.cu:44
Stores a series (sequence) of fields, logically connected to a base field, and arranged according to ...
Abstract type to add contributions to F from lagged BD terms.
Definition rhs_maker.f90:59
Abstract type to sum up contributions to kth order extrapolation scheme.
Definition rhs_maker.f90:52
Abstract type to add contributions of kth order OIFS scheme.
Definition rhs_maker.f90:66
Abstract type to compute extrapolated velocity field for the pressure equation.
Definition rhs_maker.f90:46