Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
wall_model_bc.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, intrinsic :: iso_c_binding, only : c_ptr
38 use utils, only : neko_error
39 use json_utils, only : json_get
40 use coefs, only : coef_t
41 use wall_model, only : wall_model_t, wall_model_allocator
43 use json_module, only : json_file
44 use time_state, only : time_state_t
45 implicit none
46 private
47
51 type, public, extends(shear_stress_t) :: wall_model_bc_t
53 class(wall_model_t), allocatable :: wall_model
54 contains
56 procedure, pass(this) :: init => wall_model_bc_init
58 procedure, pass(this) :: free => wall_model_bc_free
60 procedure, pass(this) :: finalize => wall_model_bc_finalize
61 procedure, pass(this) :: apply_scalar => wall_model_bc_apply_scalar
62 procedure, pass(this) :: apply_vector => wall_model_bc_apply_vector
63 procedure, pass(this) :: apply_scalar_dev => &
65 procedure, pass(this) :: apply_vector_dev => &
67 end type wall_model_bc_t
68
69contains
70
72 subroutine wall_model_bc_apply_scalar(this, x, n, time, strong)
73 class(wall_model_bc_t), intent(inout) :: this
74 integer, intent(in) :: n
75 real(kind=rp), intent(inout), dimension(n) :: x
76 type(time_state_t), intent(in), optional :: time
77 logical, intent(in), optional :: strong
78
79 call neko_error("The wall model bc is not applicable to scalar fields.")
80
81 end subroutine wall_model_bc_apply_scalar
82
90 subroutine wall_model_bc_apply_vector(this, x, y, z, n, time, strong)
91 class(wall_model_bc_t), intent(inout) :: this
92 integer, intent(in) :: n
93 real(kind=rp), intent(inout), dimension(n) :: x
94 real(kind=rp), intent(inout), dimension(n) :: y
95 real(kind=rp), intent(inout), dimension(n) :: z
96 type(time_state_t), intent(in), optional :: time
97 logical, intent(in), optional :: strong
98 integer :: i, m, k, fid
99 real(kind=rp) :: magtau
100 logical :: strong_
101
102 if (present(strong)) then
103 strong_ = strong
104 else
105 strong_ = .true.
106 end if
107
108 if (present(time) .eqv. .false.) then
109 call neko_error("wall_model_bc_apply_vector: time is required.")
110 end if
111
112 if (.not. strong_) then
113 ! Compute the wall stress using the wall model.
114 call this%wall_model%compute(time%t, time%tstep)
115
116 ! Populate the 3D wall stress field for post-processing.
117 call this%wall_model%compute_mag_field()
118
119 ! Set the computed stress for application by the underlying Neumann
120 ! boundary conditions.
121 call this%set_stress(this%wall_model%tau_x, this%wall_model%tau_y, &
122 this%wall_model%tau_z)
123 end if
124
125 ! Either add the stress to the RHS or apply the non-penetration condition
126 call this%shear_stress_t%apply_vector(x, y, z, n, time, strong_)
127
128 end subroutine wall_model_bc_apply_vector
129
132 subroutine wall_model_bc_apply_scalar_dev(this, x_d, time, strong, strm)
133 class(wall_model_bc_t), intent(inout), target :: this
134 type(c_ptr), intent(inout) :: x_d
135 type(time_state_t), intent(in), optional :: time
136 logical, intent(in), optional :: strong
137 type(c_ptr), intent(inout) :: strm
138
139 call neko_error("The wall model bc is not applicable to scalar fields.")
140
141 end subroutine wall_model_bc_apply_scalar_dev
142
145 subroutine wall_model_bc_apply_vector_dev(this, x_d, y_d, z_d, time, &
146 strong, strm)
147 class(wall_model_bc_t), intent(inout), target :: this
148 type(c_ptr), intent(inout) :: x_d
149 type(c_ptr), intent(inout) :: y_d
150 type(c_ptr), intent(inout) :: z_d
151 type(time_state_t), intent(in), optional :: time
152 logical, intent(in), optional :: strong
153 logical :: strong_
154 type(c_ptr), intent(inout) :: strm
155
156 if (present(strong)) then
157 strong_ = strong
158 else
159 strong_ = .true.
160 end if
161
162 if (present(time) .eqv. .false.) then
163 call neko_error("wall_model_bc_apply_vector_dev: time is required.")
164 end if
165
166 if (.not. strong_) then
167 ! Compute the wall stress using the wall model.
168 call this%wall_model%compute(time%t, time%tstep)
169
170 ! Populate the 3D wall stress field for post-processing.
171 call this%wall_model%compute_mag_field()
172
173 ! Set the computed stress for application by the underlying Neumann
174 ! boundary conditions.
175 call this%set_stress(this%wall_model%tau_x, this%wall_model%tau_y, &
176 this%wall_model%tau_z)
177 end if
178
179 ! Either add the stress to the RHS or apply the non-penetration condition
180 call this%shear_stress_t%apply_vector_dev(x_d, y_d, z_d, &
181 time, strong_, strm)
182
183 end subroutine wall_model_bc_apply_vector_dev
184
188 subroutine wall_model_bc_init(this, coef, json)
189 class(wall_model_bc_t), target, intent(inout) :: this
190 type(coef_t), target, intent(in) :: coef
191 type(json_file), intent(inout) :: json
192 real(kind=rp) :: value(3) = [0, 0, 0]
193 character(len=:), allocatable :: type_name
194
195 ! Initialize the shear stress base class.
196 call this%shear_stress_t%init_from_components(coef, value)
197 ! Partial initialization of the wall model by parsing the JSON.
198 call json_get(json, "model", type_name)
199 call wall_model_allocator(this%wall_model, type_name)
200
201 call this%wall_model%partial_init(coef, json)
202 end subroutine wall_model_bc_init
203
205 subroutine wall_model_bc_free(this)
206 class(wall_model_bc_t), target, intent(inout) :: this
207
208 call this%shear_stress_t%free()
209 call this%wall_model%free()
210
211 end subroutine wall_model_bc_free
212
214 subroutine wall_model_bc_finalize(this, only_facets)
215 class(wall_model_bc_t), target, intent(inout) :: this
216 logical, optional, intent(in) :: only_facets
217
218 if (present(only_facets)) then
219 if (only_facets .eqv. .false.) then
220 call neko_error("For wall_model_bc_t, only_facets has to be true.")
221 end if
222 end if
223
224 call this%shear_stress_t%finalize(.true.)
225 call this%wall_model%finalize(this%msk, this%facet)
226 end subroutine wall_model_bc_finalize
227
228end module wall_model_bc
Retrieves a parameter by name or throws an error.
Coefficients.
Definition coef.f90:34
Utilities for retrieving parameters from the case files.
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.
Module with things related to the simulation time.
Utilities.
Definition utils.f90:35
Defines the wall_model_bc_t type. Maintainer: Timofey Mukha.
subroutine wall_model_bc_init(this, coef, json)
Constructor.
subroutine wall_model_bc_apply_scalar(this, x, n, time, strong)
Apply shear stress for a scalar field x.
subroutine wall_model_bc_apply_scalar_dev(this, x_d, time, strong, strm)
Boundary condition apply for a generic wall_model_bc condition to a vector x (device version)
subroutine wall_model_bc_apply_vector_dev(this, x_d, y_d, z_d, time, strong, strm)
Boundary condition apply for a generic wall_model_bc condition to vectors x, y and z (device version)
subroutine wall_model_bc_apply_vector(this, x, y, z, n, time, strong)
Apply the boundary condition to the right-hand side.
subroutine wall_model_bc_free(this)
Destructor.
subroutine wall_model_bc_finalize(this, only_facets)
Finalize by building mask arrays and init'ing the wall model.
Implements wall_model_t.
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
Definition coef.f90:55
A shear stress boundary condition.
A struct that contains all info about the time, expand as needed.
Base abstract type for wall-stress models for wall-modelled LES.
A shear stress boundary condition, computing the stress values using a wall model.