48 type(c_ptr),
private :: usr_x_d = c_null_ptr
49 type(c_ptr),
private :: usr_y_d = c_null_ptr
50 type(c_ptr),
private :: usr_z_d = c_null_ptr
80 ix, iy, iz, ie, t, tstep)
82 real(kind=
rp),
intent(inout) :: u
83 real(kind=
rp),
intent(inout) :: v
84 real(kind=
rp),
intent(inout) :: w
85 real(kind=
rp),
intent(in) :: x
86 real(kind=
rp),
intent(in) :: y
87 real(kind=
rp),
intent(in) :: z
88 real(kind=
rp),
intent(in) :: nx
89 real(kind=
rp),
intent(in) :: ny
90 real(kind=
rp),
intent(in) :: nz
91 integer,
intent(in) :: ix
92 integer,
intent(in) :: iy
93 integer,
intent(in) :: iz
94 integer,
intent(in) :: ie
95 real(kind=
rp),
intent(in) :: t
96 integer,
intent(in) :: tstep
107 if (c_associated(this%usr_x_d))
then
111 if (c_associated(this%usr_y_d))
then
115 if (c_associated(this%usr_z_d))
then
124 integer,
intent(in) :: n
125 real(kind=
rp),
intent(inout),
dimension(n) :: x
126 real(kind=
rp),
intent(in),
optional :: t
127 integer,
intent(in),
optional :: tstep
134 real(kind=
rp),
intent(in),
optional :: t
135 integer,
intent(in),
optional :: tstep
141 integer,
intent(in) :: n
142 real(kind=
rp),
intent(inout),
dimension(n) :: x
143 real(kind=
rp),
intent(inout),
dimension(n) :: y
144 real(kind=
rp),
intent(inout),
dimension(n) :: z
145 real(kind=
rp),
intent(in),
optional :: t
146 integer,
intent(in),
optional :: tstep
147 integer :: i, m, k, idx(4), facet, tstep_
156 if (
present(tstep))
then
162 associate(xc => this%c%dof%x, yc => this%c%dof%y, zc => this%c%dof%z, &
163 nx => this%c%nx, ny => this%c%ny, nz => this%c%nz, &
168 facet = this%facet(i)
172 call this%eval(x(k), y(k), z(k), &
173 xc(idx(1), idx(2), idx(3), idx(4)), &
174 yc(idx(1), idx(2), idx(3), idx(4)), &
175 zc(idx(1), idx(2), idx(3), idx(4)), &
176 nx(idx(2), idx(3), facet, idx(4)), &
177 ny(idx(2), idx(3), facet, idx(4)), &
178 nz(idx(2), idx(3), facet, idx(4)), &
179 idx(1), idx(2), idx(3), idx(4), &
182 call this%eval(x(k), y(k), z(k), &
183 xc(idx(1), idx(2), idx(3), idx(4)), &
184 yc(idx(1), idx(2), idx(3), idx(4)), &
185 zc(idx(1), idx(2), idx(3), idx(4)), &
186 nx(idx(1), idx(3), facet, idx(4)), &
187 ny(idx(1), idx(3), facet, idx(4)), &
188 nz(idx(1), idx(3), facet, idx(4)), &
189 idx(1), idx(2), idx(3), idx(4), &
192 call this%eval(x(k), y(k), z(k), &
193 xc(idx(1), idx(2), idx(3), idx(4)), &
194 yc(idx(1), idx(2), idx(3), idx(4)), &
195 zc(idx(1), idx(2), idx(3), idx(4)), &
196 nx(idx(1), idx(2), facet, idx(4)), &
197 ny(idx(1), idx(2), facet, idx(4)), &
198 nz(idx(1), idx(2), facet, idx(4)), &
199 idx(1), idx(2), idx(3), idx(4), &
212 real(kind=rp),
intent(in),
optional :: t
213 integer,
intent(in),
optional :: tstep
214 integer :: i, m, k, idx(4), facet, tstep_
215 integer(c_size_t) :: s
217 real(kind=rp),
allocatable :: x(:)
218 real(kind=rp),
allocatable :: y(:)
219 real(kind=rp),
allocatable :: z(:)
227 if (
present(tstep))
then
233 associate(xc => this%c%dof%x, yc => this%c%dof%y, zc => this%c%dof%z, &
234 nx => this%c%nx, ny => this%c%ny, nz => this%c%nz, &
235 lx => this%c%Xh%lx, usr_x_d => this%usr_x_d, usr_y_d => this%usr_y_d, &
236 usr_z_d => this%usr_z_d)
242 if (.not. c_associated(usr_x_d))
then
243 allocate(x(m), y(m), z(m))
247 call device_alloc(usr_x_d, s)
248 call device_alloc(usr_y_d, s)
249 call device_alloc(usr_z_d, s)
251 associate(xc => this%c%dof%x, yc => this%c%dof%y, zc => this%c%dof%z, &
252 nx => this%c%nx, ny => this%c%ny, nz => this%c%nz, &
256 facet = this%facet(i)
260 call this%eval(x(i), y(i), z(i), &
261 xc(idx(1), idx(2), idx(3), idx(4)), &
262 yc(idx(1), idx(2), idx(3), idx(4)), &
263 zc(idx(1), idx(2), idx(3), idx(4)), &
264 nx(idx(2), idx(3), facet, idx(4)), &
265 ny(idx(2), idx(3), facet, idx(4)), &
266 nz(idx(2), idx(3), facet, idx(4)), &
267 idx(1), idx(2), idx(3), idx(4), &
270 call this%eval(x(i), y(i), z(i), &
271 xc(idx(1), idx(2), idx(3), idx(4)), &
272 yc(idx(1), idx(2), idx(3), idx(4)), &
273 zc(idx(1), idx(2), idx(3), idx(4)), &
274 nx(idx(1), idx(3), facet, idx(4)), &
275 ny(idx(1), idx(3), facet, idx(4)), &
276 nz(idx(1), idx(3), facet, idx(4)), &
277 idx(1), idx(2), idx(3), idx(4), &
280 call this%eval(x(i), y(i), z(i), &
281 xc(idx(1), idx(2), idx(3), idx(4)), &
282 yc(idx(1), idx(2), idx(3), idx(4)), &
283 zc(idx(1), idx(2), idx(3), idx(4)), &
284 nx(idx(1), idx(2), facet, idx(4)), &
285 ny(idx(1), idx(2), facet, idx(4)), &
286 nz(idx(1), idx(2), facet, idx(4)), &
287 idx(1), idx(2), idx(3), idx(4), &
293 call device_memcpy(x, usr_x_d, m, host_to_device, sync=.false.)
294 call device_memcpy(y, usr_y_d, m, host_to_device, sync=.false.)
295 call device_memcpy(z, usr_z_d, m, host_to_device, sync=.true.)
300 call device_inhom_dirichlet_apply_vector(this%msk_d, x_d, y_d, z_d, &
301 usr_x_d, usr_y_d, usr_z_d, m)
310 type(coef_t),
target,
intent(inout) :: c
319 this%eval => usr_eval
328 if (.not.
associated(this%c))
then
329 call neko_warning(
'Missing coefficients')
333 if (.not.
associated(this%eval))
then
334 call neko_warning(
'Missing eval function')
338 if (.not. valid)
then
339 call neko_error(
'Invalid user defined inflow condition')
__device__ void nonlinear_index(const int idx, const int lx, int *index)
Abstract interface defining a user defined inflow condition (pointwise)
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.
Defines inflow dirichlet conditions.
subroutine usr_inflow_apply_scalar(this, x, n, t, tstep)
No-op scalar apply.
subroutine usr_inflow_apply_vector_dev(this, x_d, y_d, z_d, t, tstep)
subroutine usr_inflow_apply_vector(this, x, y, z, n, t, tstep)
Apply user defined inflow conditions (vector valued)
subroutine usr_inflow_apply_scalar_dev(this, x_d, t, tstep)
No-op scalar apply (device version)
subroutine usr_inflow_validate(this)
Validate user inflow condition.
subroutine usr_inflow_free(this)
subroutine usr_inflow_set_eval(this, usr_eval)
Assign user provided eval function.
subroutine usr_inflow_set_coef(this, c)
Assign coefficients (facet normals etc)
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
Dirichlet condition for inlet (vector valued)
User defined dirichlet condition for inlet (vector valued)