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
81 ix, iy, iz, ie, t, tstep)
83 real(kind=
rp),
intent(inout) :: u
84 real(kind=
rp),
intent(inout) :: v
85 real(kind=
rp),
intent(inout) :: w
86 real(kind=
rp),
intent(in) :: x
87 real(kind=
rp),
intent(in) :: y
88 real(kind=
rp),
intent(in) :: z
89 real(kind=
rp),
intent(in) :: nx
90 real(kind=
rp),
intent(in) :: ny
91 real(kind=
rp),
intent(in) :: nz
92 integer,
intent(in) :: ix
93 integer,
intent(in) :: iy
94 integer,
intent(in) :: iz
95 integer,
intent(in) :: ie
96 real(kind=
rp),
intent(in) :: t
97 integer,
intent(in) :: tstep
108 call this%free_base()
110 if (c_associated(this%usr_x_d))
then
114 if (c_associated(this%usr_y_d))
then
118 if (c_associated(this%usr_z_d))
then
127 integer,
intent(in) :: n
128 real(kind=
rp),
intent(inout),
dimension(n) :: x
129 real(kind=
rp),
intent(in),
optional :: t
130 integer,
intent(in),
optional :: tstep
137 real(kind=
rp),
intent(in),
optional :: t
138 integer,
intent(in),
optional :: tstep
144 integer,
intent(in) :: n
145 real(kind=
rp),
intent(inout),
dimension(n) :: x
146 real(kind=
rp),
intent(inout),
dimension(n) :: y
147 real(kind=
rp),
intent(inout),
dimension(n) :: z
148 real(kind=
rp),
intent(in),
optional :: t
149 integer,
intent(in),
optional :: tstep
150 integer :: i, m, k, idx(4), facet, tstep_
159 if (
present(tstep))
then
166 xc => this%coef%dof%x, yc => this%coef%dof%y, zc => this%coef%dof%z, &
167 nx => this%coef%nx, ny => this%coef%ny, nz => this%coef%nz, &
168 lx => this%coef%Xh%lx)
172 facet = this%facet(i)
176 call this%eval(x(k), y(k), z(k), &
177 xc(idx(1), idx(2), idx(3), idx(4)), &
178 yc(idx(1), idx(2), idx(3), idx(4)), &
179 zc(idx(1), idx(2), idx(3), idx(4)), &
180 nx(idx(2), idx(3), facet, idx(4)), &
181 ny(idx(2), idx(3), facet, idx(4)), &
182 nz(idx(2), idx(3), facet, idx(4)), &
183 idx(1), idx(2), idx(3), idx(4), &
186 call this%eval(x(k), y(k), z(k), &
187 xc(idx(1), idx(2), idx(3), idx(4)), &
188 yc(idx(1), idx(2), idx(3), idx(4)), &
189 zc(idx(1), idx(2), idx(3), idx(4)), &
190 nx(idx(1), idx(3), facet, idx(4)), &
191 ny(idx(1), idx(3), facet, idx(4)), &
192 nz(idx(1), idx(3), facet, idx(4)), &
193 idx(1), idx(2), idx(3), idx(4), &
196 call this%eval(x(k), y(k), z(k), &
197 xc(idx(1), idx(2), idx(3), idx(4)), &
198 yc(idx(1), idx(2), idx(3), idx(4)), &
199 zc(idx(1), idx(2), idx(3), idx(4)), &
200 nx(idx(1), idx(2), facet, idx(4)), &
201 ny(idx(1), idx(2), facet, idx(4)), &
202 nz(idx(1), idx(2), facet, idx(4)), &
203 idx(1), idx(2), idx(3), idx(4), &
216 real(kind=rp),
intent(in),
optional :: t
217 integer,
intent(in),
optional :: tstep
218 integer :: i, m, k, idx(4), facet, tstep_
219 integer(c_size_t) :: s
221 real(kind=rp),
allocatable :: x(:)
222 real(kind=rp),
allocatable :: y(:)
223 real(kind=rp),
allocatable :: z(:)
231 if (
present(tstep))
then
238 xc => this%coef%dof%x, yc => this%coef%dof%y, zc => this%coef%dof%z, &
239 nx => this%coef%nx, ny => this%coef%ny, nz => this%coef%nz, &
240 lx => this%coef%Xh%lx, usr_x_d => this%usr_x_d, &
241 usr_y_d => this%usr_y_d, usr_z_d => this%usr_z_d)
247 if (.not. c_associated(usr_x_d))
then
248 allocate(x(m), y(m), z(m))
252 call device_alloc(usr_x_d, s)
253 call device_alloc(usr_y_d, s)
254 call device_alloc(usr_z_d, s)
256 associate(xc => this%coef%dof%x, yc => this%coef%dof%y, &
257 zc => this%coef%dof%z, &
258 nx => this%coef%nx, ny => this%coef%ny, nz => this%coef%nz, &
259 lx => this%coef%Xh%lx)
262 facet = this%facet(i)
266 call this%eval(x(i), y(i), z(i), &
267 xc(idx(1), idx(2), idx(3), idx(4)), &
268 yc(idx(1), idx(2), idx(3), idx(4)), &
269 zc(idx(1), idx(2), idx(3), idx(4)), &
270 nx(idx(2), idx(3), facet, idx(4)), &
271 ny(idx(2), idx(3), facet, idx(4)), &
272 nz(idx(2), idx(3), facet, idx(4)), &
273 idx(1), idx(2), idx(3), idx(4), &
276 call this%eval(x(i), y(i), z(i), &
277 xc(idx(1), idx(2), idx(3), idx(4)), &
278 yc(idx(1), idx(2), idx(3), idx(4)), &
279 zc(idx(1), idx(2), idx(3), idx(4)), &
280 nx(idx(1), idx(3), facet, idx(4)), &
281 ny(idx(1), idx(3), facet, idx(4)), &
282 nz(idx(1), idx(3), facet, idx(4)), &
283 idx(1), idx(2), idx(3), idx(4), &
286 call this%eval(x(i), y(i), z(i), &
287 xc(idx(1), idx(2), idx(3), idx(4)), &
288 yc(idx(1), idx(2), idx(3), idx(4)), &
289 zc(idx(1), idx(2), idx(3), idx(4)), &
290 nx(idx(1), idx(2), facet, idx(4)), &
291 ny(idx(1), idx(2), facet, idx(4)), &
292 nz(idx(1), idx(2), facet, idx(4)), &
293 idx(1), idx(2), idx(3), idx(4), &
299 call device_memcpy(x, usr_x_d, m, host_to_device, sync = .false.)
300 call device_memcpy(y, usr_y_d, m, host_to_device, sync = .false.)
301 call device_memcpy(z, usr_z_d, m, host_to_device, sync = .true.)
306 call device_inhom_dirichlet_apply_vector(this%msk_d, x_d, y_d, z_d, &
307 usr_x_d, usr_y_d, usr_z_d, m)
318 this%eval => usr_eval
327 if (.not.
associated(this%coef))
then
328 call neko_warning(
'Missing coefficients')
332 if (.not.
associated(this%eval))
then
333 call neko_warning(
'Missing eval function')
337 if (.not. valid)
then
338 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)
Defines a boundary condition.
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, public neko_warning(warning_msg)
Reports a warning to standard output.
Base type for a boundary condition.
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)