46 use,
intrinsic :: iso_c_binding, only : c_ptr, c_null_ptr
51 type,
public,
abstract ::
bc_t
53 integer,
allocatable :: msk(:)
55 integer,
allocatable :: facet(:)
67 type(c_ptr) :: msk_d = c_null_ptr
69 type(c_ptr) :: facet_d = c_null_ptr
98 class(
bc_t),
pointer :: bcp
119 class(
bc_t),
intent(inout) :: this
120 integer,
intent(in) :: n
121 real(kind=
rp),
intent(inout),
dimension(n) :: x
122 real(kind=
rp),
intent(in),
optional :: t
123 integer,
intent(in),
optional :: tstep
138 class(
bc_t),
intent(inout) :: this
139 integer,
intent(in) :: n
140 real(kind=
rp),
intent(inout),
dimension(n) :: x
141 real(kind=
rp),
intent(inout),
dimension(n) :: y
142 real(kind=
rp),
intent(inout),
dimension(n) :: z
143 real(kind=
rp),
intent(in),
optional :: t
144 integer,
intent(in),
optional :: tstep
155 class(
bc_t),
intent(inout),
target :: this
157 real(kind=
rp),
intent(in),
optional :: t
158 integer,
intent(in),
optional :: tstep
171 class(
bc_t),
intent(inout),
target :: this
175 real(kind=
rp),
intent(in),
optional :: t
176 integer,
intent(in),
optional :: tstep
192 class(
bc_t),
intent(inout) :: this
193 type(
coef_t),
target,
intent(in) :: coef
199 this%Xh => this%dof%Xh
200 this%msh => this%dof%msh
202 call this%marked_facet%init()
208 class(
bc_t),
intent(inout) :: this
210 call this%marked_facet%free()
216 if (
allocated(this%msk))
then
220 if (
allocated(this%facet))
then
221 deallocate(this%facet)
224 if (c_associated(this%msk_d))
then
226 this%msk_d = c_null_ptr
229 if (c_associated(this%facet_d))
then
231 this%facet_d = c_null_ptr
240 class(
bc_t),
intent(inout) :: this
241 integer,
intent(in) :: facet
242 integer,
intent(in) :: el
246 call this%marked_facet%push(t)
253 class(
bc_t),
intent(inout) :: this
258 fp => facet_list%array()
259 do i = 1, facet_list%size()
260 call this%marked_facet%push(fp(i))
268 class(
bc_t),
intent(inout) :: this
271 do i = 1, bc_zone%size
272 call this%marked_facet%push(bc_zone%facet_el(i))
283 class(
bc_t),
intent(inout) :: this
285 character(len=*) :: bc_key
286 character(len=100),
allocatable :: split_key(:)
287 character(len=NEKO_MSH_MAX_ZLBL_LEN) :: bc_labels(NEKO_MSH_MAX_ZLBLS)
288 integer :: i, j, k, l, msh_bc_type
291 if(trim(bc_key) .eq.
'o' .or. trim(bc_key) .eq.
'on' &
292 .or. trim(bc_key) .eq.
'o+dong' .or. trim(bc_key) .eq.
'on+dong')
then
294 else if(trim(bc_key) .eq.
'd_pres')
then
296 else if(trim(bc_key) .eq.
'w')
then
298 else if(trim(bc_key) .eq.
'v')
then
300 else if(trim(bc_key) .eq.
'd_vel_u')
then
302 else if(trim(bc_key) .eq.
'd_vel_v')
then
304 else if(trim(bc_key) .eq.
'd_vel_w')
then
306 else if(trim(bc_key) .eq.
'sym')
then
310 do i = 1, neko_msh_max_zlbls
313 if (index(trim(bc_labels(i)),
'/') .eq. 0)
then
314 if (trim(bc_key) .eq. trim(bc_labels(i)))
then
317 do j = 1,this%msh%nelv
318 do k = 1, 2 * this%msh%gdim
319 if (this%msh%facet_type(k,j) .eq. -i)
then
320 this%msh%facet_type(k,j) = msh_bc_type
327 do l = 1,
size(split_key)
328 if (trim(split_key(l)) .eq. trim(bc_key))
then
331 do j = 1,this%msh%nelv
332 do k = 1, 2 * this%msh%gdim
333 if (this%msh%facet_type(k,j) .eq. -i)
then
334 this%msh%facet_type(k,j) = msh_bc_type
349 class(
bc_t),
target,
intent(inout) :: this
352 integer :: facet_size, facet, el
353 integer :: i, j, k, l, msk_c
354 integer :: lx, ly, lz, n
364 allocate(this%msk(0:facet_size * this%marked_facet%size()))
365 allocate(this%facet(0:facet_size * this%marked_facet%size()))
368 bfp => this%marked_facet%array()
374 do i = 1, this%marked_facet%size()
376 facet = bc_facet%x(1)
384 this%facet(msk_c) = 1
392 this%facet(msk_c) = 2
400 this%facet(msk_c) = 3
408 this%facet(msk_c) = 4
416 this%facet(msk_c) = 5
424 this%facet(msk_c) = 6
431 this%facet(0) = msk_c
434 n = facet_size * this%marked_facet%size() + 1
449 type(
bc_list_t),
intent(inout),
target :: bclst
450 integer,
optional :: size
455 if (
present(size))
then
461 allocate(bclst%bc(n))
464 bclst%bc(i)%bcp => null()
478 if (
allocated(bclst%bc))
then
491 class(
bc_t),
intent(inout),
target ::
bc
492 type(
bcp_t),
allocatable :: tmp(:)
495 if(
bc%marked_facet%size() .eq. 0)
return
497 if (bclst%n .ge. bclst%size)
then
498 bclst%size = bclst%size * 2
499 allocate(tmp(bclst%size))
500 tmp(1:bclst%n) = bclst%bc
501 call move_alloc(tmp, bclst%bc)
504 bclst%n = bclst%n + 1
505 bclst%bc(bclst%n)%bcp =>
bc
516 integer,
intent(in) :: n
517 real(kind=
rp),
intent(inout),
dimension(n) :: x
518 real(kind=
rp),
intent(in),
optional :: t
519 integer,
intent(in),
optional :: tstep
525 if (
present(t) .and.
present(tstep))
then
527 call bclst%bc(i)%bcp%apply_scalar_dev(x_d, t=t, tstep=tstep)
529 else if (
present(t))
then
531 call bclst%bc(i)%bcp%apply_scalar_dev(x_d, t=t)
533 else if (
present(tstep))
then
535 call bclst%bc(i)%bcp%apply_scalar_dev(x_d, tstep=tstep)
539 call bclst%bc(i)%bcp%apply_scalar_dev(x_d)
543 if (
present(t) .and.
present(tstep))
then
545 call bclst%bc(i)%bcp%apply_scalar(x, n, t, tstep)
547 else if (
present(t))
then
549 call bclst%bc(i)%bcp%apply_scalar(x, n, t=t)
551 else if (
present(tstep))
then
553 call bclst%bc(i)%bcp%apply_scalar(x, n, tstep=tstep)
557 call bclst%bc(i)%bcp%apply_scalar(x, n)
573 integer,
intent(in) :: n
574 real(kind=
rp),
intent(inout),
dimension(n) :: x
575 real(kind=
rp),
intent(inout),
dimension(n) :: y
576 real(kind=
rp),
intent(inout),
dimension(n) :: z
577 real(kind=
rp),
intent(in),
optional :: t
578 integer,
intent(in),
optional :: tstep
588 if (
present(t) .and.
present(tstep))
then
590 call bclst%bc(i)%bcp%apply_vector_dev(x_d, y_d, z_d, t, tstep)
592 else if (
present(t))
then
594 call bclst%bc(i)%bcp%apply_vector_dev(x_d, y_d, z_d, t=t)
596 else if (
present(tstep))
then
598 call bclst%bc(i)%bcp%apply_vector_dev(x_d, y_d, z_d, tstep=tstep)
602 call bclst%bc(i)%bcp%apply_vector_dev(x_d, y_d, z_d)
606 if (
present(t) .and.
present(tstep))
then
608 call bclst%bc(i)%bcp%apply_vector(x, y, z, n, t, tstep)
610 else if (
present(t))
then
612 call bclst%bc(i)%bcp%apply_vector(x, y, z, n, t=t)
614 else if (
present(tstep))
then
616 call bclst%bc(i)%bcp%apply_vector(x, y, z, n, tstep=tstep)
620 call bclst%bc(i)%bcp%apply_vector(x, y, z, n)
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.
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, public bc_list_add(bclst, bc)
Add a condition to a list of boundary conditions.
subroutine, public bc_list_apply_vector(bclst, x, y, z, n, t, tstep)
Apply a list of boundary conditions to a vector field.
subroutine bc_finalize(this)
Finalize the construction of the bc by populting the msk and facet arrays.
subroutine, public bc_list_init(bclst, size)
Constructor for a list of boundary conditions.
subroutine bc_free(this)
Destructor.
subroutine, public bc_list_free(bclst)
Destructor for a list of boundary conditions.
subroutine bc_mark_zones_from_list(this, bc_zones, bc_key, bc_labels)
Mark all facets from a list of zones, also marks type of bc in the mesh. The facet_type in mesh is be...
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, public bc_list_apply_scalar(bclst, x, n, t, tstep)
Apply a list of boundary conditions to a scalar field.
subroutine bc_init(this, coef)
Constructor.
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.
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 split_string(string, delimiter)
Split a string based on delimiter (tokenizer) OBS: very hacky, this should really be improved,...
pure integer function 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, :)
A list of boundary conditions.
Base type for a boundary condition.
Pointer to boundary condtiion.
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
The function space for the SEM solution fields.
Integer 2-tuple based stack.