Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
field_dirichlet.f90
Go to the documentation of this file.
1! Copyright (c) 2020-2025, The Neko Authors
2! All rights reserved.
3!
4! Redistribution and use in source and binary forms, with or without
5! modification, are permitted provided that the following conditions
6! are met:
7!
8! * Redistributions of source code must retain the above copyright
9! notice, this list of conditions and the following disclaimer.
10!
11! * Redistributions in binary form must reproduce the above
12! copyright notice, this list of conditions and the following
13! disclaimer in the documentation and/or other materials provided
14! with the distribution.
15!
16! * Neither the name of the authors nor the names of its
17! contributors may be used to endorse or promote products derived
18! from this software without specific prior written permission.
19!
20! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
21! "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
22! LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
23! FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
24! COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
25! INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
26! BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
27! LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
28! CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
29! LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
30! ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
31! POSSIBILITY OF SUCH DAMAGE.
32!
35 use num_types, only: rp
36 use coefs, only: coef_t
37 use dirichlet, only: dirichlet_t
38 use bc, only: bc_t
39 use bc_list, only : bc_list_t
40 use utils, only: split_string
41 use field, only : field_t
42 use field_list, only : field_list_t
43 use math, only: masked_copy_0
45 use dofmap, only : dofmap_t
46 use utils, only: neko_error
47 use json_module, only : json_file
48 use field_list, only : field_list_t
49 use json_utils, only : json_get
50 use, intrinsic :: iso_c_binding, only : c_ptr, c_size_t
51 use time_state, only : time_state_t
52 implicit none
53 private
54
65 type, public, extends(bc_t) :: field_dirichlet_t
68 type(field_t) :: field_bc
73 procedure(field_dirichlet_update), nopass, pointer :: update => null()
74 contains
76 procedure, pass(this) :: init => field_dirichlet_init
78 procedure, pass(this) :: init_from_components => &
81 procedure, pass(this) :: free => field_dirichlet_free
83 procedure, pass(this) :: finalize => field_dirichlet_finalize
85 procedure, pass(this) :: apply_scalar => field_dirichlet_apply_scalar
87 procedure, pass(this) :: apply_vector => field_dirichlet_apply_vector
89 procedure, pass(this) :: apply_vector_dev => &
92 procedure, pass(this) :: apply_scalar_dev => &
94
95 end type field_dirichlet_t
96
105 abstract interface
106 subroutine field_dirichlet_update(fields, bc, time)
108 type(field_list_t), intent(inout) :: fields
109 type(field_dirichlet_t), intent(in) :: bc
110 type(time_state_t), intent(in) :: time
111 end subroutine field_dirichlet_update
112 end interface
113
114 public :: field_dirichlet_update
115
116contains
120 subroutine field_dirichlet_init(this, coef, json)
121 class(field_dirichlet_t), intent(inout), target :: this
122 type(coef_t), target, intent(in) :: coef
123 type(json_file), intent(inout) ::json
124 character(len=:), allocatable :: field_name
125
126 call json_get(json, "field_name", field_name)
127 call this%init_from_components(coef, field_name)
128
129 end subroutine field_dirichlet_init
130
133 subroutine field_dirichlet_init_from_components(this, coef, field_name)
134 class(field_dirichlet_t), intent(inout), target :: this
135 type(coef_t), intent(in) :: coef
136 character(len=*), intent(in) :: field_name
137
138 call this%init_base(coef)
139
140 call this%field_bc%init(this%dof, field_name)
141 call this%field_list%init(1)
142 call this%field_list%assign_to_field(1, this%field_bc)
144
147 subroutine field_dirichlet_free(this)
148 class(field_dirichlet_t), target, intent(inout) :: this
149
150 call this%free_base
151 call this%field_bc%free()
152 call this%field_list%free()
153
154 if (associated(this%update)) then
155 this%update => null()
156 end if
157
158 end subroutine field_dirichlet_free
159
164 subroutine field_dirichlet_apply_scalar(this, x, n, time, strong)
165 class(field_dirichlet_t), intent(inout) :: this
166 integer, intent(in) :: n
167 real(kind=rp), intent(inout), dimension(n) :: x
168 type(time_state_t), intent(in), optional :: time
169 logical, intent(in), optional :: strong
170 logical :: strong_
171
172 if (present(strong)) then
173 strong_ = strong
174 else
175 strong_ = .true.
176 end if
177
178 if (strong_) then
179
180 if (.not. this%updated) then
181 call this%update(this%field_list, this, time)
182 this%updated = .true.
183 end if
184
185 call masked_copy_0(x, this%field_bc%x, this%msk, n, this%msk(0))
186 end if
187
188 end subroutine field_dirichlet_apply_scalar
189
194 subroutine field_dirichlet_apply_scalar_dev(this, x_d, time, strong, strm)
195 class(field_dirichlet_t), intent(inout), target :: this
196 type(c_ptr), intent(inout) :: x_d
197 type(time_state_t), intent(in), optional :: time
198 logical, intent(in), optional :: strong
199 type(c_ptr), intent(inout) :: strm
200 logical :: strong_
201
202 if (present(strong)) then
203 strong_ = strong
204 else
205 strong_ = .true.
206 end if
207
208 if (strong_) then
209 if (.not. this%updated) then
210 call this%update(this%field_list, this, time)
211 this%updated = .true.
212 end if
213
214 if (this%msk(0) .gt. 0) then
215 call device_masked_copy_0(x_d, this%field_bc%x_d, this%msk_d, &
216 this%field_bc%dof%size(), this%msk(0), strm)
217 end if
218 end if
219
221
228 subroutine field_dirichlet_apply_vector(this, x, y, z, n, time, strong)
229 class(field_dirichlet_t), intent(inout) :: this
230 integer, intent(in) :: n
231 real(kind=rp), intent(inout), dimension(n) :: x
232 real(kind=rp), intent(inout), dimension(n) :: y
233 real(kind=rp), intent(inout), dimension(n) :: z
234 type(time_state_t), intent(in), optional :: time
235 logical, intent(in), optional :: strong
236
237 call neko_error("field_dirichlet cannot apply vector BCs.&
238 & Use field_dirichlet_vector instead!")
239
240 end subroutine field_dirichlet_apply_vector
241
248 subroutine field_dirichlet_apply_vector_dev(this, x_d, y_d, z_d, time, &
249 strong, strm)
250 class(field_dirichlet_t), intent(inout), target :: this
251 type(c_ptr), intent(inout) :: x_d
252 type(c_ptr), intent(inout) :: y_d
253 type(c_ptr), intent(inout) :: z_d
254 type(time_state_t), intent(in), optional :: time
255 logical, intent(in), optional :: strong
256 type(c_ptr), intent(inout) :: strm
257 call neko_error("field_dirichlet cannot apply vector BCs.&
258 & Use field_dirichlet_vector instead!")
259
261
263 subroutine field_dirichlet_finalize(this, only_facets)
264 class(field_dirichlet_t), target, intent(inout) :: this
265 logical, optional, intent(in) :: only_facets
266 logical :: only_facets_
267
268 if (present(only_facets)) then
269 only_facets_ = only_facets
270 else
271 only_facets_ = .false.
272 end if
273
274 call this%finalize_base(only_facets_)
275 end subroutine field_dirichlet_finalize
276end module field_dirichlet
Abstract interface defining a dirichlet condition on a list of fields.
Retrieves a parameter by name or throws an error.
Defines a list of bc_t.
Definition bc_list.f90:34
Defines a boundary condition.
Definition bc.f90:34
Coefficients.
Definition coef.f90:34
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.
Definition dirichlet.f90:34
Defines a mapping of the degrees of freedom.
Definition dofmap.f90:35
Defines user dirichlet condition for a scalar field.
subroutine field_dirichlet_finalize(this, only_facets)
Finalize.
subroutine field_dirichlet_apply_vector_dev(this, x_d, y_d, z_d, time, strong, strm)
(No-op) Apply vector (device).
subroutine field_dirichlet_apply_vector(this, x, y, z, n, time, strong)
(No-op) Apply vector.
subroutine field_dirichlet_init_from_components(this, coef, field_name)
Constructor from components.
subroutine field_dirichlet_apply_scalar(this, x, n, time, strong)
Apply scalar by performing a masked copy.
subroutine field_dirichlet_init(this, coef, json)
Constructor.
subroutine field_dirichlet_free(this)
Destructor. Currently thisfield_bc is being freed in fluid_scheme_incompressible::free
subroutine field_dirichlet_apply_scalar_dev(this, x_d, time, strong, strm)
Apply scalar (device).
Defines a field.
Definition field.f90:34
Utilities for retrieving parameters from the case files.
Definition math.f90:60
subroutine, public masked_copy_0(a, b, mask, n, n_mask)
Copy a masked vector .
Definition math.f90:274
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
Module with things related to the simulation time.
Utilities.
Definition utils.f90:35
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,...
Definition utils.f90:137
Base type for a boundary condition.
Definition bc.f90:61
A list of allocatable `bc_t`. Follows the standard interface of lists.
Definition bc_list.f90:48
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
Definition coef.f90:55
Generic Dirichlet boundary condition on .
Definition dirichlet.f90:48
User defined dirichlet condition, for which the user can work with an entire field....
field_list_t, To be able to group fields together
A struct that contains all info about the time, expand as needed.