37 use,
intrinsic :: iso_c_binding, only : c_ptr
40 use json_module,
only : json_file
52 real(kind=
rp),
allocatable :: flux_(:)
53 real(kind=
rp),
private :: init_flux_
62 procedure, pass(this) :: init_from_components => &
70 generic :: set_flux => set_flux_scalar, set_flux_array
83 class(
neumann_t),
intent(inout),
target :: this
84 type(
coef_t),
intent(in) :: coef
85 type(json_file),
intent(inout) :: json
88 call this%init_base(coef)
92 this%init_flux_ = flux
99 class(
neumann_t),
intent(inout),
target :: this
100 type(
coef_t),
intent(in) :: coef
101 real(kind=
rp),
intent(in) :: flux
103 call this%init_base(coef)
104 this%init_flux_ = flux
111 integer,
intent(in) :: n
112 real(kind=
rp),
intent(inout),
dimension(n) :: x
113 real(kind=
rp),
intent(in),
optional :: t
114 integer,
intent(in),
optional :: tstep
115 logical,
intent(in),
optional :: strong
116 integer :: i, m, k, facet
119 logical :: strong_ = .true.
121 if (
present(strong)) strong_ = strong
124 if (.not. strong_)
then
127 facet = this%facet(i)
132 x(k) = x(k) + this%flux_(i)*this%coef%area(idx(2), idx(3), facet,&
135 x(k) = x(k) + this%flux_(i)*this%coef%area(idx(1), idx(3), facet,&
138 x(k) = x(k) + this%flux_(i)*this%coef%area(idx(1), idx(2), facet,&
149 integer,
intent(in) :: n
150 real(kind=
rp),
intent(inout),
dimension(n) :: x
151 real(kind=
rp),
intent(inout),
dimension(n) :: y
152 real(kind=
rp),
intent(inout),
dimension(n) :: z
153 real(kind=
rp),
intent(in),
optional :: t
154 integer,
intent(in),
optional :: tstep
155 logical,
intent(in),
optional :: strong
157 call neko_error(
"Neumann bc not implemented for vectors")
164 class(
neumann_t),
intent(inout),
target :: this
166 real(kind=
rp),
intent(in),
optional :: t
167 integer,
intent(in),
optional :: tstep
168 logical,
intent(in),
optional :: strong
170 call neko_error(
"Neumann bc not implemented on the device")
177 class(
neumann_t),
intent(inout),
target :: this
181 real(kind=
rp),
intent(in),
optional :: t
182 integer,
intent(in),
optional :: tstep
183 logical,
intent(in),
optional :: strong
185 call neko_error(
"Neumann bc not implemented on the device")
191 class(
neumann_t),
target,
intent(inout) :: this
193 call this%free_base()
199 class(
neumann_t),
target,
intent(inout) :: this
201 call this%finalize_base()
202 allocate(this%flux_(this%msk(0)))
204 call cfill(this%flux_, this%init_flux_, this%msk(0))
210 real(kind=
rp) :: flux(this%msk(0))
218 real(kind=
rp),
intent(in) :: flux
227 real(kind=
rp),
intent(in) :: flux(this%msk(0))
229 call copy(this%flux_, flux, this%msk(0))
__device__ void nonlinear_index(const int idx, const int lx, int *index)
Retrieves a parameter by name or throws an error.
Defines a boundary condition.
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 .
Defines a Neumann boundary condition.
subroutine neumann_set_flux_array(this, flux)
Set the flux using an array of values.
subroutine neumann_init_from_components(this, coef, flux)
Constructor from components.
subroutine neumann_init(this, coef, json)
Constructor.
subroutine neumann_apply_vector_dev(this, x_d, y_d, z_d, t, tstep, strong)
Boundary condition apply for a generic Neumann condition to vectors x, y and z (device version)
subroutine neumann_apply_scalar(this, x, n, t, tstep, 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, t, tstep, strong)
Boundary condition apply for a generic Neumann condition to a vector x (device version)
subroutine neumann_free(this)
Destructor.
subroutine neumann_finalize(this)
Finalize by setting the flux.
subroutine neumann_apply_vector(this, x, y, z, n, t, tstep, strong)
Boundary condition apply for a generic Neumann condition to vectors x, y and z.
integer, parameter, public rp
Global precision used in computations.
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.