41 use,
intrinsic :: iso_fortran_env
42 use,
intrinsic :: iso_c_binding
44 use json_module,
only : json_file
54 real(kind=
rp),
dimension(3) :: uinf = [0d0, 0d0, 0d0]
55 real(kind=
rp) :: delta
57 type(c_ptr),
private :: blax_d = c_null_ptr
58 type(c_ptr),
private :: blay_d = c_null_ptr
59 type(c_ptr),
private :: blaz_d = c_null_ptr
69 procedure, pass(this) :: init_from_components => &
83 class(
blasius_t),
intent(inout),
target :: this
84 type(
coef_t),
target,
intent(in) :: coef
85 type(json_file),
intent(inout) :: json
86 real(kind=
rp) :: delta
87 real(kind=
rp),
allocatable :: uinf(:)
88 character(len=:),
allocatable :: approximation
90 call this%init_base(coef)
93 call json_get(json,
'approximation', approximation)
94 call json_get(json,
'freestream_velocity', uinf)
96 if (
size(uinf) .ne. 3)
then
97 call neko_error(
"The uinf keyword for the blasius profile should be an &
101 call this%init_from_components(coef, delta, uinf, approximation)
112 class(
blasius_t),
intent(inout),
target :: this
113 type(
coef_t),
target,
intent(in) :: coef
114 real(kind=
rp) :: delta
115 real(kind=
rp) :: uinf(3)
116 character(len=*) :: approximation
118 call this%init_base(coef)
123 select case (trim(approximation))
137 call neko_error(
'Invalid Blasius approximation')
142 class(
blasius_t),
target,
intent(inout) :: this
144 call this%free_base()
147 if (c_associated(this%blax_d))
then
151 if (c_associated(this%blay_d))
then
155 if (c_associated(this%blaz_d))
then
164 integer,
intent(in) :: n
165 real(kind=
rp),
intent(inout),
dimension(n) :: x
167 logical,
intent(in),
optional :: strong
172 class(
blasius_t),
intent(inout),
target :: this
173 type(c_ptr),
intent(inout) :: x_d
175 logical,
intent(in),
optional :: strong
176 type(c_ptr),
intent(inout) :: strm
182 integer,
intent(in) :: n
183 real(kind=
rp),
intent(inout),
dimension(n) :: x
184 real(kind=
rp),
intent(inout),
dimension(n) :: y
185 real(kind=
rp),
intent(inout),
dimension(n) :: z
187 logical,
intent(in),
optional :: strong
188 integer :: i, m, k, idx(4), facet
191 if (
present(strong))
then
197 associate(xc => this%coef%dof%x, yc => this%coef%dof%y, &
198 zc => this%coef%dof%z, nx => this%coef%nx, ny => this%coef%ny, &
199 nz => this%coef%nz, lx => this%coef%Xh%lx)
204 facet = this%facet(i)
208 x(k) = this%bla(zc(idx(1), idx(2), idx(3), idx(4)), &
209 this%delta, this%uinf(1))
214 y(k) = this%bla(xc(idx(1), idx(2), idx(3), idx(4)), &
215 this%delta, this%uinf(2))
220 z(k) = this%bla(yc(idx(1), idx(2), idx(3), idx(4)), &
221 this%delta, this%uinf(3))
230 class(
blasius_t),
intent(inout),
target :: this
231 type(c_ptr),
intent(inout) :: x_d
232 type(c_ptr),
intent(inout) :: y_d
233 type(c_ptr),
intent(inout) :: z_d
234 type(time_state_t),
intent(in),
optional :: time
235 logical,
intent(in),
optional :: strong
236 integer :: i, m, k, idx(4), facet
237 integer(c_size_t) :: s
238 real(kind=rp),
allocatable :: bla_x(:), bla_y(:), bla_z(:)
240 type(c_ptr),
intent(inout) :: strm
242 if (
present(strong))
then
248 associate(xc => this%coef%dof%x, yc => this%coef%dof%y, &
249 zc => this%coef%dof%z, nx => this%coef%nx, ny => this%coef%ny, &
250 nz => this%coef%nz, lx => this%coef%Xh%lx , &
251 blax_d => this%blax_d, blay_d => this%blay_d, &
252 blaz_d => this%blaz_d)
258 if (.not. c_associated(blax_d) .and. strong_ .and. this%msk(0) .gt. 0)
then
259 allocate(bla_x(m), bla_y(m), bla_z(m))
261 if (rp .eq. real32)
then
263 else if (rp .eq. real64)
then
267 call device_alloc(blax_d, s)
268 call device_alloc(blay_d, s)
269 call device_alloc(blaz_d, s)
273 facet = this%facet(i)
277 bla_x(i) = this%bla(zc(idx(1), idx(2), idx(3), idx(4)), &
278 this%delta, this%uinf(1))
283 bla_y(i) = this%bla(xc(idx(1), idx(2), idx(3), idx(4)), &
284 this%delta, this%uinf(2))
289 bla_z(i) = this%bla(yc(idx(1), idx(2), idx(3), idx(4)), &
290 this%delta, this%uinf(3))
294 call device_memcpy(bla_x, blax_d, m, host_to_device, sync = .false.)
295 call device_memcpy(bla_y, blay_d, m, host_to_device, sync = .false.)
296 call device_memcpy(bla_z, blaz_d, m, host_to_device, sync = .true.)
298 deallocate(bla_x, bla_y, bla_z)
301 if (strong_ .and. this%msk(0) .gt. 0)
then
302 call device_inhom_dirichlet_apply_vector(this%msk_d, x_d, y_d, z_d, &
303 blax_d, blay_d, blaz_d, m, strm)
313 real(kind=rp),
intent(in) :: uinf(3)
314 real(kind=rp),
intent(in) :: delta
315 character(len=*) :: type
319 select case (trim(type))
321 this%bla => blasius_linear
323 this%bla => blasius_quadratic
325 this%bla => blasius_cubic
327 this%bla => blasius_quartic
329 this%bla => blasius_sin
331 this%bla => blasius_tanh
333 call neko_error(
'Invalid Blasius approximation')
339 class(
blasius_t),
target,
intent(inout) :: this
340 logical,
optional,
intent(in) :: only_facets
341 logical :: only_facets_
343 if (
present(only_facets))
then
344 only_facets_ = only_facets
346 only_facets_ = .false.
349 call this%finalize_base(only_facets_)
__inline__ __device__ void nonlinear_index(const int idx, const int lx, int *index)
Copy data between host and device (or device and device)
Abstract interface for computing a Blasius flow profile.
Retrieves a parameter by name or throws an error.
Defines a boundary condition.
Defines a Blasius profile dirichlet condition.
subroutine blasius_finalize(this, only_facets)
Finalize.
subroutine blasius_apply_scalar_dev(this, x_d, time, strong, strm)
No-op scalar apply (device version)
subroutine blasius_init(this, coef, json)
Constructor.
subroutine blasius_init_from_components(this, coef, delta, uinf, approximation)
Constructor from components.
subroutine blasius_apply_vector(this, x, y, z, n, time, strong)
Apply blasius conditions (vector valued)
subroutine blasius_apply_scalar(this, x, n, time, strong)
No-op scalar apply.
subroutine blasius_apply_vector_dev(this, x_d, y_d, z_d, time, strong, strm)
Apply blasius conditions (vector valued) (device version)
subroutine blasius_free(this)
subroutine blasius_set_params(this, uinf, delta, type)
Set Blasius parameters.
Device abstraction, common interface for various accelerators.
integer, parameter, public host_to_device
subroutine, public device_free(x_d)
Deallocate memory on the device.
subroutine, public device_alloc(x_d, s)
Allocate memory on the device.
real(kind=rp) function, public blasius_quadratic(y, delta, u)
Quadratic approximate Blasius Profile .
real(kind=rp) function, public blasius_quartic(y, delta, u)
Quartic approximate Blasius Profile .
real(kind=rp) function, public blasius_sin(y, delta, u)
Sinusoidal approximate Blasius Profile .
real(kind=rp) function, public blasius_cubic(y, delta, u)
Cubic approximate Blasius Profile .
real(kind=rp) function, public blasius_tanh(y, delta, u)
Hyperbolic tangent approximate Blasius Profile from O. Savas (2012) where is the 99 percent thickne...
real(kind=rp) function, public blasius_linear(y, delta, u)
Linear approximate Blasius profile .
Utilities for retrieving parameters from the case files.
integer, parameter, public rp
Global precision used in computations.
Module with things related to the simulation time.
Base type for a boundary condition.
Blasius profile for inlet (vector valued).
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
A struct that contains all info about the time, expand as needed.