Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
dynamic_smagorinsky_device.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_list, only : field_list_t
39 use field, only : field_t
40 use operators, only : strain_rate
41 use coefs, only : coef_t
43 use gs_ops, only : gs_op_add
44 use device_math, only : device_col2
48 implicit none
49 private
50
52
53contains
54
68 subroutine dynamic_smagorinsky_compute_device(if_ext, t, tstep, coef, nut, &
69 delta, c_dyn, test_filter, mij, lij, num, den)
70 logical, intent(in) :: if_ext
71 real(kind=rp), intent(in) :: t
72 integer, intent(in) :: tstep
73 type(coef_t), intent(in) :: coef
74 type(field_t), intent(inout) :: nut
75 type(field_t), intent(in) :: delta
76 type(field_t), intent(inout) :: c_dyn
77 type(elementwise_filter_t), intent(inout) :: test_filter
78 type(field_t), intent(inout) :: mij(6), lij(6)
79 type(field_t), intent(inout) :: num, den
80
81 type(field_t), pointer :: u, v, w
82 type(field_t) :: c_dyn_curr
83 ! the strain rate tensor
84 type(field_t), pointer :: s11, s22, s33, s12, s13, s23, s_abs
85 real(kind=rp) :: alpha ! running averaging coefficient
86 integer :: temp_indices(7)
87 integer :: i
88
89 if (tstep .eq. 1) then
90 alpha = 1.0_rp
91 else
92 alpha = 0.9_rp
93 end if
94
95 if (if_ext .eqv. .true.) then
96 u => neko_field_registry%get_field_by_name("u_e")
97 v => neko_field_registry%get_field_by_name("v_e")
98 w => neko_field_registry%get_field_by_name("w_e")
99 else
100 u => neko_field_registry%get_field_by_name("u")
101 v => neko_field_registry%get_field_by_name("v")
102 w => neko_field_registry%get_field_by_name("w")
103 end if
104
105 call neko_scratch_registry%request_field(s11, temp_indices(1))
106 call neko_scratch_registry%request_field(s22, temp_indices(2))
107 call neko_scratch_registry%request_field(s33, temp_indices(3))
108 call neko_scratch_registry%request_field(s12, temp_indices(4))
109 call neko_scratch_registry%request_field(s13, temp_indices(5))
110 call neko_scratch_registry%request_field(s23, temp_indices(6))
111 call neko_scratch_registry%request_field(s_abs, temp_indices(7))
112
113 ! Compute the strain rate tensor
114 call strain_rate(s11%x, s22%x, s33%x, s12%x, s13%x, s23%x, u, v, w, coef)
115
116 call coef%gs_h%op(s11%x, s11%dof%size(), gs_op_add)
117 call coef%gs_h%op(s22%x, s11%dof%size(), gs_op_add)
118 call coef%gs_h%op(s33%x, s11%dof%size(), gs_op_add)
119 call coef%gs_h%op(s12%x, s11%dof%size(), gs_op_add)
120 call coef%gs_h%op(s13%x, s11%dof%size(), gs_op_add)
121 call coef%gs_h%op(s23%x, s11%dof%size(), gs_op_add)
122
123 call device_col2(s11%x_d, coef%mult_d, s11%dof%size())
124 call device_col2(s22%x_d, coef%mult_d, s22%dof%size())
125 call device_col2(s33%x_d, coef%mult_d, s33%dof%size())
126 call device_col2(s12%x_d, coef%mult_d, s12%dof%size())
127 call device_col2(s13%x_d, coef%mult_d, s13%dof%size())
128 call device_col2(s23%x_d, coef%mult_d, s23%dof%size())
129
130 call device_s_abs_compute(s_abs%x_d, s11%x_d, s22%x_d, s33%x_d, &
131 s12%x_d, s13%x_d, s23%x_d, &
132 s11%dof%size())
133
134 call compute_lij_device(lij, u, v, w, test_filter, u%dof%size())
135 call compute_nut_device(nut, c_dyn, num, den, lij, mij, &
136 s11, s22, s33, s12, s13, s23, &
137 s_abs, test_filter, delta, alpha, &
138 coef, u%dof%size())
139
140 call coef%gs_h%op(nut, gs_op_add)
141 call device_col2(nut%x_d, coef%mult_d, nut%dof%size())
142
143 call neko_scratch_registry%relinquish_field(temp_indices)
145
154 subroutine compute_lij_device(lij, u, v, w, test_filter, n)
155 type(field_t), intent(inout) :: lij(6)
156 type(field_t), pointer, intent(in) :: u, v, w
157 type(elementwise_filter_t), intent(inout) :: test_filter
158 integer, intent(in) :: n
159 integer :: i
161 type(field_t), pointer :: fu, fv, fw, fuu, fvv, fww, fuv, fuw, fvw
162 integer :: temp_indices(9)
163
164 call neko_scratch_registry%request_field(fu, temp_indices(1))
165 call neko_scratch_registry%request_field(fv, temp_indices(2))
166 call neko_scratch_registry%request_field(fw, temp_indices(3))
167 call neko_scratch_registry%request_field(fuu, temp_indices(4))
168 call neko_scratch_registry%request_field(fvv, temp_indices(5))
169 call neko_scratch_registry%request_field(fww, temp_indices(6))
170 call neko_scratch_registry%request_field(fuv, temp_indices(7))
171 call neko_scratch_registry%request_field(fuw, temp_indices(8))
172 call neko_scratch_registry%request_field(fvw, temp_indices(9))
173
174 ! Use test filter for the velocity fields
175 call test_filter%apply(fu, u)
176 call test_filter%apply(fv, v)
177 call test_filter%apply(fw, w)
178
179 !! ___ ___
180 !! Compute u_i*u_j and u_i*u_j
181 call device_lij_compute_part1(lij(1)%x_d, lij(2)%x_d, lij(3)%x_d, &
182 lij(4)%x_d, lij(5)%x_d, lij(6)%x_d, &
183 u%x_d, v%x_d, w%x_d, &
184 fu%x_d, fv%x_d, fw%x_d, &
185 fuu%x_d, fvv%x_d, fww%x_d, &
186 fuv%x_d, fuw%x_d, fvw%x_d, n)
187
188 !! Filter u_i*u_j by the test filter
189 call test_filter%apply(fuu, fuu)
190 call test_filter%apply(fvv, fvv)
191 call test_filter%apply(fww, fww)
192 call test_filter%apply(fuv, fuv)
193 call test_filter%apply(fuw, fuw)
194 call test_filter%apply(fvw, fvw)
195
196 !! Assember the final form
197 !! ___ ___ _______
198 !! L_ij = u_i*u_j - u_i*u_j
199 call device_lij_compute_part2(lij(1)%x_d, lij(2)%x_d, lij(3)%x_d, &
200 lij(4)%x_d, lij(5)%x_d, lij(6)%x_d, &
201 fuu%x_d, fvv%x_d, fww%x_d, &
202 fuv%x_d, fuw%x_d, fvw%x_d, n)
203 call neko_scratch_registry%relinquish_field(temp_indices)
204 end subroutine compute_lij_device
205
226 subroutine compute_nut_device(nut, c_dyn, num, den, lij, mij, &
227 s11, s22, s33, s12, s13, s23, &
228 s_abs, test_filter, delta, alpha, &
229 coef, n)
230 type(field_t), intent(inout) :: nut, c_dyn
231 type(field_t), intent(inout) :: num, den
232 type(field_t), intent(in) :: lij(6)
233 type(field_t), intent(inout) :: mij(6)
234 type(field_t), intent(inout) :: s11, s22, s33, s12, s13, s23, s_abs
235 type(elementwise_filter_t), intent(inout) :: test_filter
236 type(field_t), intent(in) :: delta
237 real(kind=rp), intent(in) :: alpha
238 type(coef_t), intent(in) :: coef
239 integer, intent(in) :: n
240
241 real(kind=rp) :: delta_ratio2
242 integer :: temp_indices(13)
243 type(field_t), pointer :: fs11, fs22, fs33, fs12, fs13, fs23, fs_abs, &
244 fsabss11, fsabss22, fsabss33, &
245 fsabss12, fsabss13, fsabss23
246
247 delta_ratio2 = ((test_filter%nx-1.0_rp)/(test_filter%nt-1.0_rp))**2
248
249 call neko_scratch_registry%request_field(fs11, temp_indices(1))
250 call neko_scratch_registry%request_field(fs22, temp_indices(2))
251 call neko_scratch_registry%request_field(fs33, temp_indices(3))
252 call neko_scratch_registry%request_field(fs12, temp_indices(4))
253 call neko_scratch_registry%request_field(fs13, temp_indices(5))
254 call neko_scratch_registry%request_field(fs23, temp_indices(6))
255 call neko_scratch_registry%request_field(fsabss11, temp_indices(7))
256 call neko_scratch_registry%request_field(fsabss22, temp_indices(8))
257 call neko_scratch_registry%request_field(fsabss33, temp_indices(9))
258 call neko_scratch_registry%request_field(fsabss12, temp_indices(10))
259 call neko_scratch_registry%request_field(fsabss13, temp_indices(11))
260 call neko_scratch_registry%request_field(fsabss23, temp_indices(12))
261 call neko_scratch_registry%request_field(fs_abs, temp_indices(13))
262
263 !! Compute M_ij
264 !! _____ ____
265 !! Compute s_abs and s_ij
266 call test_filter%apply(fs_abs, s_abs)
267 call test_filter%apply(fs11, s11)
268 call test_filter%apply(fs22, s22)
269 call test_filter%apply(fs33, s33)
270 call test_filter%apply(fs12, s12)
271 call test_filter%apply(fs13, s13)
272 call test_filter%apply(fs23, s23)
273
274 !! _____ ____
275 !! Compute (delta_test/delta)^2 s_abs*s_ij and s_abs*s_ij
276 call device_mij_compute_part1(mij(1)%x_d, mij(2)%x_d, mij(3)%x_d, &
277 mij(4)%x_d, mij(5)%x_d, mij(6)%x_d, &
278 s_abs%x_d, s11%x_d, s22%x_d, s33%x_d, &
279 s12%x_d, s13%x_d, s23%x_d, &
280 fs_abs%x_d, fs11%x_d, fs22%x_d, fs33%x_d, &
281 fs12%x_d, fs13%x_d, fs23%x_d, &
282 fsabss11%x_d, fsabss22%x_d, fsabss33%x_d, &
283 fsabss12%x_d, fsabss13%x_d, fsabss23%x_d, &
284 delta_ratio2, n)
285
286 !! Filter s_abs*s_ij by the test filter
287 call test_filter%apply(fsabss11, fsabss11)
288 call test_filter%apply(fsabss22, fsabss22)
289 call test_filter%apply(fsabss33, fsabss33)
290 call test_filter%apply(fsabss12, fsabss12)
291 call test_filter%apply(fsabss13, fsabss13)
292 call test_filter%apply(fsabss23, fsabss23)
293
294 !! Finalise the compute of Mij and nut
295 call device_mij_nut_compute_part2(mij(1)%x_d, mij(2)%x_d, mij(3)%x_d, &
296 mij(4)%x_d, mij(5)%x_d, mij(6)%x_d, &
297 lij(1)%x_d, lij(2)%x_d, lij(3)%x_d, &
298 lij(4)%x_d, lij(5)%x_d, lij(6)%x_d, &
299 fsabss11%x_d, fsabss22%x_d, fsabss33%x_d, &
300 fsabss12%x_d, fsabss13%x_d, fsabss23%x_d, &
301 num%x_d, den%x_d, c_dyn%x_d, delta%x_d, &
302 s_abs%x_d, nut%x_d, alpha, n)
303 call neko_scratch_registry%relinquish_field(temp_indices)
304 end subroutine compute_nut_device
305
307
Coefficients.
Definition coef.f90:34
subroutine, public device_lij_compute_part1(l11_d, l22_d, l33_d, l12_d, l13_d, l23_d, u_d, v_d, w_d, fu_d, fv_d, fw_d, fuu_d, fvv_d, fww_d, fuv_d, fuw_d, fvw_d, n)
part 1 of the computing of the lij field
subroutine, public device_mij_nut_compute_part2(m11_d, m22_d, m33_d, m12_d, m13_d, m23_d, l11_d, l22_d, l33_d, l12_d, l13_d, l23_d, fsabss11_d, fsabss22_d, fsabss33_d, fsabss12_d, fsabss13_d, fsabss23_d, num_d, den_d, c_dyn_d, delta_d, s_abs_d, nut_d, alpha, n)
part 1 of the computing of the mij field
subroutine, public device_mij_compute_part1(m11_d, m22_d, m33_d, m12_d, m13_d, m23_d, s_abs_d, s11_d, s22_d, s33_d, s12_d, s13_d, s23_d, fs_abs_d, fs11_d, fs22_d, fs33_d, fs12_d, fs13_d, fs23_d, fsabss11_d, fsabss22_d, fsabss33_d, fsabss12_d, fsabss13_d, fsabss23_d, delta_ratio2, n)
part 1 of the computing of the mij field
subroutine, public device_s_abs_compute(s_abs_d, s11_d, s22_d, s33_d, s12_d, s13_d, s23_d, n)
Compute the s_abs field for the Sigma model indevice.
subroutine, public device_lij_compute_part2(l11_d, l22_d, l33_d, l12_d, l13_d, l23_d, fuu_d, fvv_d, fww_d, fuv_d, fuw_d, fvw_d, n)
part 2 of the computing of the lij field
subroutine, public device_col2(a_d, b_d, n, strm)
Vector multiplication .
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.
subroutine compute_nut_device(nut, c_dyn, num, den, lij, mij, s11, s22, s33, s12, s13, s23, s_abs, test_filter, delta, alpha, coef, n)
Compute M_ij and nut on the device.
subroutine compute_lij_device(lij, u, v, w, test_filter, n)
Compute Germano Identity on the device.
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
Defines Gather-scatter operations.
Definition gs_ops.f90:34
integer, parameter, public gs_op_add
Definition gs_ops.f90:36
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
Operators.
Definition operators.f90:34
subroutine, public strain_rate(s11, s22, s33, s12, s13, s23, u, v, w, coef)
Compute the strain rate tensor, i.e 0.5 * du_i/dx_j + du_j/dx_i.
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.
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
Definition coef.f90:55
Implements the elementwise filter for SEM.
field_list_t, To be able to group fields together