42 use,
intrinsic :: iso_fortran_env
43 use,
intrinsic :: iso_c_binding
50 real(kind=
rp) :: delta
52 type(c_ptr),
private :: blax_d = c_null_ptr
53 type(c_ptr),
private :: blay_d = c_null_ptr
54 type(c_ptr),
private :: blaz_d = c_null_ptr
71 if (c_associated(this%blax_d))
then
75 if (c_associated(this%blay_d))
then
79 if (c_associated(this%blaz_d))
then
88 integer,
intent(in) :: n
89 real(kind=
rp),
intent(inout),
dimension(n) :: x
90 real(kind=
rp),
intent(in),
optional :: t
91 integer,
intent(in),
optional :: tstep
96 class(
blasius_t),
intent(inout),
target :: this
98 real(kind=
rp),
intent(in),
optional :: t
99 integer,
intent(in),
optional :: tstep
105 integer,
intent(in) :: n
106 real(kind=
rp),
intent(inout),
dimension(n) :: x
107 real(kind=
rp),
intent(inout),
dimension(n) :: y
108 real(kind=
rp),
intent(inout),
dimension(n) :: z
109 real(kind=
rp),
intent(in),
optional :: t
110 integer,
intent(in),
optional :: tstep
111 integer :: i, m, k, idx(4), facet
113 associate(xc => this%c%dof%x, yc => this%c%dof%y, zc => this%c%dof%z, &
114 nx => this%c%nx, ny => this%c%ny, nz => this%c%nz, &
119 facet = this%facet(i)
123 x(k) = this%bla(zc(idx(1), idx(2), idx(3), idx(4)), &
124 this%delta, this%x(1))
129 y(k) = this%bla(xc(idx(1), idx(2), idx(3), idx(4)), &
130 this%delta, this%x(2))
135 z(k) = this%bla(yc(idx(1), idx(2), idx(3), idx(4)), &
136 this%delta, this%x(3))
144 class(
blasius_t),
intent(inout),
target :: this
148 real(kind=rp),
intent(in),
optional :: t
149 integer,
intent(in),
optional :: tstep
150 integer :: i, m, k, idx(4), facet
151 integer(c_size_t) :: s
152 real(kind=rp),
allocatable :: bla_x(:), bla_y(:), bla_z(:)
154 associate(xc => this%c%dof%x, yc => this%c%dof%y, zc => this%c%dof%z, &
155 nx => this%c%nx, ny => this%c%ny, nz => this%c%nz, &
156 lx => this%c%Xh%lx, blax_d => this%blax_d, blay_d => this%blay_d, &
157 blaz_d => this%blaz_d)
163 if (.not. c_associated(blax_d))
then
164 allocate(bla_x(m), bla_y(m), bla_z(m))
166 if (rp .eq. real32)
then
168 else if (rp .eq. real64)
then
172 call device_alloc(blax_d, s)
173 call device_alloc(blay_d, s)
174 call device_alloc(blaz_d, s)
178 facet = this%facet(i)
182 bla_x(i) = this%bla(zc(idx(1), idx(2), idx(3), idx(4)), &
183 this%delta, this%x(1))
188 bla_y(i) = this%bla(xc(idx(1), idx(2), idx(3), idx(4)), &
189 this%delta, this%x(2))
194 bla_z(i) = this%bla(yc(idx(1), idx(2), idx(3), idx(4)), &
195 this%delta, this%x(3))
199 call device_memcpy(bla_x, blax_d, m, host_to_device, sync=.false.)
200 call device_memcpy(bla_y, blay_d, m, host_to_device, sync=.false.)
201 call device_memcpy(bla_z, blaz_d, m, host_to_device, sync=.true.)
203 deallocate(bla_x, bla_y, bla_z)
206 call device_inhom_dirichlet_apply_vector(this%msk_d, x_d, y_d, z_d, &
207 blax_d, blay_d, blaz_d, m)
216 real(kind=rp) :: delta
217 character(len=*) :: type
220 select case(trim(type))
222 this%bla => blasius_linear
224 this%bla => blasius_quadratic
226 this%bla => blasius_cubic
228 this%bla => blasius_quartic
230 this%bla => blasius_sin
232 call neko_error(
'Invalid Blasius approximation')
239 type(coef_t),
target,
intent(inout) :: c
__device__ void nonlinear_index(const int idx, const int lx, int *index)
Abstract interface for computing a Blasius flow profile.
Defines a Blasius profile dirichlet condition.
subroutine blasius_set_coef(this, c)
Assign coefficients (facet normals etc)
subroutine blasius_apply_vector(this, x, y, z, n, t, tstep)
Apply blasius conditions (vector valued)
subroutine blasius_apply_vector_dev(this, x_d, y_d, z_d, t, tstep)
Apply blasius conditions (vector valued) (device version)
subroutine blasius_apply_scalar_dev(this, x_d, t, tstep)
No-op scalar apply (device version)
subroutine blasius_apply_scalar(this, x, n, t, tstep)
No-op scalar apply.
subroutine blasius_free(this)
subroutine blasius_set_params(this, delta, type)
Set Blasius parameters.
Device abstraction, common interface for various accelerators.
subroutine, public device_free(x_d)
Deallocate memory on the device.
Defines inflow dirichlet conditions.
integer, parameter, public rp
Global precision used in computations.
Blasius profile for inlet (vector valued)
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
Dirichlet condition for inlet (vector valued)