Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
spalding.f90
Go to the documentation of this file.
1! Copyright (c) 2024, 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!
33!
36 use field, only: field_t
37 use num_types, only : rp
38 use json_module, only : json_file
39 use coefs, only : coef_t
41 use wall_model, only : wall_model_t
46 use field_math, only: field_invcol3
47 use vector, only : vector_t
48 use math, only: masked_gather_copy_0
51
52 implicit none
53 private
54
57 type, public, extends(wall_model_t) :: spalding_t
59 real(kind=rp) :: kappa = 0.41_rp
61 real(kind=rp) :: b = 5.2_rp
63 type(vector_t) :: nu
64 contains
66 procedure, pass(this) :: init => spalding_init
69 procedure, pass(this) :: partial_init => spalding_partial_init
71 procedure, pass(this) :: finalize => spalding_finalize
73 procedure, pass(this) :: init_from_components => &
76 procedure, pass(this) :: free => spalding_free
78 procedure, pass(this) :: compute_nu => spalding_compute_nu
80 procedure, pass(this) :: compute => spalding_compute
81 end type spalding_t
82
83contains
91 subroutine spalding_init(this, scheme_name, coef, msk, facet, h_index, json)
92 class(spalding_t), intent(inout) :: this
93 character(len=*), intent(in) :: scheme_name
94 type(coef_t), intent(in) :: coef
95 integer, intent(in) :: msk(:)
96 integer, intent(in) :: facet(:)
97 integer, intent(in) :: h_index
98 type(json_file), intent(inout) :: json
99 real(kind=rp) :: kappa, b
100
101 call json_get_or_default(json, "kappa", kappa, 0.41_rp)
102 call json_get_or_default(json, "B", b, 5.2_rp)
103
104 call this%init_from_components(scheme_name, coef, msk, facet, h_index, &
105 kappa, b)
106 end subroutine spalding_init
107
111 subroutine spalding_partial_init(this, coef, json)
112 class(spalding_t), intent(inout) :: this
113 type(coef_t), intent(in) :: coef
114 type(json_file), intent(inout) :: json
115
116 call this%partial_init_base(coef, json)
117 call json_get_or_default(json, "kappa", this%kappa, 0.41_rp)
118 call json_get_or_default(json, "B", this%B, 5.2_rp)
119
120 end subroutine spalding_partial_init
121
125 subroutine spalding_finalize(this, msk, facet)
126 class(spalding_t), intent(inout) :: this
127 integer, intent(in) :: msk(:)
128 integer, intent(in) :: facet(:)
129
130 call this%finalize_base(msk, facet)
131 call this%nu%init(this%n_nodes)
132 end subroutine spalding_finalize
133
142 subroutine spalding_init_from_components(this, scheme_name, coef, msk, &
143 facet, h_index, kappa, B)
144 class(spalding_t), intent(inout) :: this
145 character(len=*), intent(in) :: scheme_name
146 type(coef_t), intent(in) :: coef
147 integer, intent(in) :: msk(:)
148 integer, intent(in) :: facet(:)
149 integer, intent(in) :: h_index
150 real(kind=rp), intent(in) :: kappa
151 real(kind=rp), intent(in) :: b
152
153 call this%free()
154 call this%init_base(scheme_name, coef, msk, facet, h_index)
155
156 this%kappa = kappa
157 this%B = b
158
159 call this%nu%init(this%n_nodes)
160 end subroutine spalding_init_from_components
161
163 subroutine spalding_compute_nu(this)
164 class(spalding_t), intent(inout) :: this
165 type(field_t), pointer :: temp
166 integer :: idx
167
168 call neko_scratch_registry%request_field(temp, idx)
169 call field_invcol3(temp, this%mu, this%rho)
170
171 if (neko_bcknd_device .eq. 1) then
172 call device_masked_gather_copy_0(this%nu%x_d, temp%x_d, this%msk_d, &
173 temp%size(), this%nu%size())
174 else
175 call masked_gather_copy_0(this%nu%x, temp%x, this%msk, temp%size(), &
176 this%nu%size())
177 end if
178
179 call neko_scratch_registry%relinquish_field(idx)
180 end subroutine spalding_compute_nu
181
183 subroutine spalding_free(this)
184 class(spalding_t), intent(inout) :: this
185
186 call this%free_base()
187
188 end subroutine spalding_free
189
193 subroutine spalding_compute(this, t, tstep)
194 class(spalding_t), intent(inout) :: this
195 real(kind=rp), intent(in) :: t
196 integer, intent(in) :: tstep
197 type(field_t), pointer :: u
198 type(field_t), pointer :: v
199 type(field_t), pointer :: w
200 integer :: i
201 real(kind=rp) :: ui, vi, wi, magu, utau, normu, guess
202
203 call this%compute_nu()
204
205 u => neko_field_registry%get_field("u")
206 v => neko_field_registry%get_field("v")
207 w => neko_field_registry%get_field("w")
208
209 if (neko_bcknd_device .eq. 1) then
210 call spalding_compute_device(u%x_d, v%x_d, w%x_d, this%ind_r_d, &
211 this%ind_s_d, this%ind_t_d, this%ind_e_d, &
212 this%n_x%x_d, this%n_y%x_d, this%n_z%x_d, &
213 this%nu%x_d, this%h%x_d, &
214 this%tau_x%x_d, this%tau_y%x_d, this%tau_z%x_d, &
215 this%n_nodes, u%Xh%lx, &
216 this%kappa, this%B, tstep)
217 else
218 call spalding_compute_cpu(u%x, v%x, w%x, &
219 this%ind_r, this%ind_s, this%ind_t, this%ind_e, &
220 this%n_x%x, this%n_y%x, this%n_z%x, &
221 this%nu%x, this%h%x, &
222 this%tau_x%x, this%tau_y%x, this%tau_z%x, &
223 this%n_nodes, u%Xh%lx, u%msh%nelv, &
224 this%kappa, this%B, tstep)
225 end if
226
227 end subroutine spalding_compute
228end module spalding
__global__ void spalding_compute(const T *__restrict__ u_d, const T *__restrict__ v_d, const T *__restrict__ w_d, const int *__restrict__ ind_r_d, const int *__restrict__ ind_s_d, const int *__restrict__ ind_t_d, const int *__restrict__ ind_e_d, const T *__restrict__ n_x_d, const T *__restrict__ n_y_d, const T *__restrict__ n_z_d, const T *__restrict__ nu_d, const T *__restrict__ h_d, T *__restrict__ tau_x_d, T *__restrict__ tau_y_d, T *__restrict__ tau_z_d, const int n_nodes, const int lx, const T kappa, const T B, const int tstep)
Retrieves a parameter by name or assigns a provided default value. In the latter case also adds the m...
Coefficients.
Definition coef.f90:34
subroutine, public device_masked_gather_copy_0(a_d, b_d, mask_d, n, n_mask, strm)
Gather a masked vector .
subroutine, public field_invcol3(a, b, c, n)
Invert a vector .
Defines a registry for storing solution fields.
type(field_registry_t), target, public neko_field_registry
Global field registry.
Defines a field.
Definition field.f90:34
Utilities for retrieving parameters from the case files.
Definition math.f90:60
subroutine, public masked_gather_copy_0(a, b, mask, n, n_mask)
Gather a masked vector to reduced contigous vector .
Definition math.f90:318
Build configurations.
integer, parameter neko_bcknd_device
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
Defines a registry for storing and requesting temporary fields This can be used when you have a funct...
type(scratch_registry_t), target, public neko_scratch_registry
Global scratch registry.
Implements the CPU kernel for the spalding_t type.
subroutine, public spalding_compute_cpu(u, v, w, ind_r, ind_s, ind_t, ind_e, n_x, n_y, n_z, nu, h, tau_x, tau_y, tau_z, n_nodes, lx, nelv, kappa, b, tstep)
Compute the wall shear stress on cpu using Spalding's model.
Implements the device kernel for the spalding_t type.
subroutine, public spalding_compute_device(u_d, v_d, w_d, ind_r_d, ind_s_d, ind_t_d, ind_e_d, n_x_d, n_y_d, n_z_d, nu_d, h_d, tau_x_d, tau_y_d, tau_z_d, n_nodes, lx, kappa, b, tstep)
Compute the wall shear stress on device using Spalding's model.
Implements spalding_t.
Definition spalding.f90:35
subroutine spalding_init_from_components(this, scheme_name, coef, msk, facet, h_index, kappa, b)
Constructor from components.
Definition spalding.f90:144
subroutine spalding_partial_init(this, coef, json)
Constructor from JSON.
Definition spalding.f90:112
subroutine spalding_compute_nu(this)
Compute the kinematic viscosity vector.
Definition spalding.f90:164
subroutine spalding_free(this)
Destructor for the spalding_t (base) class.
Definition spalding.f90:184
subroutine spalding_finalize(this, msk, facet)
Finalize the construction using the mask and facet arrays of the bc.
Definition spalding.f90:126
subroutine spalding_init(this, scheme_name, coef, msk, facet, h_index, json)
Constructor from JSON.
Definition spalding.f90:92
Defines a vector.
Definition vector.f90:34
Implements wall_model_t.
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
Definition coef.f90:55
Wall model based on Spalding's law of the wall. Reference: http://dx.doi.org/10.1115/1....
Definition spalding.f90:57
Base abstract type for wall-stress models for wall-modelled LES.