39 use iso_c_binding,
only: c_associated
48 use,
intrinsic :: iso_c_binding, only : c_ptr, c_null_ptr
49 use json_module,
only : json_file
54 type,
public,
abstract ::
bc_t
56 integer,
allocatable :: msk(:)
58 integer,
allocatable :: facet(:)
70 type(c_ptr) :: msk_d = c_null_ptr
72 type(c_ptr) :: facet_d = c_null_ptr
77 logical :: strong = .true.
119 class(
bc_t),
pointer :: ptr => null()
124 class(
bc_t),
allocatable :: obj
132 class(
bc_t),
intent(inout),
target :: this
133 type(
coef_t),
intent(in) :: coef
134 type(json_file),
intent(inout) :: json
142 class(
bc_t),
intent(inout),
target :: this
150 class(
bc_t),
intent(inout),
target :: this
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
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
206 class(
bc_t),
intent(inout),
target :: this
208 real(kind=
rp),
intent(in),
optional :: t
209 integer,
intent(in),
optional :: tstep
210 logical,
intent(in),
optional :: strong
226 class(
bc_t),
intent(inout),
target :: this
230 real(kind=
rp),
intent(in),
optional :: t
231 integer,
intent(in),
optional :: tstep
232 logical,
intent(in),
optional :: strong
241 class(
bc_t),
intent(inout) :: this
242 type(
coef_t),
target,
intent(in) :: coef
248 this%Xh => this%dof%Xh
249 this%msh => this%dof%msh
251 call this%marked_facet%init()
257 class(
bc_t),
intent(inout) :: this
259 call this%marked_facet%free()
266 if (
allocated(this%msk))
then
270 if (
allocated(this%facet))
then
271 deallocate(this%facet)
274 if (c_associated(this%msk_d))
then
276 this%msk_d = c_null_ptr
279 if (c_associated(this%facet_d))
then
281 this%facet_d = c_null_ptr
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
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)
319 call this%apply_vector_dev(x_d, y_d, z_d)
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)
329 call this%apply_vector(x, y, z, n)
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
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)
362 call this%apply_scalar_dev(x_d)
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)
372 call this%apply_scalar(x, n)
382 class(
bc_t),
intent(inout) :: this
383 integer,
intent(in) :: facet
384 integer,
intent(in) :: el
388 call this%marked_facet%push(t)
395 class(
bc_t),
intent(inout) :: this
400 fp => facet_list%array()
401 do i = 1, facet_list%size()
402 call this%marked_facet%push(fp(i))
410 class(
bc_t),
intent(inout) :: this
413 do i = 1, bc_zone%size
414 call this%marked_facet%push(bc_zone%facet_el(i))
422 class(
bc_t),
target,
intent(inout) :: this
425 integer :: facet_size, facet, el
426 integer :: i, j, k, l, msk_c
427 integer :: lx, ly, lz, n
437 allocate(this%msk(0:facet_size * this%marked_facet%size()))
438 allocate(this%facet(0:facet_size * this%marked_facet%size()))
441 bfp => this%marked_facet%array()
447 do i = 1, this%marked_facet%size()
449 facet = bc_facet%x(1)
457 this%facet(msk_c) = 1
464 this%msk(msk_c) =
linear_index(lx, k, l, el, lx, ly, lz)
465 this%facet(msk_c) = 2
473 this%facet(msk_c) = 3
480 this%msk(msk_c) =
linear_index(j, ly, l, el, lx, ly, lz)
481 this%facet(msk_c) = 4
489 this%facet(msk_c) = 5
496 this%msk(msk_c) =
linear_index(j, k, lz, el, lx, ly, lz)
497 this%facet(msk_c) = 6
504 this%facet(0) = msk_c
507 n = facet_size * this%marked_facet%size() + 1
526 class(
bc_t),
intent(inout) :: this
527 character(len=*),
intent(in) :: file_name
532 call bdry_field%init(this%coef%dof,
'bdry')
536 bdry_field%x(k,1,1,1) = 1.0_rp
539 dump_file =
file_t(file_name)
540 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.
Return the device pointer for an associated Fortran array.
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_finalize_base(this)
Finalize the construction of the bc by populting the msk and facet arrays.
subroutine bc_free_base(this)
Destructor for the base type, bc_t.
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.
subroutine bc_init_base(this, coef)
Constructor.
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.
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.
Device abstraction, common interface for various accelerators.
integer, parameter, public host_to_device
subroutine, public device_free(x_d)
Deallocate memory on the device.
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.
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.
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.