Neko 0.9.99
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
usr_inflow.f90
Go to the documentation of this file.
1! Copyright (c) 2020-2023, 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!
35 use num_types, only : rp
36 use coefs, only : coef_t
37 use inflow, only : inflow_t
38 use device
41 use bc, only : bc_t
42 use json_module, only : json_file
43 implicit none
44 private
45
47 type, public, extends(bc_t) :: usr_inflow_t
48 procedure(usr_inflow_eval), nopass, pointer :: eval => null()
49 type(c_ptr), private :: usr_x_d = c_null_ptr
50 type(c_ptr), private :: usr_y_d = c_null_ptr
51 type(c_ptr), private :: usr_z_d = c_null_ptr
52 contains
53 procedure, pass(this) :: apply_scalar => usr_inflow_apply_scalar
54 procedure, pass(this) :: apply_vector => usr_inflow_apply_vector
55 procedure, pass(this) :: validate => usr_inflow_validate
56 procedure, pass(this) :: set_eval => usr_inflow_set_eval
57 procedure, pass(this) :: apply_vector_dev => usr_inflow_apply_vector_dev
58 procedure, pass(this) :: apply_scalar_dev => usr_inflow_apply_scalar_dev
60 procedure, pass(this) :: init => usr_inflow_init
62 procedure, pass(this) :: free => usr_inflow_free
64 procedure, pass(this) :: finalize => usr_inflow_finalize
65 end type usr_inflow_t
66
67 abstract interface
68
85 subroutine usr_inflow_eval(u, v, w, x, y, z, nx, ny, nz, &
86 ix, iy, iz, ie, t, tstep)
87 import rp
88 real(kind=rp), intent(inout) :: u
89 real(kind=rp), intent(inout) :: v
90 real(kind=rp), intent(inout) :: w
91 real(kind=rp), intent(in) :: x
92 real(kind=rp), intent(in) :: y
93 real(kind=rp), intent(in) :: z
94 real(kind=rp), intent(in) :: nx
95 real(kind=rp), intent(in) :: ny
96 real(kind=rp), intent(in) :: nz
97 integer, intent(in) :: ix
98 integer, intent(in) :: iy
99 integer, intent(in) :: iz
100 integer, intent(in) :: ie
101 real(kind=rp), intent(in) :: t
102 integer, intent(in) :: tstep
103 end subroutine usr_inflow_eval
104 end interface
105
106 public :: usr_inflow_eval
107
108contains
109
113 subroutine usr_inflow_init(this, coef, json)
114 class(usr_inflow_t), intent(inout), target :: this
115 type(coef_t), intent(in) :: coef
116 type(json_file), intent(inout) ::json
117
118 call this%init_base(coef)
119 end subroutine usr_inflow_init
120
121 subroutine usr_inflow_free(this)
122 class(usr_inflow_t), target, intent(inout) :: this
123
124 call this%free_base()
125
126 if (c_associated(this%usr_x_d)) then
127 call device_free(this%usr_x_d)
128 end if
129
130 if (c_associated(this%usr_y_d)) then
131 call device_free(this%usr_y_d)
132 end if
133
134 if (c_associated(this%usr_z_d)) then
135 call device_free(this%usr_z_d)
136 end if
137
138 end subroutine usr_inflow_free
139
141 subroutine usr_inflow_apply_scalar(this, x, n, t, tstep, strong)
142 class(usr_inflow_t), intent(inout) :: this
143 integer, intent(in) :: n
144 real(kind=rp), intent(inout), dimension(n) :: x
145 real(kind=rp), intent(in), optional :: t
146 integer, intent(in), optional :: tstep
147 logical, intent(in), optional :: strong
148 end subroutine usr_inflow_apply_scalar
149
151 subroutine usr_inflow_apply_scalar_dev(this, x_d, t, tstep, strong)
152 class(usr_inflow_t), intent(inout), target :: this
153 type(c_ptr) :: x_d
154 real(kind=rp), intent(in), optional :: t
155 integer, intent(in), optional :: tstep
156 logical, intent(in), optional :: strong
157 end subroutine usr_inflow_apply_scalar_dev
158
160 subroutine usr_inflow_apply_vector(this, x, y, z, n, t, tstep, strong)
161 class(usr_inflow_t), intent(inout) :: this
162 integer, intent(in) :: n
163 real(kind=rp), intent(inout), dimension(n) :: x
164 real(kind=rp), intent(inout), dimension(n) :: y
165 real(kind=rp), intent(inout), dimension(n) :: z
166 real(kind=rp), intent(in), optional :: t
167 integer, intent(in), optional :: tstep
168 logical, intent(in), optional :: strong
169 integer :: i, m, k, idx(4), facet, tstep_
170 real(kind=rp) :: t_
171 logical :: strong_ = .true.
172
173 if (present(strong)) strong_ = strong
174
175 if (present(t)) then
176 t_ = t
177 else
178 t_ = 0.0_rp
179 end if
180
181 if (present(tstep)) then
182 tstep_ = tstep
183 else
184 tstep_ = 1
185 end if
186
187 associate(&
188 xc => this%coef%dof%x, yc => this%coef%dof%y, zc => this%coef%dof%z, &
189 nx => this%coef%nx, ny => this%coef%ny, nz => this%coef%nz, &
190 lx => this%coef%Xh%lx)
191 m = this%msk(0)
192 if (strong_) then
193 do i = 1, m
194 k = this%msk(i)
195 facet = this%facet(i)
196 idx = nonlinear_index(k, lx, lx, lx)
197 select case (facet)
198 case (1,2)
199 call this%eval(x(k), y(k), z(k), &
200 xc(idx(1), idx(2), idx(3), idx(4)), &
201 yc(idx(1), idx(2), idx(3), idx(4)), &
202 zc(idx(1), idx(2), idx(3), idx(4)), &
203 nx(idx(2), idx(3), facet, idx(4)), &
204 ny(idx(2), idx(3), facet, idx(4)), &
205 nz(idx(2), idx(3), facet, idx(4)), &
206 idx(1), idx(2), idx(3), idx(4), &
207 t_, tstep_)
208 case (3,4)
209 call this%eval(x(k), y(k), z(k), &
210 xc(idx(1), idx(2), idx(3), idx(4)), &
211 yc(idx(1), idx(2), idx(3), idx(4)), &
212 zc(idx(1), idx(2), idx(3), idx(4)), &
213 nx(idx(1), idx(3), facet, idx(4)), &
214 ny(idx(1), idx(3), facet, idx(4)), &
215 nz(idx(1), idx(3), facet, idx(4)), &
216 idx(1), idx(2), idx(3), idx(4), &
217 t_, tstep_)
218 case (5,6)
219 call this%eval(x(k), y(k), z(k), &
220 xc(idx(1), idx(2), idx(3), idx(4)), &
221 yc(idx(1), idx(2), idx(3), idx(4)), &
222 zc(idx(1), idx(2), idx(3), idx(4)), &
223 nx(idx(1), idx(2), facet, idx(4)), &
224 ny(idx(1), idx(2), facet, idx(4)), &
225 nz(idx(1), idx(2), facet, idx(4)), &
226 idx(1), idx(2), idx(3), idx(4), &
227 t_, tstep_)
228 end select
229 end do
230 end if
231 end associate
232
233 end subroutine usr_inflow_apply_vector
234
235 subroutine usr_inflow_apply_vector_dev(this, x_d, y_d, z_d, t, tstep, strong)
236 class(usr_inflow_t), intent(inout), target :: this
237 type(c_ptr) :: x_d
238 type(c_ptr) :: y_d
239 type(c_ptr) :: z_d
240 real(kind=rp), intent(in), optional :: t
241 integer, intent(in), optional :: tstep
242 logical, intent(in), optional :: strong
243 integer :: i, m, k, idx(4), facet, tstep_
244 integer(c_size_t) :: s
245 real(kind=rp) :: t_
246 real(kind=rp), allocatable :: x(:)
247 real(kind=rp), allocatable :: y(:)
248 real(kind=rp), allocatable :: z(:)
249 logical :: strong_ = .true.
250
251 if (present(strong)) strong_ = strong
252
253 if (present(t)) then
254 t_ = t
255 else
256 t_ = 0.0_rp
257 end if
258
259 if (present(tstep)) then
260 tstep_ = tstep
261 else
262 tstep_ = 1
263 end if
264
265 associate(&
266 xc => this%coef%dof%x, yc => this%coef%dof%y, zc => this%coef%dof%z, &
267 nx => this%coef%nx, ny => this%coef%ny, nz => this%coef%nz, &
268 lx => this%coef%Xh%lx, usr_x_d => this%usr_x_d, &
269 usr_y_d => this%usr_y_d, usr_z_d => this%usr_z_d)
270
271 m = this%msk(0)
272
273
274 ! Pretabulate values during first call to apply
275 if (.not. c_associated(usr_x_d) .and. strong_ .and. &
276 (this%msk(0) .gt. 0)) then
277 allocate(x(m), y(m), z(m)) ! Temp arrays
278
279 s = m*rp
280
281 call device_alloc(usr_x_d, s)
282 call device_alloc(usr_y_d, s)
283 call device_alloc(usr_z_d, s)
284
285 associate(xc => this%coef%dof%x, yc => this%coef%dof%y, &
286 zc => this%coef%dof%z, &
287 nx => this%coef%nx, ny => this%coef%ny, nz => this%coef%nz, &
288 lx => this%coef%Xh%lx)
289 do i = 1, m
290 k = this%msk(i)
291 facet = this%facet(i)
292 idx = nonlinear_index(k, lx, lx, lx)
293 select case (facet)
294 case (1,2)
295 call this%eval(x(i), y(i), z(i), &
296 xc(idx(1), idx(2), idx(3), idx(4)), &
297 yc(idx(1), idx(2), idx(3), idx(4)), &
298 zc(idx(1), idx(2), idx(3), idx(4)), &
299 nx(idx(2), idx(3), facet, idx(4)), &
300 ny(idx(2), idx(3), facet, idx(4)), &
301 nz(idx(2), idx(3), facet, idx(4)), &
302 idx(1), idx(2), idx(3), idx(4), &
303 t_, tstep_)
304 case (3,4)
305 call this%eval(x(i), y(i), z(i), &
306 xc(idx(1), idx(2), idx(3), idx(4)), &
307 yc(idx(1), idx(2), idx(3), idx(4)), &
308 zc(idx(1), idx(2), idx(3), idx(4)), &
309 nx(idx(1), idx(3), facet, idx(4)), &
310 ny(idx(1), idx(3), facet, idx(4)), &
311 nz(idx(1), idx(3), facet, idx(4)), &
312 idx(1), idx(2), idx(3), idx(4), &
313 t_, tstep_)
314 case (5,6)
315 call this%eval(x(i), y(i), z(i), &
316 xc(idx(1), idx(2), idx(3), idx(4)), &
317 yc(idx(1), idx(2), idx(3), idx(4)), &
318 zc(idx(1), idx(2), idx(3), idx(4)), &
319 nx(idx(1), idx(2), facet, idx(4)), &
320 ny(idx(1), idx(2), facet, idx(4)), &
321 nz(idx(1), idx(2), facet, idx(4)), &
322 idx(1), idx(2), idx(3), idx(4), &
323 t_, tstep_)
324 end select
325 end do
326 end associate
327
328 call device_memcpy(x, usr_x_d, m, host_to_device, sync = .false.)
329 call device_memcpy(y, usr_y_d, m, host_to_device, sync = .false.)
330 call device_memcpy(z, usr_z_d, m, host_to_device, sync = .true.)
331
332 deallocate(x, y, z)
333 end if
334
335 if (strong_ .and. (this%msk(0) .gt. 0)) then
336 call device_inhom_dirichlet_apply_vector(this%msk_d, x_d, y_d, z_d, &
337 usr_x_d, usr_y_d, usr_z_d, m)
338 end if
339
340 end associate
341
342 end subroutine usr_inflow_apply_vector_dev
343
346 subroutine usr_inflow_set_eval(this, usr_eval)
347 class(usr_inflow_t), intent(inout) :: this
348 procedure(usr_inflow_eval) :: usr_eval
349 this%eval => usr_eval
350 end subroutine usr_inflow_set_eval
351
353 subroutine usr_inflow_validate(this)
354 class(usr_inflow_t), intent(inout) :: this
355 logical :: valid
356
357 valid = .true. ! Assert it's going to be ok...
358 if (.not. associated(this%coef)) then
359 call neko_warning('Missing coefficients')
360 valid = .false.
361 end if
362
363 if (.not. associated(this%eval)) then
364 call neko_warning('Missing eval function')
365 valid = .false.
366 end if
367
368 if (.not. valid) then
369 call neko_error('Invalid user defined inflow condition')
370 end if
371
372 end subroutine usr_inflow_validate
373
375 subroutine usr_inflow_finalize(this)
376 class(usr_inflow_t), target, intent(inout) :: this
377
378 call this%finalize_base()
379 end subroutine usr_inflow_finalize
380end module usr_inflow
__device__ void nonlinear_index(const int idx, const int lx, int *index)
Abstract interface defining a user defined inflow condition (pointwise)
Defines a boundary condition.
Definition bc.f90:34
Coefficients.
Definition coef.f90:34
Device abstraction, common interface for various accelerators.
Definition device.F90:34
subroutine, public device_free(x_d)
Deallocate memory on the device.
Definition device.F90:197
Defines inflow dirichlet conditions.
Definition inflow.f90:34
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
Defines inflow dirichlet conditions.
subroutine usr_inflow_apply_scalar(this, x, n, t, tstep, strong)
No-op scalar apply.
subroutine usr_inflow_finalize(this)
Finalize.
subroutine usr_inflow_apply_vector_dev(this, x_d, y_d, z_d, t, tstep, strong)
subroutine usr_inflow_apply_scalar_dev(this, x_d, t, tstep, strong)
No-op scalar apply (device version)
subroutine usr_inflow_apply_vector(this, x, y, z, n, t, tstep, strong)
Apply user defined inflow conditions (vector valued)
subroutine usr_inflow_init(this, coef, json)
Constructor.
subroutine usr_inflow_validate(this)
Validate user inflow condition.
subroutine usr_inflow_free(this)
subroutine usr_inflow_set_eval(this, usr_eval)
Assign user provided eval function.
Utilities.
Definition utils.f90:35
subroutine, public neko_warning(warning_msg)
Reports a warning to standard output.
Definition utils.f90:266
Base type for a boundary condition.
Definition bc.f90:54
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
Definition coef.f90:55
Dirichlet condition for inlet (vector valued)
Definition inflow.f90:46
User defined dirichlet condition for velocity.