38 use json_module,
only : json_file
51 use,
intrinsic :: iso_c_binding, only : c_ptr, c_null_ptr, c_associated
61 type(
coef_t),
pointer :: coef => null()
71 character(len=:),
allocatable :: scheme_name
73 integer,
pointer :: msk(:) => null()
74 type(c_ptr) :: msk_d = c_null_ptr
76 integer,
pointer :: facet(:) => null()
90 integer,
allocatable :: ind_r(:)
91 type(c_ptr) :: ind_r_d = c_null_ptr
93 integer,
allocatable :: ind_s(:)
94 type(c_ptr) :: ind_s_d = c_null_ptr
96 integer,
allocatable :: ind_t(:)
97 type(c_ptr) :: ind_t_d = c_null_ptr
99 integer,
allocatable :: ind_e(:)
100 type(c_ptr) :: ind_e_d = c_null_ptr
104 integer :: h_index = 0
106 integer :: n_nodes = 0
147 real(kind=
rp),
intent(in) :: t
148 integer,
intent(in) :: tstep
165 character(len=*),
intent(in) :: scheme_name
166 type(
coef_t),
intent(in) :: coef
167 integer,
intent(in) :: msk(:)
168 integer,
intent(in) :: facet(:)
169 integer,
intent(in) :: h_index
170 type(json_file),
intent(inout) :: json
182 type(
coef_t),
intent(in) :: coef
183 type(json_file),
intent(inout) :: json
194 integer,
intent(in) :: msk(:)
195 integer,
intent(in) :: facet(:)
217 module subroutine wall_model_factory(object, scheme_name, coef, msk, &
219 class(wall_model_t),
allocatable,
intent(inout) :: object
220 character(len=*),
intent(in) :: scheme_name
221 type(
coef_t),
intent(in) :: coef
222 integer,
intent(in) :: msk(:)
223 integer,
intent(in) :: facet(:)
224 type(json_file),
intent(inout) :: json
225 end subroutine wall_model_factory
232 module subroutine wall_model_allocator(object, type_name)
233 class(wall_model_t),
allocatable,
intent(inout) :: object
234 character(len=*),
intent(in) :: type_name
235 end subroutine wall_model_allocator
246 subroutine wall_model_allocate(obj)
248 class(wall_model_t),
allocatable,
intent(inout) :: obj
249 end subroutine wall_model_allocate
254 module subroutine register_wall_model(type_name, allocator)
255 character(len=*),
intent(in) :: type_name
256 procedure(wall_model_allocate),
pointer,
intent(in) :: allocator
257 end subroutine register_wall_model
263 character(len=20) :: type_name
264 procedure(wall_model_allocate),
pointer,
nopass :: allocator
265 end type allocator_entry
268 type(allocator_entry),
allocatable :: wall_model_registry(:)
271 integer :: wall_model_registry_size = 0
273 public :: wall_model_factory, wall_model_allocator, register_wall_model, &
283 subroutine wall_model_init_base(this, scheme_name, coef, msk, facet, index)
284 class(wall_model_t),
intent(inout) :: this
285 type(
coef_t),
target,
intent(in) :: coef
286 integer,
target,
intent(in) :: msk(0:)
287 integer,
target,
intent(in) :: facet(0:)
288 character(len=*) :: scheme_name
289 integer,
intent(in) :: index
296 this%scheme_name = trim(scheme_name)
302 ignore_existing = .true.)
305 call this%finalize_base(msk, facet)
306 end subroutine wall_model_init_base
312 subroutine wall_model_partial_init_base(this, coef, json)
313 class(wall_model_t),
intent(inout) :: this
314 type(
coef_t),
target,
intent(in) :: coef
315 type(json_file),
intent(inout) :: json
317 call this%free_base()
321 call json_get(json,
"h_index", this%h_index)
322 call json_get(json,
"scheme_name", this%scheme_name)
329 ignore_existing = .true.)
331 end subroutine wall_model_partial_init_base
333 subroutine wall_model_finalize_base(this, msk, facet)
334 class(wall_model_t),
intent(inout) :: this
335 integer,
target,
intent(in) :: msk(0:)
336 integer,
target,
intent(in) :: facet(:)
338 this%msk(0:msk(0)) => msk
340 this%facet(0:msk(0)) => facet
342 call this%tau_x%init(this%msk(0))
343 call this%tau_y%init(this%msk(0))
344 call this%tau_z%init(this%msk(0))
346 allocate(this%ind_r(this%msk(0)))
347 allocate(this%ind_s(this%msk(0)))
348 allocate(this%ind_t(this%msk(0)))
349 allocate(this%ind_e(this%msk(0)))
351 call this%h%init(this%msk(0))
352 call this%n_x%init(this%msk(0))
353 call this%n_y%init(this%msk(0))
354 call this%n_z%init(this%msk(0))
356 call this%find_points()
360 call device_map(this%ind_r, this%ind_r_d, this%n_nodes)
361 call device_map(this%ind_s, this%ind_s_d, this%n_nodes)
362 call device_map(this%ind_t, this%ind_t_d, this%n_nodes)
363 call device_map(this%ind_e, this%ind_e_d, this%n_nodes)
374 end subroutine wall_model_finalize_base
377 subroutine wall_model_free_base(this)
378 class(wall_model_t),
intent(inout) :: this
383 nullify(this%tau_field)
387 call this%tau_x%free()
388 call this%tau_y%free()
389 call this%tau_z%free()
391 if (
allocated(this%ind_r))
then
392 deallocate(this%ind_r)
394 if (
allocated(this%ind_s))
then
395 deallocate(this%ind_s)
397 if (
allocated(this%ind_t))
then
398 deallocate(this%ind_t)
401 if (c_associated(this%msk_d))
then
404 if (c_associated(this%ind_r_d))
then
407 if (c_associated(this%ind_s_d))
then
410 if (c_associated(this%ind_t_d))
then
413 if (c_associated(this%ind_e_d))
then
421 end subroutine wall_model_free_base
424 subroutine wall_model_find_points(this)
425 class(wall_model_t),
intent(inout) :: this
426 integer :: n_nodes, fid, idx(4), i, linear
427 real(kind=
rp) :: normal(3), p(3), x, y, z, xw, yw, zw, magp
428 real(kind=
rp) :: hmin, hmax
429 type(
field_t),
pointer :: h_field
431 character(len=LOG_SIZE),
allocatable :: log_msg
433 n_nodes = this%msk(0)
434 this%n_nodes = n_nodes
437 ignore_existing=.true.)
446 normal = this%coef%get_normal(idx(1), idx(2), idx(3), idx(4), fid)
448 this%n_x%x(i) = normal(1)
449 this%n_y%x(i) = normal(2)
450 this%n_z%x(i) = normal(3)
457 this%ind_r(i) = idx(1) + this%h_index
458 this%ind_s(i) = idx(2)
459 this%ind_t(i) = idx(3)
461 this%ind_r(i) = idx(1) - this%h_index
462 this%ind_s(i) = idx(2)
463 this%ind_t(i) = idx(3)
465 this%ind_r(i) = idx(1)
466 this%ind_s(i) = idx(2) + this%h_index
467 this%ind_t(i) = idx(3)
469 this%ind_r(i) = idx(1)
470 this%ind_s(i) = idx(2) - this%h_index
471 this%ind_t(i) = idx(3)
473 this%ind_r(i) = idx(1)
474 this%ind_s(i) = idx(2)
475 this%ind_t(i) = idx(3) + this%h_index
477 this%ind_r(i) = idx(1)
478 this%ind_s(i) = idx(2)
479 this%ind_t(i) = idx(3) - this%h_index
481 call neko_error(
"The face index is not correct ")
483 this%ind_e(i) = idx(4)
486 xw = this%dof%x(idx(1), idx(2), idx(3), idx(4))
487 yw = this%dof%y(idx(1), idx(2), idx(3), idx(4))
488 zw = this%dof%z(idx(1), idx(2), idx(3), idx(4))
491 x = this%dof%x(this%ind_r(i), this%ind_s(i), this%ind_t(i), &
493 y = this%dof%y(this%ind_r(i), this%ind_s(i), this%ind_t(i), &
495 z = this%dof%z(this%ind_r(i), this%ind_s(i), this%ind_t(i), &
504 magp = sqrt(p(1)**2 + p(2)**2 + p(3)**2)
507 this%h%x(i) = p(1)*normal(1) + p(2)*normal(2) + p(3)*normal(3)
509 h_field%x(linear,1,1,1) = this%h%x(i)
513 if ((this%h%x(i) - magp) / magp > 0.1)
then
514 write(log_msg,*)
"Significant misalignment between wall normal and"
516 write(log_msg,*)
"sampling point direction at wall node", xw, yw, zw
540 call h_file%init(
"sampling_height.fld")
541 call h_file%write(h_field)
542 end subroutine wall_model_find_points
545 class(wall_model_t),
intent(inout) :: this
547 real(kind=
rp) :: magtau
555 this%tau_field%x_d, &
559 magtau = sqrt(this%tau_x%x(i)**2 + &
560 this%tau_y%x(i)**2 + &
562 this%tau_field%x(this%msk(i),1,1,1) = magtau
__inline__ __device__ void nonlinear_index(const int idx, const int lx, int *index)
__global__ void wall_model_compute_mag_field(const T *__restrict__ tau_x_d, const T *__restrict__ tau_y_d, const T *__restrict__ tau_z_d, T *__restrict__ tau_field_d, const int *__restrict__ msk_d, const int m)
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)
Retrieves a parameter by name or throws an error.
Compute wall shear stress.
Finilzation of partial construction, similar to bc_t
Partial constructor from JSON, meant to work as the first stage of initialization before the finalize...
integer, public pe_rank
MPI rank.
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 registry for storing solution fields.
type(field_registry_t), target, public neko_field_registry
Global field registry.
Module for file I/O operations.
Utilities for retrieving parameters from the case files.
integer, parameter, public neko_log_debug
Debug log level.
type(log_t), public neko_log
Global log stream.
integer, parameter, public log_size
real(kind=rp) function, public glmax(a, n)
Max of a vector of length n.
real(kind=rp) function, public glmin(a, n)
Min of a vector of length n.
integer, parameter neko_bcknd_device
integer, parameter, public rp
Global precision used in computations.
Implements the device kernel for the wall_model_t type.
subroutine, public wall_model_compute_mag_field_device(tau_x_d, tau_y_d, tau_z_d, tau_field_d, msk_d, m)
Compute the wall shear stress's magnitude on device.
subroutine wall_model_finalize_base(this, msk, facet)
subroutine wall_model_init_base(this, scheme_name, coef, msk, facet, index)
Constructor for the wall_model_t (base) class.
subroutine wall_model_partial_init_base(this, coef, json)
Partial initialization based on JSON, prior to knowing the mask and facets.
subroutine wall_model_find_points(this)
Find sampling points based on the requested index.
subroutine wall_model_free_base(this)
Destructor for the wall_model_t (base) class.
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...
Base abstract type for wall-stress models for wall-modelled LES.