Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
shear_stress.f90
Go to the documentation of this file.
1! Copyright (c) 2024-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!
36 use num_types, only : rp
37 use bc, only : bc_t
38 use, intrinsic :: iso_c_binding, only : c_ptr
39 use utils, only : neko_error
40 use coefs, only : coef_t
41 use symmetry, only : symmetry_t
42 use neumann, only : neumann_t
43 use json_module, only : json_file
44 use json_utils, only : json_get
45 use vector, only : vector_t
46 use time_state, only : time_state_t
47 implicit none
48 private
49
52 type, public, extends(bc_t) :: shear_stress_t
53 ! This bc takes care of setting the wall-normal component to zero.
54 ! It can be passed to associated bc lists, which take care of masking
55 ! changes in residuals and solution increments.
57
59 type(neumann_t) :: neumann_x
61 type(neumann_t) :: neumann_y
63 type(neumann_t) :: neumann_z
64 contains
65 procedure, pass(this) :: apply_scalar => shear_stress_apply_scalar
66 procedure, pass(this) :: apply_vector => shear_stress_apply_vector
67 procedure, pass(this) :: apply_scalar_dev => shear_stress_apply_scalar_dev
68 procedure, pass(this) :: apply_vector_dev => shear_stress_apply_vector_dev
70 procedure, pass(this) :: init => shear_stress_init
72 procedure, pass(this) :: init_from_components => &
74 procedure, pass(this) :: set_stress_scalar => &
76 procedure, pass(this) :: set_stress_array => &
79 generic :: set_stress => set_stress_scalar, set_stress_array
81 procedure, pass(this) :: free => shear_stress_free
83 procedure, pass(this) :: finalize => shear_stress_finalize
84 end type shear_stress_t
85
86contains
87
89 subroutine shear_stress_apply_scalar(this, x, n, time, strong)
90 class(shear_stress_t), intent(inout) :: this
91 integer, intent(in) :: n
92 real(kind=rp), intent(inout), dimension(n) :: x
93 type(time_state_t), intent(in), optional :: time
94 logical, intent(in), optional :: strong
95 integer :: i, m, k, facet
96 ! Store non-linear index
97 integer :: idx(4)
98
99 call neko_error("The shear stress bc is not applicable to scalar fields.")
100
101 end subroutine shear_stress_apply_scalar
102
105 subroutine shear_stress_apply_vector(this, x, y, z, n, time, strong)
106 class(shear_stress_t), intent(inout) :: this
107 integer, intent(in) :: n
108 real(kind=rp), intent(inout), dimension(n) :: x
109 real(kind=rp), intent(inout), dimension(n) :: y
110 real(kind=rp), intent(inout), dimension(n) :: z
111 type(time_state_t), intent(in), optional :: time
112 logical, intent(in), optional :: strong
113 logical :: strong_
114
115 if (present(strong)) then
116 strong_ = strong
117 else
118 strong_ = .true.
119 end if
120
121 if (strong_) then
122 call this%symmetry%apply_vector(x, y, z, n, strong = .true.)
123 else
124 call this%neumann_x%apply_scalar(x, n, strong = .false.)
125 call this%neumann_y%apply_scalar(y, n, strong = .false.)
126 call this%neumann_z%apply_scalar(z, n, strong = .false.)
127 end if
128
129 end subroutine shear_stress_apply_vector
130
133 subroutine shear_stress_apply_scalar_dev(this, x_d, time, strong, strm)
134 class(shear_stress_t), intent(inout), target :: this
135 type(c_ptr), intent(inout) :: x_d
136 type(time_state_t), intent(in), optional :: time
137 logical, intent(in), optional :: strong
138 type(c_ptr), intent(inout) :: strm
139
140 call neko_error("The shear stress bc is not applicable to scalar fields.")
141
142 end subroutine shear_stress_apply_scalar_dev
143
146 subroutine shear_stress_apply_vector_dev(this, x_d, y_d, z_d, time, &
147 strong, strm)
148 class(shear_stress_t), intent(inout), target :: this
149 type(c_ptr), intent(inout) :: x_d
150 type(c_ptr), intent(inout) :: y_d
151 type(c_ptr), intent(inout) :: z_d
152 type(time_state_t), intent(in), optional :: time
153 logical, intent(in), optional :: strong
154 type(c_ptr), intent(inout) :: strm
155 logical :: strong_
156
157 if (present(strong)) then
158 strong_ = strong
159 else
160 strong_ = .true.
161 end if
162
163 if (strong_) then
164 call this%symmetry%apply_vector_dev(x_d, y_d, z_d, strong = .true., &
165 strm = strm)
166 else
167 call this%neumann_x%apply_scalar_dev(x_d, strong = .false., strm = strm)
168 call this%neumann_y%apply_scalar_dev(y_d, strong = .false., strm = strm)
169 call this%neumann_z%apply_scalar_dev(z_d, strong = .false., strm = strm)
170 end if
171
172 end subroutine shear_stress_apply_vector_dev
173
177 subroutine shear_stress_init(this, coef, json)
178 class(shear_stress_t), target, intent(inout) :: this
179 type(coef_t), target, intent(in) :: coef
180 type(json_file), intent(inout) ::json
181 real(kind=rp), allocatable :: value(:)
182
183 call json_get(json, 'value', value)
184
185 if (size(value) .ne. 3) then
186 call neko_error ("The shear stress vector provided for the shear stress &
187 & boundary condition should have 3 components.")
188 end if
189
190 call this%init_from_components(coef, value)
191 end subroutine shear_stress_init
192
196 subroutine shear_stress_init_from_components(this, coef, value)
197 class(shear_stress_t), target, intent(inout) :: this
198 type(coef_t), intent(in) :: coef
199 real(kind=rp), intent(in) :: value(3)
200
201 call this%init_base(coef)
202 this%strong = .false.
203
204 call this%symmetry%free()
205 call this%symmetry%init_from_components(this%coef)
206
207 call this%neumann_x%free()
208 call this%neumann_y%free()
209 call this%neumann_z%free()
210
211 call this%neumann_x%init_from_components(this%coef, value(1))
212 call this%neumann_y%init_from_components(this%coef, value(2))
213 call this%neumann_z%init_from_components(this%coef, value(3))
214
216
217 subroutine shear_stress_finalize(this, only_facets)
218 class(shear_stress_t), target, intent(inout) :: this
219 logical, optional, intent(in) :: only_facets
220 logical :: only_facets_
221
222 if (present(only_facets)) then
223 only_facets_ = only_facets
224 else
225 only_facets_ = .false.
226 end if
227
228 call this%finalize_base(only_facets_)
229
230 call this%symmetry%mark_facets(this%marked_facet)
231 call this%symmetry%finalize()
232
233
234 call this%neumann_x%mark_facets(this%marked_facet)
235 call this%neumann_y%mark_facets(this%marked_facet)
236 call this%neumann_z%mark_facets(this%marked_facet)
237
238 call this%neumann_x%finalize(only_facets_)
239 call this%neumann_y%finalize(only_facets_)
240 call this%neumann_z%finalize(only_facets_)
241
242 end subroutine shear_stress_finalize
243
245 subroutine shear_stress_set_stress_scalar(this, tau_x, tau_y, tau_z)
246 class(shear_stress_t), intent(inout) :: this
247 real(kind=rp), intent(in) :: tau_x
248 real(kind=rp), intent(in) :: tau_y
249 real(kind=rp), intent(in) :: tau_z
250
251 ! Calls finalize and allocates the flux arrays
252 call this%neumann_x%set_flux(tau_x)
253 call this%neumann_y%set_flux(tau_y)
254 call this%neumann_z%set_flux(tau_z)
255
256
257 end subroutine shear_stress_set_stress_scalar
258
263 subroutine shear_stress_set_stress_array(this, tau_x, tau_y, tau_z)
264 class(shear_stress_t), intent(inout) :: this
265 type(vector_t), intent(in) :: tau_x
266 type(vector_t), intent(in) :: tau_y
267 type(vector_t), intent(in) :: tau_z
268
269 call this%neumann_x%set_flux(tau_x)
270 call this%neumann_y%set_flux(tau_y)
271 call this%neumann_z%set_flux(tau_z)
272
273 end subroutine shear_stress_set_stress_array
274
276 subroutine shear_stress_free(this)
277 class(shear_stress_t), target, intent(inout) :: this
278 call this%free_base
279 call this%symmetry%free
280
281 call this%neumann_x%free
282 call this%neumann_y%free
283 call this%neumann_z%free
284
285 end subroutine shear_stress_free
286end module shear_stress
Retrieves a parameter by name or throws an error.
Defines a boundary condition.
Definition bc.f90:34
Coefficients.
Definition coef.f90:34
Utilities for retrieving parameters from the case files.
Defines a Neumann boundary condition.
Definition neumann.f90:34
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
Defines a shear stress boundary condition for a vector field. Maintainer: Timofey Mukha.
subroutine shear_stress_init_from_components(this, coef, value)
Constructor from components.
subroutine shear_stress_free(this)
Destructor.
subroutine shear_stress_apply_scalar(this, x, n, time, strong)
Apply shear stress for a scalar field x.
subroutine shear_stress_set_stress_array(this, tau_x, tau_y, tau_z)
Set the shear stress components.
subroutine shear_stress_apply_scalar_dev(this, x_d, time, strong, strm)
Boundary condition apply for a generic shear_stress condition to a vector x (device version)
subroutine shear_stress_set_stress_scalar(this, tau_x, tau_y, tau_z)
Set the value of the shear stress vector using 3 scalars.
subroutine shear_stress_finalize(this, only_facets)
subroutine shear_stress_init(this, coef, json)
Constructor.
subroutine shear_stress_apply_vector_dev(this, x_d, y_d, z_d, time, strong, strm)
Boundary condition apply for a generic shear_stress condition to vectors x, y and z (device version)
subroutine shear_stress_apply_vector(this, x, y, z, n, time, strong)
Boundary condition apply for a generic shear_stress condition to vectors x, y and z.
Mixed Dirichlet-Neumann axis aligned symmetry plane.
Definition symmetry.f90:34
Module with things related to the simulation time.
Utilities.
Definition utils.f90:35
Defines a vector.
Definition vector.f90:34
Base type for a boundary condition.
Definition bc.f90:61
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
Definition coef.f90:55
A Neumann boundary condition for scalar fields. Sets the flux of the field to the chosen value.
Definition neumann.f90:56
A shear stress boundary condition.
Mixed Dirichlet-Neumann symmetry plane condition.
Definition symmetry.f90:50
A struct that contains all info about the time, expand as needed.