48 use json_module,
only : json_file
50 use,
intrinsic :: iso_c_binding, only : c_ptr, c_size_t
73 procedure, pass(this) :: init_from_components => &
80 procedure, pass(this) :: apply_scalar => &
83 procedure, pass(this) :: apply_vector => &
86 procedure, pass(this) :: apply_vector_dev => &
89 procedure, pass(this) :: apply_scalar_dev => &
100 type(
coef_t),
target,
intent(in) :: coef
101 type(json_file),
intent(inout) ::json
103 call this%init_from_components(coef)
111 type(
coef_t),
intent(in) :: coef
113 call this%init_base(coef)
115 call this%bc_u%init_from_components(coef,
"u")
116 call this%bc_v%init_from_components(coef,
"v")
117 call this%bc_w%init_from_components(coef,
"w")
121 call this%field_list%init(3)
122 call this%field_list%assign_to_field(1, this%bc_u%field_bc)
123 call this%field_list%assign_to_field(2, this%bc_v%field_bc)
124 call this%field_list%assign_to_field(3, this%bc_w%field_bc)
133 call this%bc_u%free()
134 call this%bc_v%free()
135 call this%bc_w%free()
137 call this%field_list%free()
139 if (
associated(this%update))
then
150 integer,
intent(in) :: n
151 real(kind=
rp),
intent(inout),
dimension(n) :: x
153 logical,
intent(in),
optional :: strong
155 call neko_error(
"field_dirichlet_vector cannot apply scalar BCs.&
156 & Use field_dirichlet instead!")
166 type(c_ptr),
intent(inout) :: x_d
168 logical,
intent(in),
optional :: strong
169 type(c_ptr),
intent(inout) :: strm
171 call neko_error(
"field_dirichlet_vector cannot apply scalar BCs.&
172 & Use field_dirichlet instead!")
185 integer,
intent(in) :: n
186 real(kind=
rp),
intent(inout),
dimension(n) :: x
187 real(kind=
rp),
intent(inout),
dimension(n) :: y
188 real(kind=
rp),
intent(inout),
dimension(n) :: z
190 logical,
intent(in),
optional :: strong
193 if (
present(strong))
then
203 if (.not. this%updated)
then
204 call this%update(this%field_list, this%bc_u, time)
205 this%updated = .true.
208 call masked_copy_0(x, this%bc_u%field_bc%x, this%msk, n, this%msk(0))
209 call masked_copy_0(y, this%bc_v%field_bc%x, this%msk, n, this%msk(0))
210 call masked_copy_0(z, this%bc_w%field_bc%x, this%msk, n, this%msk(0))
224 type(c_ptr),
intent(inout) :: x_d
225 type(c_ptr),
intent(inout) :: y_d
226 type(c_ptr),
intent(inout) :: z_d
228 logical,
intent(in),
optional :: strong
229 type(c_ptr),
intent(inout) :: strm
232 if (
present(strong))
then
239 if (.not. this%updated)
then
240 call this%update(this%field_list, this%bc_u, time)
241 this%updated = .true.
244 if (this%msk(0) .gt. 0)
then
246 this%bc_u%dof%size(), this%msk(0), strm)
248 this%bc_v%dof%size(), this%msk(0), strm)
250 this%bc_w%dof%size(), this%msk(0), strm)
259 logical,
optional,
intent(in) :: only_facets
260 logical :: only_facets_
262 if (
present(only_facets))
then
263 only_facets_ = only_facets
265 only_facets_ = .false.
268 call this%finalize_base(only_facets_)
270 call this%bc_u%mark_facets(this%marked_facet)
271 call this%bc_v%mark_facets(this%marked_facet)
272 call this%bc_w%mark_facets(this%marked_facet)
274 call this%bc_u%finalize(only_facets_)
275 call this%bc_v%finalize(only_facets_)
276 call this%bc_w%finalize(only_facets_)
Abstract interface defining a dirichlet condition on a list of fields.
Defines a boundary condition.
subroutine, public device_masked_copy_0(a_d, b_d, mask_d, n, n_mask, strm)
Copy a masked vector .
Defines a dirichlet boundary condition.
Defines a mapping of the degrees of freedom.
Defines inflow dirichlet conditions.
subroutine field_dirichlet_vector_apply_vector_dev(this, x_d, y_d, z_d, time, strong, strm)
Apply the boundary condition to a vector field on the device.
subroutine field_dirichlet_vector_finalize(this, only_facets)
Finalize by building the mask arrays and propagating to underlying bcs.
subroutine field_dirichlet_vector_init_from_components(this, coef)
Constructor from components.
subroutine field_dirichlet_vector_apply_scalar_dev(this, x_d, time, strong, strm)
No-op apply scalar (device).
subroutine field_dirichlet_vector_init(this, coef, json)
Constructor.
subroutine field_dirichlet_vector_apply_vector(this, x, y, z, n, time, strong)
Apply the boundary condition to a vector field.
subroutine field_dirichlet_vector_apply_scalar(this, x, n, time, strong)
No-op apply scalar.
subroutine field_dirichlet_vector_free(this)
Destructor. Currently unused as is, all field_dirichlet attributes are freed in fluid_scheme_incompre...
Defines user dirichlet condition for a scalar field.
subroutine, public masked_copy_0(a, b, mask, n, n_mask)
Copy a masked vector .
integer, parameter, public rp
Global precision used in computations.
Module with things related to the simulation time.
character(len=100) function, dimension(:), allocatable, public split_string(string, delimiter)
Split a string based on delimiter (tokenizer) OBS: very hacky, this should really be improved,...
Base type for a boundary condition.
A list of allocatable `bc_t`. Follows the standard interface of lists.
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
Generic Dirichlet boundary condition on .
User defined dirichlet condition, for which the user can work with an entire field....
Extension of the user defined dirichlet condition field_dirichlet
field_list_t, To be able to group fields together
A struct that contains all info about the time, expand as needed.