65 c_dyn, test_filter, mij, lij, num, den)
66 real(kind=
rp),
intent(in) :: t
67 integer,
intent(in) :: tstep
68 type(
coef_t),
intent(in) :: coef
69 type(
field_t),
intent(inout) :: nut
70 type(
field_t),
intent(in) :: delta
71 type(
field_t),
intent(inout) :: c_dyn
73 type(
field_t),
intent(inout) :: mij(6), lij(6)
74 type(
field_t),
intent(inout) :: num, den
76 type(
field_t),
pointer :: u, v, w
79 type(
field_t),
pointer :: s11, s22, s33, s12, s13, s23, s_abs
80 real(kind=
rp) :: alpha
81 integer :: temp_indices(7)
84 if (tstep .eq. 1)
then
103 call strain_rate(s11%x, s22%x, s33%x, s12%x, s13%x, s23%x, u, v, w, coef)
105 call coef%gs_h%op(s11%x, s11%dof%size(),
gs_op_add)
106 call coef%gs_h%op(s22%x, s11%dof%size(),
gs_op_add)
107 call coef%gs_h%op(s33%x, s11%dof%size(),
gs_op_add)
108 call coef%gs_h%op(s12%x, s11%dof%size(),
gs_op_add)
109 call coef%gs_h%op(s13%x, s11%dof%size(),
gs_op_add)
110 call coef%gs_h%op(s23%x, s11%dof%size(),
gs_op_add)
112 do concurrent(i = 1:s11%dof%size())
113 s11%x(i,1,1,1) = s11%x(i,1,1,1) * coef%mult(i,1,1,1)
114 s22%x(i,1,1,1) = s22%x(i,1,1,1) * coef%mult(i,1,1,1)
115 s33%x(i,1,1,1) = s33%x(i,1,1,1) * coef%mult(i,1,1,1)
116 s12%x(i,1,1,1) = s12%x(i,1,1,1) * coef%mult(i,1,1,1)
117 s13%x(i,1,1,1) = s13%x(i,1,1,1) * coef%mult(i,1,1,1)
118 s23%x(i,1,1,1) = s23%x(i,1,1,1) * coef%mult(i,1,1,1)
121 do concurrent(i = 1:u%dof%size())
122 s_abs%x(i,1,1,1) = sqrt(2.0_rp * (s11%x(i,1,1,1)*s11%x(i,1,1,1) + &
123 s22%x(i,1,1,1)*s22%x(i,1,1,1) + &
124 s33%x(i,1,1,1)*s33%x(i,1,1,1)) + &
125 4.0_rp * (s12%x(i,1,1,1)*s12%x(i,1,1,1) + &
126 s13%x(i,1,1,1)*s13%x(i,1,1,1) + &
127 s23%x(i,1,1,1)*s23%x(i,1,1,1)))
130 call compute_lij_cpu(lij, u, v, w, test_filter, u%dof%size(), u%msh%nelv)
132 s_abs, test_filter, delta, u%dof%size(), u%msh%nelv)
135 do concurrent(i =1:u%dof%size())
136 if (den%x(i,1,1,1) .gt. 0.0_rp)
then
137 c_dyn%x(i,1,1,1) = 0.5_rp * (num%x(i,1,1,1)/den%x(i,1,1,1))
139 c_dyn%x(i,1,1,1) = 0.0_rp
141 c_dyn%x(i,1,1,1) =
max(c_dyn%x(i,1,1,1),0.0_rp)
142 nut%x(i,1,1,1) = c_dyn%x(i,1,1,1) * delta%x(i,1,1,1)**2 * s_abs%x(i,1,1,1)
158 type(
field_t),
intent(inout) :: lij(6)
159 type(
field_t),
pointer,
intent(in) :: u, v, w
161 integer,
intent(in) :: n
162 integer,
intent(inout) :: nelv
165 real(kind=
rp),
dimension(u%dof%size()) :: fu, fv, fw
168 call test_filter%filter_3d(fu, u%x, nelv)
169 call test_filter%filter_3d(fv, v%x, nelv)
170 call test_filter%filter_3d(fw, w%x, nelv)
173 do concurrent(i = 1:n)
174 lij(1)%x(i,1,1,1) = fu(i) * fu(i)
175 lij(2)%x(i,1,1,1) = fv(i) * fv(i)
176 lij(3)%x(i,1,1,1) = fw(i) * fw(i)
177 lij(4)%x(i,1,1,1) = fu(i) * fv(i)
178 lij(5)%x(i,1,1,1) = fu(i) * fw(i)
179 lij(6)%x(i,1,1,1) = fv(i) * fw(i)
185 call col3(fu, u%x, u%x, n)
186 call test_filter%filter_3d(fv, fu, nelv)
187 call sub2(lij(1)%x, fv, n)
189 call col3(fu, v%x, v%x, n)
190 call test_filter%filter_3d(fv, fu, nelv)
191 call sub2(lij(2)%x, fv, n)
193 call col3(fu, w%x, w%x, n)
194 call test_filter%filter_3d(fv, fu, nelv)
195 call sub2(lij(3)%x, fv, n)
197 call col3(fu, u%x, v%x, n)
198 call test_filter%filter_3d(fv, fu, nelv)
199 call sub2(lij(4)%x, fv, n)
201 call col3(fu, u%x, w%x, n)
202 call test_filter%filter_3d(fv, fu, nelv)
203 call sub2(lij(5)%x, fv, n)
205 call col3(fu, v%x, w%x, n)
206 call test_filter%filter_3d(fv, fu, nelv)
207 call sub2(lij(6)%x, fv, n)
220 s_abs, test_filter, delta, n, nelv)
221 type(
field_t),
intent(inout) :: mij(6)
222 type(
field_t),
intent(inout) :: s11, s22, s33, s12, s13, s23, s_abs
224 type(
field_t),
intent(in) :: delta
225 integer,
intent(in) :: n
226 integer,
intent(inout) :: nelv
228 real(kind=
rp),
dimension(n) :: fs11, fs22, fs33, fs12, fs13, fs23, fs_abs
229 real(kind=
rp) :: delta_ratio2
231 real(kind=
rp) :: delta2
233 delta_ratio2 = ((test_filter%nx-1)/(test_filter%nt-1))**2
238 call test_filter%filter_3d(fs_abs, s_abs%x, nelv)
240 call test_filter%filter_3d(fs11, s11%x, nelv)
241 call col3(mij(1)%x, fs_abs, fs11, n)
242 call cmult(mij(1)%x, delta_ratio2, n)
244 call test_filter%filter_3d(fs22, s22%x, nelv)
245 call col3(mij(2)%x, fs_abs, fs11, n)
246 call cmult(mij(2)%x, delta_ratio2, n)
248 call test_filter%filter_3d(fs33, s33%x, nelv)
249 call col3(mij(3)%x, fs_abs, fs11, n)
250 call cmult(mij(3)%x, delta_ratio2, n)
252 call test_filter%filter_3d(fs12, s12%x, nelv)
253 call col3(mij(4)%x, fs_abs, fs11, n)
254 call cmult(mij(4)%x, delta_ratio2, n)
256 call test_filter%filter_3d(fs13, s13%x, nelv)
257 call col3(mij(5)%x, fs_abs, fs11, n)
258 call cmult(mij(5)%x, delta_ratio2, n)
260 call test_filter%filter_3d(fs23, s23%x, nelv)
261 call col3(mij(6)%x, fs_abs, fs11, n)
262 call cmult(mij(6)%x, delta_ratio2, n)
268 call col3(fs11, s_abs%x, s11%x, n)
269 call test_filter%filter_3d(fs22, fs11, nelv)
270 call sub2(mij(1)%x, fs22, n)
272 call col3(fs11, s_abs%x, s22%x, n)
273 call test_filter%filter_3d(fs22, fs11, nelv)
274 call sub2(mij(2)%x, fs22, n)
276 call col3(fs11, s_abs%x, s33%x, n)
277 call test_filter%filter_3d(fs22, fs11, nelv)
278 call sub2(mij(3)%x, fs22, n)
280 call col3(fs11, s_abs%x, s12%x, n)
281 call test_filter%filter_3d(fs22, fs11, nelv)
282 call sub2(mij(4)%x, fs22, n)
284 call col3(fs11, s_abs%x, s13%x, n)
285 call test_filter%filter_3d(fs22, fs11, nelv)
286 call sub2(mij(5)%x, fs22, n)
288 call col3(fs11, s_abs%x, s23%x, n)
289 call test_filter%filter_3d(fs22, fs11, nelv)
290 call sub2(mij(6)%x, fs22, n)
293 do concurrent(i = 1:n)
294 delta2 = delta%x(i,1,1,1)**2
295 mij(1)%x(i,1,1,1) = mij(1)%x(i,1,1,1) * delta2
296 mij(2)%x(i,1,1,1) = mij(2)%x(i,1,1,1) * delta2
297 mij(3)%x(i,1,1,1) = mij(3)%x(i,1,1,1) * delta2
298 mij(4)%x(i,1,1,1) = mij(4)%x(i,1,1,1) * delta2
299 mij(5)%x(i,1,1,1) = mij(5)%x(i,1,1,1) * delta2
300 mij(6)%x(i,1,1,1) = mij(6)%x(i,1,1,1) * delta2
312 type(
field_t),
intent(inout) :: num, den
313 type(
field_t),
intent(in) :: lij(6), mij(6)
314 real(kind=
rp),
intent(in) :: alpha
315 integer,
intent(in) :: n
317 real(kind=
rp),
dimension(n) :: num_curr, den_curr
320 do concurrent(i = 1:n)
321 num_curr(i) = mij(1)%x(i,1,1,1)*lij(1)%x(i,1,1,1) + &
322 mij(2)%x(i,1,1,1)*lij(2)%x(i,1,1,1) + &
323 mij(3)%x(i,1,1,1)*lij(3)%x(i,1,1,1) + &
324 2.0_rp*(mij(4)%x(i,1,1,1)*lij(4)%x(i,1,1,1) + &
325 mij(5)%x(i,1,1,1)*lij(5)%x(i,1,1,1) + &
326 mij(6)%x(i,1,1,1)*lij(6)%x(i,1,1,1))
327 den_curr(i) = mij(1)%x(i,1,1,1)*mij(1)%x(i,1,1,1) + &
328 mij(2)%x(i,1,1,1)*mij(2)%x(i,1,1,1) + &
329 mij(3)%x(i,1,1,1)*mij(3)%x(i,1,1,1) + &
330 2.0_rp*(mij(4)%x(i,1,1,1)*mij(4)%x(i,1,1,1) + &
331 mij(5)%x(i,1,1,1)*mij(5)%x(i,1,1,1) + &
332 mij(6)%x(i,1,1,1)*mij(6)%x(i,1,1,1))
336 do concurrent(i = 1:n)
337 num%x(i,1,1,1) = alpha * num%x(i,1,1,1) + (1.0_rp - alpha) * num_curr(i)
338 den%x(i,1,1,1) = alpha * den%x(i,1,1,1) + (1.0_rp - alpha) * den_curr(i)
Implements the CPU kernel for the smagorinsky_t type.
subroutine compute_num_den_cpu(num, den, lij, mij, alpha, n)
Compute numerator and denominator for c_dyn on the CPU.
subroutine compute_mij_cpu(mij, s11, s22, s33, s12, s13, s23, s_abs, test_filter, delta, n, nelv)
Compute M_ij on the CPU.
subroutine compute_lij_cpu(lij, u, v, w, test_filter, n, nelv)
Compute Germano Identity on the CPU.
subroutine, public dynamic_smagorinsky_compute_cpu(t, tstep, coef, nut, delta, c_dyn, test_filter, mij, lij, num, den)
Compute eddy viscosity on the CPU.
Implements explicit_filter_t.
Defines a registry for storing solution fields.
type(field_registry_t), target, public neko_field_registry
Global field registry.
Defines Gather-scatter operations.
integer, parameter, public gs_op_add
subroutine, public cmult(a, c, n)
Multiplication by constant c .
subroutine, public cadd(a, s, n)
Add a scalar to vector .
subroutine, public col2(a, b, n)
Vector multiplication .
subroutine, public col3(a, b, c, n)
Vector multiplication with 3 vectors .
real(kind=rp), parameter, public neko_eps
Machine epsilon .
subroutine, public sub2(a, b, n)
Vector substraction .
integer, parameter, public rp
Global precision used in computations.
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,...
Implements the explicit filter for SEM.
field_list_t, To be able to group fields together