39 use iso_c_binding,
only: c_associated
51 use,
intrinsic :: iso_c_binding, only : c_ptr, c_null_ptr
52 use json_module,
only : json_file
61 type,
public,
abstract ::
bc_t
63 integer,
allocatable :: msk(:)
65 integer,
allocatable :: facet(:)
69 type(
coef_t),
pointer :: coef => null()
71 type(
mesh_t),
pointer :: msh => null()
77 type(c_ptr) :: msk_d = c_null_ptr
79 type(c_ptr) :: facet_d = c_null_ptr
84 logical :: strong = .true.
87 logical :: updated = .false.
129 class(
bc_t),
pointer :: ptr => null()
134 class(
bc_t),
allocatable :: obj
142 class(
bc_t),
intent(inout),
target :: this
143 type(
coef_t),
target,
intent(in) :: coef
144 type(json_file),
intent(inout) :: json
152 class(
bc_t),
intent(inout),
target :: this
160 class(
bc_t),
intent(inout),
target :: this
161 logical,
optional,
intent(in) :: only_facets
174 class(
bc_t),
intent(inout) :: this
175 integer,
intent(in) :: n
176 real(kind=
rp),
intent(inout),
dimension(n) :: x
178 logical,
intent(in),
optional :: strong
194 class(
bc_t),
intent(inout) :: this
195 integer,
intent(in) :: n
196 real(kind=
rp),
intent(inout),
dimension(n) :: x
197 real(kind=
rp),
intent(inout),
dimension(n) :: y
198 real(kind=
rp),
intent(inout),
dimension(n) :: z
200 logical,
intent(in),
optional :: strong
214 class(
bc_t),
intent(inout),
target :: this
215 type(c_ptr),
intent(inout) :: x_d
217 logical,
intent(in),
optional :: strong
218 type(c_ptr),
intent(inout) :: strm
233 class(
bc_t),
intent(inout),
target :: this
234 type(c_ptr),
intent(inout) :: x_d
235 type(c_ptr),
intent(inout) :: y_d
236 type(c_ptr),
intent(inout) :: z_d
238 logical,
intent(in),
optional :: strong
239 type(c_ptr),
intent(inout) :: strm
248 class(
bc_t),
intent(inout) :: this
249 type(
coef_t),
target,
intent(in) :: coef
255 this%Xh => this%dof%Xh
256 this%msh => this%dof%msh
258 call this%marked_facet%init()
264 class(
bc_t),
intent(inout) :: this
266 call this%marked_facet%free()
273 if (
allocated(this%msk))
then
277 if (
allocated(this%facet))
then
278 deallocate(this%facet)
281 if (c_associated(this%msk_d))
then
283 this%msk_d = c_null_ptr
286 if (c_associated(this%facet_d))
then
288 this%facet_d = c_null_ptr
302 class(
bc_t),
intent(inout) :: this
303 type(
field_t),
intent(inout) :: x
304 type(
field_t),
intent(inout) :: y
305 type(
field_t),
intent(inout) :: z
307 logical,
intent(in),
optional :: strong
308 type(c_ptr),
intent(inout),
optional :: strm
311 character(len=256) :: msg
317 if (y%size() .ne. n .or. z%size() .ne. n)
then
318 msg =
"Fields x, y, z must have the same size in " // &
319 "bc_list_apply_vector_field"
325 if (
present(strm))
then
331 call this%apply_vector_dev(x%x_d, y%x_d, z%x_d, time = time, &
332 strong = strong, strm = strm_)
334 call this%apply_vector(x%x, y%x, z%x, n, time = time, strong = strong)
346 class(
bc_t),
intent(inout) :: this
347 type(
field_t),
intent(inout) :: x
349 logical,
intent(in),
optional :: strong
350 type(c_ptr),
intent(inout),
optional :: strm
359 if (
present(strm))
then
365 call this%apply_scalar_dev(x%x_d, time = time, strong = strong, &
368 call this%apply_scalar(x%x, n, time = time)
377 class(
bc_t),
intent(inout) :: this
378 integer,
intent(in) :: facet
379 integer,
intent(in) :: el
383 call this%marked_facet%push(t)
390 class(
bc_t),
intent(inout) :: this
395 fp => facet_list%array()
396 do i = 1, facet_list%size()
397 call this%marked_facet%push(fp(i))
405 class(
bc_t),
intent(inout) :: this
408 do i = 1, bc_zone%size
409 call this%marked_facet%push(bc_zone%facet_el(i))
421 class(
bc_t),
target,
intent(inout) :: this
422 logical,
optional,
intent(in) :: only_facets
426 integer :: facet_size, facet, el
427 logical :: only_facet = .false.
428 integer :: i, j, k, l, msk_c
429 integer :: lx, ly, lz, n
433 if (
present(only_facets))
then
434 only_facet = only_facets
442 allocate(this%msk(0:facet_size * this%marked_facet%size()))
443 allocate(this%facet(0:facet_size * this%marked_facet%size()))
446 bfp => this%marked_facet%array()
452 do i = 1, this%marked_facet%size()
454 facet = bc_facet%x(1)
462 this%facet(msk_c) = 1
469 this%msk(msk_c) =
linear_index(lx, k, l, el, lx, ly, lz)
470 this%facet(msk_c) = 2
478 this%facet(msk_c) = 3
485 this%msk(msk_c) =
linear_index(j, ly, l, el, lx, ly, lz)
486 this%facet(msk_c) = 4
494 this%facet(msk_c) = 5
501 this%msk(msk_c) =
linear_index(j, k, lz, el, lx, ly, lz)
502 this%facet(msk_c) = 6
507 this%facet(0) = msk_c
515 if ( .not. only_facet)
then
517 call test_field%init(this%dof)
519 n = test_field%size()
520 test_field%x = 0.0_rp
523 test_field%x(this%msk(i),1,1,1) = 1.0
530 call this%coef%gs_h%op(test_field,
gs_op_add)
536 do i = 1, this%dof%size()
537 if (test_field%x(i,1,1,1) .gt. 0.5)
then
543 allocate(this%msk(0:msk_c))
545 do i = 1, this%dof%size()
546 if (test_field%x(i,1,1,1) .gt. 0.5)
then
552 call test_field%free()
570 class(
bc_t),
intent(inout) :: this
571 character(len=*),
intent(in) :: file_name
576 call bdry_field%init(this%coef%dof,
'bdry')
580 bdry_field%x(k,1,1,1) = 1.0_rp
582 call dump_file%init(file_name)
583 call dump_file%write(bdry_field)
Apply the boundary condition to a scalar field on the device.
Apply the boundary condition to a scalar field.
Apply the boundary condition to a vector field on the device.
Apply the boundary condition to a vector field.
Finalize by building the mask and facet arrays.
Map a Fortran array to a device (allocate and associate)
Copy data between host and device (or device and device)
Defines a boundary condition.
subroutine bc_mark_zone(this, bc_zone)
Mark all facets from a zone.
subroutine bc_free_base(this)
Destructor for the base type, bc_t.
subroutine bc_init_base(this, coef)
Constructor.
subroutine bc_apply_scalar_generic(this, x, time, strong, strm)
Apply the boundary condition to a scalar field. Dispatches to the CPU or the device version.
subroutine bc_finalize_base(this, only_facets)
Finalize the construction of the bc by populting the msk and facet arrays.
subroutine bc_apply_vector_generic(this, x, y, z, time, strong, strm)
Apply the boundary condition to a vector field. Dispatches to the CPU or the device version.
subroutine bc_debug_mask(this, file_name)
Write a field showing the mask of the bc.
subroutine bc_mark_facet(this, facet, el)
Mark facet on element el as part of the boundary condition.
subroutine bc_mark_facets(this, facet_list)
Mark all facets from a (facet, el) tuple list.
Device abstraction, common interface for various accelerators.
integer, parameter, public host_to_device
subroutine, public device_free(x_d)
Deallocate memory on the device.
integer, parameter, public device_to_host
type(c_ptr), bind(C), public glb_cmd_queue
Global command queue.
Defines a mapping of the degrees of freedom.
Defines a zone as a subset of facets in a mesh.
Module for file I/O operations.
Defines Gather-scatter operations.
integer, parameter, public gs_op_add
integer, parameter, public neko_msh_max_zlbls
Max num. zone labels.
integer, parameter, public neko_msh_max_zlbl_len
Max length of a zone label.
integer, parameter neko_bcknd_device
integer, parameter, public rp
Global precision used in computations.
Defines a function space.
Implements a dynamic stack ADT.
Module with things related to the simulation time.
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,...
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, :)
Base type for a boundary condition.
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
A wrapper around a polymorphic generic_file_t that handles its init. This is essentially a factory fo...
The function space for the SEM solution fields.
Integer 2-tuple based stack.
A struct that contains all info about the time, expand as needed.