Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
dynamic_smagorinsky.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!
35 use num_types, only : rp
36 use field, only : field_t
38 use les_model, only : les_model_t
40 use json_module, only : json_file
41 use utils, only : neko_error
45 use logger, only : log_size, neko_log
48 implicit none
49 private
50
53 type, public, extends(les_model_t) :: dynamic_smagorinsky_t
55 type(field_t) :: c_dyn
57 type(elementwise_filter_t) :: test_filter
61 type(field_t) :: mij(6)
63 type(field_t) :: lij(6)
65 type(field_t) :: num
67 type(field_t) :: den
68 contains
70 procedure, pass(this) :: init => dynamic_smagorinsky_init
72 procedure, pass(this) :: free => dynamic_smagorinsky_free
74 procedure, pass(this) :: compute => dynamic_smagorinsky_compute
75
77
78contains
82 subroutine dynamic_smagorinsky_init(this, fluid, json)
83 class(dynamic_smagorinsky_t), intent(inout) :: this
84 class(fluid_scheme_base_t), intent(inout), target :: fluid
85 type(json_file), intent(inout) :: json
86 type(json_file) :: json_subdict
87 character(len=:), allocatable :: nut_name
88 integer :: i
89 character(len=:), allocatable :: delta_type
90 logical :: if_ext
91 character(len=:), allocatable :: filter_type
92 character(len=LOG_SIZE) :: log_buf
93
94 associate(dofmap => fluid%dm_Xh, &
95 coef => fluid%c_Xh)
96
97 call json_get_or_default(json, "nut_field", nut_name, "nut")
98 call json_get_or_default(json, "delta_type", delta_type, "pointwise")
99 call json_get_or_default(json, "extrapolation", if_ext, .false.)
100
101 call this%free()
102 call this%init_base(fluid, nut_name, delta_type, if_ext)
103 call json_extract_object(json, "test_filter", json_subdict)
104 call this%test_filter%init(json_subdict, coef)
105 if (json%valid_path('filter.transfer_function')) then
106 call neko_error("Dynamic Smagorinsky model does not support " // &
107 "transfer function specified in the json file. " // &
108 "Please hard-code it in " // &
109 "subroutine set_ds_filt() in " // &
110 "src/les/dynamic_smagorisnky.f90")
111 end if
112 if (json%valid_path('filter.type')) then
113 call json_get(json, "filter.type", filter_type)
114 if (trim(filter_type) .ne. "elementwise") then
115 call neko_error("Currently only elementwise filter is supported &
116 &for dynamic smagorinsky model.")
117 end if
118 end if
119
120 call set_ds_filt(this%test_filter)
121
122 call neko_log%section('LES model')
123 write(log_buf, '(A)') 'Model : Dynamic Smagorinsky'
124 call neko_log%message(log_buf)
125 write(log_buf, '(A, A)') 'Delta evaluation : ', delta_type
126 call neko_log%message(log_buf)
127 write(log_buf, '(A, A)') 'Test filter type : ', &
128 this%test_filter%filter_type
129 call neko_log%message(log_buf)
130 write(log_buf, '(A, L1)') 'extrapolation : ', if_ext
131 call neko_log%message(log_buf)
132 call neko_log%end_section()
133
134 call this%c_dyn%init(dofmap, "ds_c_dyn")
135 call this%num%init(dofmap, "ds_num")
136 call this%den%init(dofmap, "ds_den")
137
138 do i = 1, 6
139 call this%mij(i)%init(dofmap)
140 call this%lij(i)%init(dofmap)
141 end do
142
143 end associate
144
145 end subroutine dynamic_smagorinsky_init
146
149 class(dynamic_smagorinsky_t), intent(inout) :: this
150 integer :: i
151
152 call this%c_dyn%free()
153 do i = 1, 6
154 call this%mij(i)%free()
155 call this%lij(i)%free()
156 end do
157 call this%num%free()
158 call this%den%free()
159 call this%test_filter%free()
160 call this%free_base()
161
162 end subroutine dynamic_smagorinsky_free
163
167 subroutine dynamic_smagorinsky_compute(this, t, tstep)
168 class(dynamic_smagorinsky_t), intent(inout) :: this
169 real(kind=rp), intent(in) :: t
170 integer, intent(in) :: tstep
171
172 type(field_t), pointer :: u, v, w, u_e, v_e, w_e
173
174 if (this%if_ext .eqv. .true.) then
175 ! Extrapolate the velocity fields
176 associate(ulag => this%ulag, vlag => this%vlag, &
177 wlag => this%wlag, ext_bdf => this%ext_bdf)
178
179 u => neko_field_registry%get_field_by_name("u")
180 v => neko_field_registry%get_field_by_name("v")
181 w => neko_field_registry%get_field_by_name("w")
182 u_e => neko_field_registry%get_field_by_name("u_e")
183 v_e => neko_field_registry%get_field_by_name("v_e")
184 w_e => neko_field_registry%get_field_by_name("w_e")
185
186 call this%sumab%compute_fluid(u_e, v_e, w_e, u, v, w, &
187 ulag, vlag, wlag, ext_bdf%advection_coeffs, ext_bdf%nadv)
188
189 end associate
190 end if
191
192 ! Compute the eddy viscosity field
193 if (neko_bcknd_device .eq. 1) then
194 call dynamic_smagorinsky_compute_device(this%if_ext, t, tstep, &
195 this%coef, this%nut, &
196 this%delta, this%c_dyn, this%test_filter, &
197 this%mij, this%lij, this%num, this%den)
198 else
199 call dynamic_smagorinsky_compute_cpu(this%if_ext, t, tstep, &
200 this%coef, this%nut, &
201 this%delta, this%c_dyn, this%test_filter, &
202 this%mij, this%lij, this%num, this%den)
203 end if
204
205 end subroutine dynamic_smagorinsky_compute
206
208 subroutine set_ds_filt(filter_1d)
209 type(elementwise_filter_t), intent(inout) :: filter_1d
210 integer :: i
211
212 if (filter_1d%nx .le. 2) then
213 call neko_error("Dynamic Smagorinsky model error: test filter is not &
214 &defined for the current polynomial order")
215 end if
216 if (mod(filter_1d%nx,2) .eq. 0) then ! number of grid spacing is odd
217 ! cutoff at polynomial order int((filter_1d%nx)/2)
218 filter_1d%transfer(int((filter_1d%nx)/2)) = 0.95_rp
219 filter_1d%transfer(int((filter_1d%nx)/2) + 1) = 0.50_rp
220 filter_1d%transfer(int((filter_1d%nx)/2) + 2) = 0.05_rp
221 if ((int((filter_1d%nx)/2) + 2) .lt. filter_1d%nx) then
222 do i = int((filter_1d%nx)/2) + 3, filter_1d%nx
223 filter_1d%transfer(i) = 0.0_rp
224 end do
225 end if
226 ! make delta_ratio = (nx-1)/(nt-1) as close to 2
227 filter_1d%nt = int(filter_1d%nx/2) + 1
228 else ! number of grid spacing is even
229 ! cutoff at polynomial order int((filter_1d%nx-1)/2)
230 filter_1d%transfer(int((filter_1d%nx-1)/2)) = 0.95_rp
231 filter_1d%transfer(int((filter_1d%nx-1)/2) + 1) = 0.50_rp
232 filter_1d%transfer(int((filter_1d%nx-1)/2) + 2) = 0.05_rp
233 if ((int((filter_1d%nx-1)/2) + 2) .lt. filter_1d%nx) then
234 do i = int((filter_1d%nx-1)/2) + 3, filter_1d%nx
235 filter_1d%transfer(i) = 0.0_rp
236 end do
237 end if
238 ! make delta_ratio = (nx-1)/(nt-1) = 2
239 filter_1d%nt = int((filter_1d%nx-1)/2) + 1
240 end if
241
242 call filter_1d%build_1d()
243
244 end subroutine set_ds_filt
245
246end module dynamic_smagorinsky
Retrieves a parameter by name or assigns a provided default value. In the latter case also adds the m...
Retrieves a parameter by name or throws an error.
Defines a mapping of the degrees of freedom.
Definition dofmap.f90:35
Implements the CPU kernel for the smagorinsky_t type.
subroutine, public dynamic_smagorinsky_compute_cpu(if_ext, t, tstep, coef, nut, delta, c_dyn, test_filter, mij, lij, num, den)
Compute eddy viscosity on the CPU.
Implements the device kernel for the smagorinsky_t type.
subroutine, public dynamic_smagorinsky_compute_device(if_ext, t, tstep, coef, nut, delta, c_dyn, test_filter, mij, lij, num, den)
Compute eddy viscosity on the device.
Implements dynamic_smagorinsky_t.
subroutine set_ds_filt(filter_1d)
Set up the test filter.
subroutine dynamic_smagorinsky_free(this)
Destructor for the les_model_t (base) class.
subroutine dynamic_smagorinsky_compute(this, t, tstep)
Compute eddy viscosity.
subroutine dynamic_smagorinsky_init(this, fluid, json)
Constructor.
Implements elementwise_filter_t.
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.
subroutine, public json_extract_object(json, name, object)
Extract object as a separate JSON dictionary.
Implements les_model_t.
Definition les_model.f90:35
Logging routines.
Definition log.f90:34
type(log_t), public neko_log
Global log stream.
Definition log.f90:70
integer, parameter, public log_size
Definition log.f90:46
Build configurations.
integer, parameter neko_bcknd_device
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
Utilities.
Definition utils.f90:35
Implements the dynamic Smagorinsky LES model.
Implements the elementwise filter for SEM.
Base type of all fluid formulations.
Base abstract type for LES models based on the Boussinesq approximation.
Definition les_model.f90:64