Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
device_zero_dirichlet.F90
Go to the documentation of this file.
1! Copyright (c) 2021-2025, 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 utils, only : neko_error
35 use, intrinsic :: iso_c_binding, only : c_ptr, c_int
36 private
37
38#ifdef HAVE_HIP
39 interface
40 subroutine hip_zero_dirichlet_apply_scalar(msk, x, m, strm) &
41 bind(c, name = 'hip_zero_dirichlet_apply_scalar')
42 use, intrinsic :: iso_c_binding
43 implicit none
44 integer(c_int) :: m
45 type(c_ptr), value :: msk, x, strm
47 end interface
48
49 interface
50 subroutine hip_zero_dirichlet_apply_vector(msk, x, y, z, m, strm) &
51 bind(c, name = 'hip_zero_dirichlet_apply_vector')
52 use, intrinsic :: iso_c_binding
53 implicit none
54 integer(c_int) :: m
55 type(c_ptr), value :: msk, x, y, z, strm
57 end interface
58#elif HAVE_CUDA
59 interface
60 subroutine cuda_zero_dirichlet_apply_scalar(msk, x, m, strm) &
61 bind(c, name = 'cuda_zero_dirichlet_apply_scalar')
62 use, intrinsic :: iso_c_binding
63 implicit none
64 integer(c_int) :: m
65 type(c_ptr), value :: msk, x, strm
67 end interface
68
69 interface
70 subroutine cuda_zero_dirichlet_apply_vector(msk, x, y, z, m, strm) &
71 bind(c, name = 'cuda_zero_dirichlet_apply_vector')
72 use, intrinsic :: iso_c_binding
73 implicit none
74 integer(c_int) :: m
75 type(c_ptr), value :: msk, x, y, z, strm
77 end interface
78#elif HAVE_OPENCL
79 interface
80 subroutine opencl_zero_dirichlet_apply_scalar(msk, x, m, strm) &
81 bind(c, name = 'opencl_zero_dirichlet_apply_scalar')
82 use, intrinsic :: iso_c_binding
83 implicit none
84 integer(c_int) :: m
85 type(c_ptr), value :: msk, x, strm
87 end interface
88
89 interface
90 subroutine opencl_zero_dirichlet_apply_vector(msk, x, y, z, m, strm) &
91 bind(c, name = 'opencl_zero_dirichlet_apply_vector')
92 use, intrinsic :: iso_c_binding
93 implicit none
94 integer(c_int) :: m
95 type(c_ptr), value :: msk, x, y, z, strm
97 end interface
98#endif
99
102
103contains
104
105 subroutine device_zero_dirichlet_apply_scalar(msk, x, m, strm)
106 integer, intent(in) :: m
107 type(c_ptr) :: msk, x, strm
108
109#ifdef HAVE_HIP
110 call hip_zero_dirichlet_apply_scalar(msk, x, m, strm)
111#elif HAVE_CUDA
112 call cuda_zero_dirichlet_apply_scalar(msk, x, m, strm)
113#elif HAVE_OPENCL
114 call opencl_zero_dirichlet_apply_scalar(msk, x, m, strm)
115#else
116 call neko_error('No device backend configured')
117#endif
118
120
121 subroutine device_zero_dirichlet_apply_vector(msk, x, y, z, m, strm)
122 integer, intent(in) :: m
123 type(c_ptr) :: msk, x, y, z, strm
124
125#ifdef HAVE_HIP
126 call hip_zero_dirichlet_apply_vector(msk, x, y, z, m, strm)
127#elif HAVE_CUDA
128 call cuda_zero_dirichlet_apply_vector(msk, x, y, z, m, strm)
129#elif HAVE_OPENCL
130 call opencl_zero_dirichlet_apply_vector(msk, x, y, z, m, strm)
131#else
132 call neko_error('No device backend configured')
133#endif
134
136
137end module device_zero_dirichlet
subroutine, public device_zero_dirichlet_apply_vector(msk, x, y, z, m, strm)
subroutine, public device_zero_dirichlet_apply_scalar(msk, x, m, strm)
Utilities.
Definition utils.f90:35
void opencl_zero_dirichlet_apply_vector(void *msk, void *x, void *y, void *z, int *m, cl_command_queue cmd_queue)
void opencl_zero_dirichlet_apply_scalar(void *msk, void *x, int *m, cl_command_queue cmd_queue)
void cuda_zero_dirichlet_apply_scalar(void *msk, void *x, int *m, cudaStream_t strm)
void cuda_zero_dirichlet_apply_vector(void *msk, void *x, void *y, void *z, int *m, cudaStream_t strm)