Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
field_dirichlet_vector.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
47 use utils, only: neko_error
48 use json_module, only : json_file
49 use field_list, only : field_list_t
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
56 ! for the application on a vector field.
57 type, public, extends(bc_t) :: field_dirichlet_vector_t
58 ! The bc for the first compoent.
59 type(field_dirichlet_t) :: bc_u
60 ! The bc for the second compoent.
61 type(field_dirichlet_t) :: bc_v
62 ! The bc for the third compoent.
63 type(field_dirichlet_t) :: bc_w
68 procedure(field_dirichlet_update), nopass, pointer :: update => null()
69 contains
71 procedure, pass(this) :: init => field_dirichlet_vector_init
73 procedure, pass(this) :: init_from_components => &
76 procedure, pass(this) :: free => field_dirichlet_vector_free
78 procedure, pass(this) :: finalize => field_dirichlet_vector_finalize
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 => &
92
93contains
94
98 subroutine field_dirichlet_vector_init(this, coef, json)
99 class(field_dirichlet_vector_t), intent(inout), target :: this
100 type(coef_t), target, intent(in) :: coef
101 type(json_file), intent(inout) ::json
102
103 call this%init_from_components(coef)
104
105 end subroutine field_dirichlet_vector_init
106
110 class(field_dirichlet_vector_t), intent(inout), target :: this
111 type(coef_t), intent(in) :: coef
112
113 call this%init_base(coef)
114
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")
118
119 ! TODO set to u v w values
120
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)
125
127
131 class(field_dirichlet_vector_t), target, intent(inout) :: this
132
133 call this%bc_u%free()
134 call this%bc_v%free()
135 call this%bc_w%free()
136
137 call this%field_list%free()
138
139 if (associated(this%update)) then
140 nullify(this%update)
141 end if
142 end subroutine field_dirichlet_vector_free
143
148 subroutine field_dirichlet_vector_apply_scalar(this, x, n, time, strong)
149 class(field_dirichlet_vector_t), intent(inout) :: this
150 integer, intent(in) :: n
151 real(kind=rp), intent(inout), dimension(n) :: x
152 type(time_state_t), intent(in), optional :: time
153 logical, intent(in), optional :: strong
154
155 call neko_error("field_dirichlet_vector cannot apply scalar BCs.&
156 & Use field_dirichlet instead!")
157
159
163 subroutine field_dirichlet_vector_apply_scalar_dev(this, x_d, time, &
164 strong, strm)
165 class(field_dirichlet_vector_t), intent(inout), target :: this
166 type(c_ptr), intent(inout) :: x_d
167 type(time_state_t), intent(in), optional :: time
168 logical, intent(in), optional :: strong
169 type(c_ptr), intent(inout) :: strm
170
171 call neko_error("field_dirichlet_vector cannot apply scalar BCs.&
172 & Use field_dirichlet instead!")
173
175
182 subroutine field_dirichlet_vector_apply_vector(this, x, y, z, n, time, &
183 strong)
184 class(field_dirichlet_vector_t), intent(inout) :: this
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
189 type(time_state_t), intent(in), optional :: time
190 logical, intent(in), optional :: strong
191 logical :: strong_
192
193 if (present(strong)) then
194 strong_ = strong
195 else
196 strong_ = .true.
197 end if
198
199 if (strong_) then
200
201 ! We can send any of the 3 bcs we have as argument, since they are all
202 ! the same boundary.
203 if (.not. this%updated) then
204 call this%update(this%field_list, this%bc_u, time)
205 this%updated = .true.
206 end if
207
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))
211 end if
212
214
221 subroutine field_dirichlet_vector_apply_vector_dev(this, x_d, y_d, z_d, &
222 time, strong, strm)
223 class(field_dirichlet_vector_t), intent(inout), target :: this
224 type(c_ptr), intent(inout) :: x_d
225 type(c_ptr), intent(inout) :: y_d
226 type(c_ptr), intent(inout) :: z_d
227 type(time_state_t), intent(in), optional :: time
228 logical, intent(in), optional :: strong
229 type(c_ptr), intent(inout) :: strm
230 logical :: strong_
231
232 if (present(strong)) then
233 strong_ = strong
234 else
235 strong_ = .true.
236 end if
237
238 if (strong_) then
239 if (.not. this%updated) then
240 call this%update(this%field_list, this%bc_u, time)
241 this%updated = .true.
242 end if
243
244 if (this%msk(0) .gt. 0) then
245 call device_masked_copy_0(x_d, this%bc_u%field_bc%x_d, this%bc_u%msk_d,&
246 this%bc_u%dof%size(), this%msk(0), strm)
247 call device_masked_copy_0(y_d, this%bc_v%field_bc%x_d, this%bc_v%msk_d,&
248 this%bc_v%dof%size(), this%msk(0), strm)
249 call device_masked_copy_0(z_d, this%bc_w%field_bc%x_d, this%bc_w%msk_d,&
250 this%bc_w%dof%size(), this%msk(0), strm)
251 end if
252 end if
253
255
257 subroutine field_dirichlet_vector_finalize(this, only_facets)
258 class(field_dirichlet_vector_t), target, intent(inout) :: this
259 logical, optional, intent(in) :: only_facets
260 logical :: only_facets_
261
262 if (present(only_facets)) then
263 only_facets_ = only_facets
264 else
265 only_facets_ = .false.
266 end if
267
268 call this%finalize_base(only_facets_)
269
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)
273
274 call this%bc_u%finalize(only_facets_)
275 call this%bc_v%finalize(only_facets_)
276 call this%bc_w%finalize(only_facets_)
277
279
280end module field_dirichlet_vector
Abstract interface defining a dirichlet condition on a list of fields.
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 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.
Defines a field.
Definition field.f90:34
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....
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.