47 type(c_ptr),
private :: usr_x_d = c_null_ptr
78 ix, iy, iz, ie, t, tstep)
80 real(kind=
rp),
intent(inout) :: s
81 real(kind=
rp),
intent(in) :: x
82 real(kind=
rp),
intent(in) :: y
83 real(kind=
rp),
intent(in) :: z
84 real(kind=
rp),
intent(in) :: nx
85 real(kind=
rp),
intent(in) :: ny
86 real(kind=
rp),
intent(in) :: nz
87 integer,
intent(in) :: ix
88 integer,
intent(in) :: iy
89 integer,
intent(in) :: iz
90 integer,
intent(in) :: ie
91 real(kind=
rp),
intent(in) :: t
92 integer,
intent(in) :: tstep
106 if (c_associated(this%usr_x_d))
then
119 integer,
intent(in) :: n
120 real(kind=
rp),
intent(inout),
dimension(n) :: x
121 real(kind=
rp),
intent(in),
optional :: t
122 integer,
intent(in),
optional :: tstep
123 integer :: i, m, k, idx(4), facet, tstep_
132 if (
present(tstep))
then
138 associate(xc => this%coef%dof%x, yc => this%coef%dof%y, &
139 zc => this%coef%dof%z, nx => this%coef%nx, ny => this%coef%ny, &
140 nz => this%coef%nz, lx => this%coef%Xh%lx)
144 facet = this%facet(i)
148 call this%eval(x(k), &
149 xc(idx(1), idx(2), idx(3), idx(4)), &
150 yc(idx(1), idx(2), idx(3), idx(4)), &
151 zc(idx(1), idx(2), idx(3), idx(4)), &
152 nx(idx(2), idx(3), facet, idx(4)), &
153 ny(idx(2), idx(3), facet, idx(4)), &
154 nz(idx(2), idx(3), facet, idx(4)), &
155 idx(1), idx(2), idx(3), idx(4), &
158 call this%eval(x(k), &
159 xc(idx(1), idx(2), idx(3), idx(4)), &
160 yc(idx(1), idx(2), idx(3), idx(4)), &
161 zc(idx(1), idx(2), idx(3), idx(4)), &
162 nx(idx(1), idx(3), facet, idx(4)), &
163 ny(idx(1), idx(3), facet, idx(4)), &
164 nz(idx(1), idx(3), facet, idx(4)), &
165 idx(1), idx(2), idx(3), idx(4), &
168 call this%eval(x(k), &
169 xc(idx(1), idx(2), idx(3), idx(4)), &
170 yc(idx(1), idx(2), idx(3), idx(4)), &
171 zc(idx(1), idx(2), idx(3), idx(4)), &
172 nx(idx(1), idx(2), facet, idx(4)), &
173 ny(idx(1), idx(2), facet, idx(4)), &
174 nz(idx(1), idx(2), facet, idx(4)), &
175 idx(1), idx(2), idx(3), idx(4), &
190 real(kind=rp),
intent(in),
optional :: t
191 integer,
intent(in),
optional :: tstep
192 integer :: i, m, k, idx(4), facet, tstep_
194 integer(c_size_t) :: s
195 real(kind=rp),
allocatable :: x(:)
204 if (
present(tstep))
then
210 associate(xc => this%coef%dof%x, yc => this%coef%dof%y, &
211 zc => this%coef%dof%z, nx => this%coef%nx, ny => this%coef%ny, &
212 nz => this%coef%nz, lx => this%coef%Xh%lx, usr_x_d => this%usr_x_d)
216 if (.not. c_associated(usr_x_d))
then
220 call device_alloc(this%usr_x_d, s)
224 facet = this%facet(i)
228 call this%eval(x(i), &
229 xc(idx(1), idx(2), idx(3), idx(4)), &
230 yc(idx(1), idx(2), idx(3), idx(4)), &
231 zc(idx(1), idx(2), idx(3), idx(4)), &
232 nx(idx(2), idx(3), facet, idx(4)), &
233 ny(idx(2), idx(3), facet, idx(4)), &
234 nz(idx(2), idx(3), facet, idx(4)), &
235 idx(1), idx(2), idx(3), idx(4), &
238 call this%eval(x(i), &
239 xc(idx(1), idx(2), idx(3), idx(4)), &
240 yc(idx(1), idx(2), idx(3), idx(4)), &
241 zc(idx(1), idx(2), idx(3), idx(4)), &
242 nx(idx(1), idx(3), facet, idx(4)), &
243 ny(idx(1), idx(3), facet, idx(4)), &
244 nz(idx(1), idx(3), facet, idx(4)), &
245 idx(1), idx(2), idx(3), idx(4), &
248 call this%eval(x(i), &
249 xc(idx(1), idx(2), idx(3), idx(4)), &
250 yc(idx(1), idx(2), idx(3), idx(4)), &
251 zc(idx(1), idx(2), idx(3), idx(4)), &
252 nx(idx(1), idx(2), facet, idx(4)), &
253 ny(idx(1), idx(2), facet, idx(4)), &
254 nz(idx(1), idx(2), facet, idx(4)), &
255 idx(1), idx(2), idx(3), idx(4), &
260 call device_memcpy(x, this%usr_x_d, m, host_to_device, sync=.true.)
265 call device_inhom_dirichlet_apply_scalar(this%msk_d, x_d, &
275 integer,
intent(in) :: n
276 real(kind=rp),
intent(inout),
dimension(n) :: x
277 real(kind=rp),
intent(inout),
dimension(n) :: y
278 real(kind=rp),
intent(inout),
dimension(n) :: z
279 real(kind=rp),
intent(in),
optional :: t
280 integer,
intent(in),
optional :: tstep
290 real(kind=rp),
intent(in),
optional :: t
291 integer,
intent(in),
optional :: tstep
300 this%eval => user_scalar_bc
309 if (.not.
associated(this%coef))
then
310 call neko_warning(
'Missing coefficients')
314 if (.not.
associated(this%eval))
then
315 call neko_warning(
'Missing eval function')
319 if (.not. valid)
then
320 call neko_error(
'Invalid user defined scalar condition')
__device__ void nonlinear_index(const int idx, const int lx, int *index)
Abstract interface defining a user defined scalar boundary condition (pointwise) Just imitating inflo...
Defines a boundary condition.
Device abstraction, common interface for various accelerators.
subroutine, public device_free(x_d)
Deallocate memory on the device.
integer, parameter, public rp
Global precision used in computations.
Defines dirichlet conditions for scalars.
subroutine usr_scalar_apply_scalar_dev(this, x_d, t, tstep)
Scalar apply (device version) Just imitating inflow for now, but we should look this over Applies bou...
subroutine usr_scalar_apply_vector_dev(this, x_d, y_d, z_d, t, tstep)
No-op vector apply (device version)
subroutine usr_scalar_apply_vector(this, x, y, z, n, t, tstep)
No-op vector apply.
subroutine usr_scalar_set_eval(this, user_scalar_bc)
Assign user provided eval function.
subroutine usr_scalar_apply_scalar(this, x, n, t, tstep)
Scalar apply Just imitating inflow for now, but we should look this over Applies boundary conditions ...
subroutine usr_scalar_validate(this)
Validate user scalar condition.
subroutine usr_scalar_free(this)
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,...
User defined dirichlet condition for scalars.