37 use,
intrinsic :: iso_c_binding, only : c_ptr, c_null_ptr
40 use json_module,
only : json_file
58 real(kind=
rp),
private :: init_flux_
59 logical :: uniform_0 = .false.
68 procedure, pass(this) :: init_from_components => &
76 generic :: set_flux => set_flux_scalar, set_flux_array
89 class(
neumann_t),
intent(inout),
target :: this
90 type(
coef_t),
target,
intent(in) :: coef
91 type(json_file),
intent(inout) :: json
94 call this%init_base(coef)
98 this%init_flux_ = flux
105 class(
neumann_t),
intent(inout),
target :: this
106 type(
coef_t),
intent(in) :: coef
107 real(kind=
rp),
intent(in) :: flux
109 call this%init_base(coef)
110 this%init_flux_ = flux
117 integer,
intent(in) :: n
118 real(kind=
rp),
intent(inout),
dimension(n) :: x
120 logical,
intent(in),
optional :: strong
121 integer :: i, m, k, facet
126 if (
present(strong))
then
133 if (.not. strong_)
then
136 facet = this%facet(i)
142 this%flux_%x(i)*this%coef%area(idx(2), idx(3), facet, idx(4))
145 this%flux_%x(i)*this%coef%area(idx(1), idx(3), facet, idx(4))
148 this%flux_%x(i)*this%coef%area(idx(1), idx(2), facet, idx(4))
158 integer,
intent(in) :: n
159 real(kind=
rp),
intent(inout),
dimension(n) :: x
160 real(kind=
rp),
intent(inout),
dimension(n) :: y
161 real(kind=
rp),
intent(inout),
dimension(n) :: z
163 logical,
intent(in),
optional :: strong
165 if (.not. this%uniform_0 .and. this%msk(0) .gt. 0)
then
166 call neko_error(
"Neumann bc not implemented for vectors")
173 class(
neumann_t),
intent(inout),
target :: this
174 type(c_ptr),
intent(inout) :: x_d
176 logical,
intent(in),
optional :: strong
177 type(c_ptr),
intent(inout) :: strm
180 if (
present(strong))
then
186 if (.not. this%uniform_0 .and. this%msk(0) .gt. 0 .and. &
189 this%flux_%x_d, this%coef%area_d, this%coef%Xh%lx, &
190 size(this%msk), strm)
198 class(
neumann_t),
intent(inout),
target :: this
199 type(c_ptr),
intent(inout) :: x_d
200 type(c_ptr),
intent(inout) :: y_d
201 type(c_ptr),
intent(inout) :: z_d
203 logical,
intent(in),
optional :: strong
204 type(c_ptr),
intent(inout) :: strm
206 if (.not. this%uniform_0 .and. this%msk(0) .gt. 0)
then
207 call neko_error(
"Neumann bc not implemented for vectors.")
214 class(
neumann_t),
target,
intent(inout) :: this
216 call this%free_base()
222 class(
neumann_t),
target,
intent(inout) :: this
223 logical,
optional,
intent(in) :: only_facets
226 if (
present(only_facets))
then
227 if (only_facets .eqv. .false.)
then
228 call neko_error(
"For neumann_t, only_facets has to be true.")
232 call this%finalize_base(.true.)
234 call this%flux_%init(this%msk(0))
236 this%flux_ = this%init_flux_
237 this%uniform_0 = .true.
240 this%uniform_0 =
abscmp(this%init_flux_, 0.0_rp) .and. this%uniform_0
248 real(kind=
rp) :: flux(this%msk(0))
257 real(kind=
rp),
intent(in) :: flux
260 this%uniform_0 =
abscmp(flux, 0.0_rp)
277 this%uniform_0 = .false.
__inline__ __device__ void nonlinear_index(const int idx, const int lx, int *index)
Copy data between host and device (or device and device)
Retrieves a parameter by name or throws an error.
Defines a boundary condition.
subroutine, public device_copy(a_d, b_d, n, strm)
Copy a vector .
subroutine, public device_cfill(a_d, c, n, strm)
Set all elements to a constant c .
subroutine, public device_neumann_apply_scalar(msk, facet, x, flux, area, lx, m, strm)
Device abstraction, common interface for various accelerators.
integer, parameter, public device_to_host
Utilities for retrieving parameters from the case files.
subroutine, public cfill(a, c, n)
Set all elements to a constant c .
subroutine, public copy(a, b, n)
Copy a vector .
integer, parameter neko_bcknd_device
Defines a Neumann boundary condition.
subroutine neumann_set_flux_array(this, flux)
Set the flux using an array of values.
subroutine neumann_apply_vector_dev(this, x_d, y_d, z_d, time, strong, strm)
Boundary condition apply for a generic Neumann condition to vectors x, y and z (device version)
subroutine neumann_init_from_components(this, coef, flux)
Constructor from components.
subroutine neumann_finalize(this, only_facets)
Finalize by setting the flux.
subroutine neumann_init(this, coef, json)
Constructor.
subroutine neumann_apply_vector(this, x, y, z, n, time, strong)
Boundary condition apply for a generic Neumann condition to vectors x, y and z.
subroutine neumann_apply_scalar(this, x, n, time, strong)
Boundary condition apply for a generic Neumann condition to a vector x.
subroutine neumann_set_flux_scalar(this, flux)
Set the flux using a scalar.
pure real(kind=rp) function, dimension(this%msk(0)) neumann_flux(this)
Get the flux.
subroutine neumann_apply_scalar_dev(this, x_d, time, strong, strm)
Boundary condition apply for a generic Neumann condition to a vector x (device version)
subroutine neumann_free(this)
Destructor.
integer, parameter, public rp
Global precision used in computations.
Module with things related to the simulation time.
Base type for a boundary condition.
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
A Neumann boundary condition for scalar fields. Sets the flux of the field to the chosen value.
A struct that contains all info about the time, expand as needed.