Neko 0.9.99
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
bc.f90
Go to the documentation of this file.
1! Copyright (c) 2020-2024, 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!
34module bc
36 use num_types, only : rp
39 use iso_c_binding, only: c_associated
40 use dofmap, only : dofmap_t
41 use coefs, only : coef_t
42 use space, only : space_t
44 use facet_zone, only : facet_zone_t
45 use stack, only : stack_i4t2_t
46 use tuple, only : tuple_i4_t
48 use, intrinsic :: iso_c_binding, only : c_ptr, c_null_ptr
49 use json_module, only : json_file
50 implicit none
51 private
52
54 type, public, abstract :: bc_t
56 integer, allocatable :: msk(:)
58 integer, allocatable :: facet(:)
60 type(dofmap_t), pointer :: dof
62 type(coef_t), pointer :: coef
64 type(mesh_t), pointer :: msh
66 type(space_t), pointer :: xh
68 type(stack_i4t2_t) :: marked_facet
70 type(c_ptr) :: msk_d = c_null_ptr
72 type(c_ptr) :: facet_d = c_null_ptr
77 logical :: strong = .true.
78 contains
80 procedure, pass(this) :: init_base => bc_init_base
82 procedure, pass(this) :: free_base => bc_free_base
84 procedure, pass(this) :: mark_facet => bc_mark_facet
86 procedure, pass(this) :: mark_facets => bc_mark_facets
88 procedure, pass(this) :: mark_zone => bc_mark_zone
91 procedure, pass(this) :: finalize_base => bc_finalize_base
92
95 procedure, pass(this) :: apply_scalar_generic => bc_apply_scalar_generic
98 procedure, pass(this) :: apply_vector_generic => bc_apply_vector_generic
100 procedure, pass(this) :: debug_mask_ => bc_debug_mask
102 procedure(bc_apply_scalar), pass(this), deferred :: apply_scalar
104 procedure(bc_apply_vector), pass(this), deferred :: apply_vector
106 procedure(bc_apply_scalar_dev), pass(this), deferred :: apply_scalar_dev
108 procedure(bc_apply_vector_dev), pass(this), deferred :: apply_vector_dev
110 procedure(bc_destructor), pass(this), deferred :: free
112 procedure(bc_constructor), pass(this), deferred :: init
114 procedure(bc_finalize), pass(this), deferred :: finalize
115 end type bc_t
116
118 type, public :: bc_ptr_t
119 class(bc_t), pointer :: ptr => null()
120 end type bc_ptr_t
121
122 ! Helper type to have an array of polymorphic bc_t objects.
123 type, public :: bc_alloc_t
124 class(bc_t), allocatable :: obj
125 end type
126
127
128 abstract interface
129
130 subroutine bc_constructor(this, coef, json)
131 import :: bc_t, coef_t, json_file
132 class(bc_t), intent(inout), target :: this
133 type(coef_t), intent(in) :: coef
134 type(json_file), intent(inout) :: json
135 end subroutine bc_constructor
136 end interface
137
138 abstract interface
139
140 subroutine bc_destructor(this)
141 import :: bc_t
142 class(bc_t), intent(inout), target :: this
143 end subroutine bc_destructor
144 end interface
145
146 abstract interface
147
148 subroutine bc_finalize(this)
149 import :: bc_t
150 class(bc_t), intent(inout), target :: this
151 end subroutine bc_finalize
152 end interface
153
154 abstract interface
155
161 subroutine bc_apply_scalar(this, x, n, t, tstep, strong)
162 import :: bc_t
163 import :: rp
164 class(bc_t), intent(inout) :: this
165 integer, intent(in) :: n
166 real(kind=rp), intent(inout), dimension(n) :: x
167 real(kind=rp), intent(in), optional :: t
168 integer, intent(in), optional :: tstep
169 logical, intent(in), optional :: strong
170 end subroutine bc_apply_scalar
171 end interface
172
173 abstract interface
174
182 subroutine bc_apply_vector(this, x, y, z, n, t, tstep, strong)
183 import :: bc_t
184 import :: rp
185 class(bc_t), intent(inout) :: this
186 integer, intent(in) :: n
187 real(kind=rp), intent(inout), dimension(n) :: x
188 real(kind=rp), intent(inout), dimension(n) :: y
189 real(kind=rp), intent(inout), dimension(n) :: z
190 real(kind=rp), intent(in), optional :: t
191 integer, intent(in), optional :: tstep
192 logical, intent(in), optional :: strong
193 end subroutine bc_apply_vector
194 end interface
195
196 abstract interface
197
202 subroutine bc_apply_scalar_dev(this, x_d, t, tstep, strong)
203 import :: c_ptr
204 import :: bc_t
205 import :: rp
206 class(bc_t), intent(inout), target :: this
207 type(c_ptr) :: x_d
208 real(kind=rp), intent(in), optional :: t
209 integer, intent(in), optional :: tstep
210 logical, intent(in), optional :: strong
211 end subroutine bc_apply_scalar_dev
212 end interface
213
214 abstract interface
215
222 subroutine bc_apply_vector_dev(this, x_d, y_d, z_d, t, tstep, strong)
223 import :: c_ptr
224 import :: bc_t
225 import :: rp
226 class(bc_t), intent(inout), target :: this
227 type(c_ptr) :: x_d
228 type(c_ptr) :: y_d
229 type(c_ptr) :: z_d
230 real(kind=rp), intent(in), optional :: t
231 integer, intent(in), optional :: tstep
232 logical, intent(in), optional :: strong
233 end subroutine bc_apply_vector_dev
234 end interface
235
236contains
237
240 subroutine bc_init_base(this, coef)
241 class(bc_t), intent(inout) :: this
242 type(coef_t), target, intent(in) :: coef
243
244 call this%free_base
245
246 this%dof => coef%dof
247 this%coef => coef
248 this%Xh => this%dof%Xh
249 this%msh => this%dof%msh
250
251 call this%marked_facet%init()
252
253 end subroutine bc_init_base
254
256 subroutine bc_free_base(this)
257 class(bc_t), intent(inout) :: this
258
259 call this%marked_facet%free()
260
261 nullify(this%Xh)
262 nullify(this%msh)
263 nullify(this%dof)
264 nullify(this%coef)
265
266 if (allocated(this%msk)) then
267 deallocate(this%msk)
268 end if
269
270 if (allocated(this%facet)) then
271 deallocate(this%facet)
272 end if
273
274 if (c_associated(this%msk_d)) then
275 call device_free(this%msk_d)
276 this%msk_d = c_null_ptr
277 end if
278
279 if (c_associated(this%facet_d)) then
280 call device_free(this%facet_d)
281 this%facet_d = c_null_ptr
282 end if
283
284 end subroutine bc_free_base
285
294 subroutine bc_apply_vector_generic(this, x, y, z, n, t, tstep)
295 class(bc_t), intent(inout) :: this
296 integer, intent(in) :: n
297 real(kind=rp), intent(inout), dimension(n) :: x
298 real(kind=rp), intent(inout), dimension(n) :: y
299 real(kind=rp), intent(inout), dimension(n) :: z
300 real(kind=rp), intent(in), optional :: t
301 integer, intent(in), optional :: tstep
302 type(c_ptr) :: x_d
303 type(c_ptr) :: y_d
304 type(c_ptr) :: z_d
305 integer :: i
306
307
308 if (neko_bcknd_device .eq. 1) then
309 x_d = device_get_ptr(x)
310 y_d = device_get_ptr(y)
311 z_d = device_get_ptr(z)
312 if (present(t) .and. present(tstep)) then
313 call this%apply_vector_dev(x_d, y_d, z_d, t = t, tstep = tstep)
314 else if (present(t)) then
315 call this%apply_vector_dev(x_d, y_d, z_d, t = t)
316 else if (present(tstep)) then
317 call this%apply_vector_dev(x_d, y_d, z_d, tstep = tstep)
318 else
319 call this%apply_vector_dev(x_d, y_d, z_d)
320 end if
321 else
322 if (present(t) .and. present(tstep)) then
323 call this%apply_vector(x, y, z, n, t = t, tstep = tstep)
324 else if (present(t)) then
325 call this%apply_vector(x, y, z, n, t = t)
326 else if (present(tstep)) then
327 call this%apply_vector(x, y, z, n, tstep = tstep)
328 else
329 call this%apply_vector(x, y, z, n)
330 end if
331 end if
332
333 end subroutine bc_apply_vector_generic
334
343 subroutine bc_apply_scalar_generic(this, x, n, t, tstep)
344 class(bc_t), intent(inout) :: this
345 integer, intent(in) :: n
346 real(kind=rp), intent(inout), dimension(n) :: x
347 real(kind=rp), intent(in), optional :: t
348 integer, intent(in), optional :: tstep
349 type(c_ptr) :: x_d
350 integer :: i
351
352
353 if (neko_bcknd_device .eq. 1) then
354 x_d = device_get_ptr(x)
355 if (present(t) .and. present(tstep)) then
356 call this%apply_scalar_dev(x_d, t = t, tstep = tstep)
357 else if (present(t)) then
358 call this%apply_scalar_dev(x_d, t = t)
359 else if (present(tstep)) then
360 call this%apply_scalar_dev(x_d, tstep = tstep)
361 else
362 call this%apply_scalar_dev(x_d)
363 end if
364 else
365 if (present(t) .and. present(tstep)) then
366 call this%apply_scalar(x, n, t = t, tstep = tstep)
367 else if (present(t)) then
368 call this%apply_scalar(x, n, t = t)
369 else if (present(tstep)) then
370 call this%apply_scalar(x, n, tstep = tstep)
371 else
372 call this%apply_scalar(x, n)
373 end if
374 end if
375
376 end subroutine bc_apply_scalar_generic
377
381 subroutine bc_mark_facet(this, facet, el)
382 class(bc_t), intent(inout) :: this
383 integer, intent(in) :: facet
384 integer, intent(in) :: el
385 type(tuple_i4_t) :: t
386
387 t%x = [facet, el]
388 call this%marked_facet%push(t)
389
390 end subroutine bc_mark_facet
391
394 subroutine bc_mark_facets(this, facet_list)
395 class(bc_t), intent(inout) :: this
396 type(stack_i4t2_t), intent(inout) :: facet_list
397 type(tuple_i4_t), pointer :: fp(:)
398 integer :: i
399
400 fp => facet_list%array()
401 do i = 1, facet_list%size()
402 call this%marked_facet%push(fp(i))
403 end do
404
405 end subroutine bc_mark_facets
406
409 subroutine bc_mark_zone(this, bc_zone)
410 class(bc_t), intent(inout) :: this
411 class(facet_zone_t), intent(in) :: bc_zone
412 integer :: i
413 do i = 1, bc_zone%size
414 call this%marked_facet%push(bc_zone%facet_el(i))
415 end do
416 end subroutine bc_mark_zone
417
421 subroutine bc_finalize_base(this)
422 class(bc_t), target, intent(inout) :: this
423 type(tuple_i4_t), pointer :: bfp(:)
424 type(tuple_i4_t) :: bc_facet
425 integer :: facet_size, facet, el
426 integer :: i, j, k, l, msk_c
427 integer :: lx, ly, lz, n
428
429 lx = this%Xh%lx
430 ly = this%Xh%ly
431 lz = this%Xh%lz
432
434
435 ! Note we assume that lx = ly = lz
436 facet_size = lx**2
437 allocate(this%msk(0:facet_size * this%marked_facet%size()))
438 allocate(this%facet(0:facet_size * this%marked_facet%size()))
439
440 msk_c = 0
441 bfp => this%marked_facet%array()
442
443 ! Loop through each (facet, element) id tuple
444 ! Then loop over all the nodes of the face and compute their linear index
445 ! This index goes into this%msk, whereas the corresponding face id goes into
446 ! this%facet
447 do i = 1, this%marked_facet%size()
448 bc_facet = bfp(i)
449 facet = bc_facet%x(1)
450 el = bc_facet%x(2)
451 select case (facet)
452 case (1)
453 do l = 1, lz
454 do k = 1, ly
455 msk_c = msk_c + 1
456 this%msk(msk_c) = linear_index(1, k, l, el, lx, ly, lz)
457 this%facet(msk_c) = 1
458 end do
459 end do
460 case (2)
461 do l = 1, lz
462 do k = 1, ly
463 msk_c = msk_c + 1
464 this%msk(msk_c) = linear_index(lx, k, l, el, lx, ly, lz)
465 this%facet(msk_c) = 2
466 end do
467 end do
468 case (3)
469 do l = 1, lz
470 do j = 1, lx
471 msk_c = msk_c + 1
472 this%msk(msk_c) = linear_index(j, 1, l, el, lx, ly, lz)
473 this%facet(msk_c) = 3
474 end do
475 end do
476 case (4)
477 do l = 1, lz
478 do j = 1, lx
479 msk_c = msk_c + 1
480 this%msk(msk_c) = linear_index(j, ly, l, el, lx, ly, lz)
481 this%facet(msk_c) = 4
482 end do
483 end do
484 case (5)
485 do k = 1, ly
486 do j = 1, lx
487 msk_c = msk_c + 1
488 this%msk(msk_c) = linear_index(j, k, 1, el, lx, ly, lz)
489 this%facet(msk_c) = 5
490 end do
491 end do
492 case (6)
493 do k = 1, ly
494 do j = 1, lx
495 msk_c = msk_c + 1
496 this%msk(msk_c) = linear_index(j, k, lz, el, lx, ly, lz)
497 this%facet(msk_c) = 6
498 end do
499 end do
500 end select
501 end do
502
503 this%msk(0) = msk_c
504 this%facet(0) = msk_c
505
506 if (neko_bcknd_device .eq. 1) then
507 n = facet_size * this%marked_facet%size() + 1
508 call device_map(this%msk, this%msk_d, n)
509 call device_map(this%facet, this%facet_d, n)
510
511 call device_memcpy(this%msk, this%msk_d, n, &
512 host_to_device, sync = .false.)
513 call device_memcpy(this%facet, this%facet_d, n, &
514 host_to_device, sync = .false.)
515 end if
516
517 end subroutine bc_finalize_base
518
522 subroutine bc_debug_mask(this, file_name)
523 use field, only : field_t
524 use file, only : file_t
525
526 class(bc_t), intent(inout) :: this
527 character(len=*), intent(in) :: file_name
528 type(field_t) :: bdry_field
529 integer:: i, m, k
530 type(file_t) :: dump_file
531
532 call bdry_field%init(this%coef%dof, 'bdry')
533 m = this%msk(0)
534 do i = 1, m
535 k = this%msk(i)
536 bdry_field%x(k,1,1,1) = 1.0_rp
537 end do
538
539 dump_file = file_t(file_name)
540 call dump_file%write(bdry_field)
541
542
543 end subroutine bc_debug_mask
544end module bc
Apply the boundary condition to a scalar field on the device.
Definition bc.f90:202
Apply the boundary condition to a scalar field.
Definition bc.f90:161
Apply the boundary condition to a vector field on the device.
Definition bc.f90:222
Apply the boundary condition to a vector field.
Definition bc.f90:182
Constructor.
Definition bc.f90:130
Destructor.
Definition bc.f90:140
Finalize by building the mask and facet arrays.
Definition bc.f90:148
Return the device pointer for an associated Fortran array.
Definition device.F90:92
Map a Fortran array to a device (allocate and associate)
Definition device.F90:68
Copy data between host and device (or device and device)
Definition device.F90:62
Defines a boundary condition.
Definition bc.f90:34
subroutine bc_mark_zone(this, bc_zone)
Mark all facets from a zone.
Definition bc.f90:410
subroutine bc_finalize_base(this)
Finalize the construction of the bc by populting the msk and facet arrays.
Definition bc.f90:422
subroutine bc_free_base(this)
Destructor for the base type, bc_t.
Definition bc.f90:257
subroutine bc_apply_scalar_generic(this, x, n, t, tstep)
Apply the boundary condition to a scalar field. Dispatches to the CPU or the device version.
Definition bc.f90:344
subroutine bc_init_base(this, coef)
Constructor.
Definition bc.f90:241
subroutine bc_debug_mask(this, file_name)
Write a field showing the mask of the bc.
Definition bc.f90:523
subroutine bc_mark_facet(this, facet, el)
Mark facet on element el as part of the boundary condition.
Definition bc.f90:382
subroutine bc_mark_facets(this, facet_list)
Mark all facets from a (facet, el) tuple list.
Definition bc.f90:395
subroutine bc_apply_vector_generic(this, x, y, z, n, t, tstep)
Apply the boundary condition to a vector field. Dispatches to the CPU or the device version.
Definition bc.f90:295
Coefficients.
Definition coef.f90:34
Device abstraction, common interface for various accelerators.
Definition device.F90:34
integer, parameter, public host_to_device
Definition device.F90:46
subroutine, public device_free(x_d)
Deallocate memory on the device.
Definition device.F90:197
Defines a mapping of the degrees of freedom.
Definition dofmap.f90:35
Defines a zone as a subset of facets in a mesh.
Defines a field.
Definition field.f90:34
Module for file I/O operations.
Definition file.f90:34
Defines a mesh.
Definition mesh.f90:34
integer, parameter, public neko_msh_max_zlbls
Max num. zone labels.
Definition mesh.f90:56
integer, parameter, public neko_msh_max_zlbl_len
Max length of a zone label.
Definition mesh.f90:58
Build configurations.
integer, parameter neko_bcknd_device
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
Defines a function space.
Definition space.f90:34
Implements a dynamic stack ADT.
Definition stack.f90:35
Implements a n-tuple.
Definition tuple.f90:34
Utilities.
Definition utils.f90:35
character(len=100) function, dimension(:), allocatable, public split_string(string, delimiter)
Split a string based on delimiter (tokenizer) OBS: very hacky, this should really be improved,...
Definition utils.f90:136
pure integer function, public linear_index(i, j, k, l, lx, ly, lz)
Compute the address of a (i,j,k,l) array with sizes (1:lx, 1:ly, 1:lz, :)
Definition utils.f90:174
Pointer to a `bc_t`.
Definition bc.f90:118
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
A wrapper around a polymorphic generic_file_t that handles its init. This is essentially a factory fo...
Definition file.f90:54
The function space for the SEM solution fields.
Definition space.f90:62
Integer 2-tuple based stack.
Definition stack.f90:84
Integer based 2-tuple.
Definition tuple.f90:51