60 logical,
intent(in) :: if_ext
61 real(kind=
rp),
intent(in) :: t
62 integer,
intent(in) :: tstep
63 type(
coef_t),
intent(in) :: coef
64 type(
field_t),
intent(inout) :: nut
65 type(
field_t),
intent(in) :: delta
66 real(kind=
rp),
intent(in) :: c_w
67 type(
field_t),
pointer :: u, v, w
69 type(
field_t),
pointer :: g11, g12, g13, g21, g22, g23, g31, g32, g33
71 type(
field_t),
pointer :: s11, s22, s33, s12, s13, s23
73 real(kind=
rp) :: gsqr_11, gsqr_12, gsqr_13, gsqr_21, gsqr_22, gsqr_23, gsqr_31, gsqr_32, gsqr_33
74 real(kind=
rp) :: sd11, sd22, sd33, sd12, sd13, sd23
75 real(kind=
rp) :: sdij_sdij
76 real(kind=
rp) :: sij_sij
77 real(kind=
rp) :: op_wale
78 integer :: temp_indices(15)
82 if (if_ext .eqv. .true.)
then
111 call dudxyz(g11%x, u%x, coef%drdx, coef%dsdx, coef%dtdx, coef)
112 call dudxyz(g12%x, u%x, coef%drdy, coef%dsdy, coef%dtdy, coef)
113 call dudxyz(g13%x, u%x, coef%drdz, coef%dsdz, coef%dtdz, coef)
115 call dudxyz(g21%x, v%x, coef%drdx, coef%dsdx, coef%dtdx, coef)
116 call dudxyz(g22%x, v%x, coef%drdy, coef%dsdy, coef%dtdy, coef)
117 call dudxyz(g23%x, v%x, coef%drdz, coef%dsdz, coef%dtdz, coef)
119 call dudxyz(g31%x, w%x, coef%drdx, coef%dsdx, coef%dtdx, coef)
120 call dudxyz(g32%x, w%x, coef%drdy, coef%dsdy, coef%dtdy, coef)
121 call dudxyz(g33%x, w%x, coef%drdz, coef%dsdz, coef%dtdz, coef)
134 call strain_rate(s11%x, s22%x, s33%x, s12%x, s13%x, s23%x, u, v, w, coef)
144 do concurrent(e = 1:coef%msh%nelv)
145 do concurrent(i = 1:coef%Xh%lxyz)
147 gsqr_11 = g11%x(i,1,1,e)*g11%x(i,1,1,e) + &
148 g12%x(i,1,1,e)*g21%x(i,1,1,e) + &
149 g13%x(i,1,1,e)*g31%x(i,1,1,e)
150 gsqr_12 = g11%x(i,1,1,e)*g12%x(i,1,1,e) + &
151 g12%x(i,1,1,e)*g22%x(i,1,1,e) + &
152 g13%x(i,1,1,e)*g32%x(i,1,1,e)
153 gsqr_13 = g11%x(i,1,1,e)*g13%x(i,1,1,e) + &
154 g12%x(i,1,1,e)*g23%x(i,1,1,e) + &
155 g13%x(i,1,1,e)*g33%x(i,1,1,e)
156 gsqr_21 = g21%x(i,1,1,e)*g11%x(i,1,1,e) + &
157 g22%x(i,1,1,e)*g21%x(i,1,1,e) + &
158 g23%x(i,1,1,e)*g31%x(i,1,1,e)
159 gsqr_22 = g21%x(i,1,1,e)*g12%x(i,1,1,e) + &
160 g22%x(i,1,1,e)*g22%x(i,1,1,e) + &
161 g23%x(i,1,1,e)*g32%x(i,1,1,e)
162 gsqr_23 = g21%x(i,1,1,e)*g13%x(i,1,1,e) + &
163 g22%x(i,1,1,e)*g23%x(i,1,1,e) + &
164 g23%x(i,1,1,e)*g33%x(i,1,1,e)
165 gsqr_31 = g31%x(i,1,1,e)*g11%x(i,1,1,e) + &
166 g32%x(i,1,1,e)*g21%x(i,1,1,e) + &
167 g33%x(i,1,1,e)*g31%x(i,1,1,e)
168 gsqr_32 = g31%x(i,1,1,e)*g12%x(i,1,1,e) + &
169 g32%x(i,1,1,e)*g22%x(i,1,1,e) + &
170 g33%x(i,1,1,e)*g32%x(i,1,1,e)
171 gsqr_33 = g31%x(i,1,1,e)*g13%x(i,1,1,e) + &
172 g32%x(i,1,1,e)*g23%x(i,1,1,e) + &
173 g33%x(i,1,1,e)*g33%x(i,1,1,e)
176 sd11 = gsqr_11 - ( (gsqr_11 + gsqr_22 + gsqr_33) / 3.0_rp)
177 sd22 = gsqr_22 - ( (gsqr_11 + gsqr_22 + gsqr_33) / 3.0_rp)
178 sd33 = gsqr_33 - ( (gsqr_11 + gsqr_22 + gsqr_33) / 3.0_rp)
179 sd12 = 0.5_rp * (gsqr_12 + gsqr_21)
180 sd13 = 0.5_rp * (gsqr_13 + gsqr_31)
181 sd23 = 0.5_rp * (gsqr_23 + gsqr_32)
184 sdij_sdij = sd11*sd11 + sd22*sd22 + sd33*sd33 + &
185 2.0_rp * (sd12*sd12 + sd13*sd13 + sd23*sd23)
187 sij_sij = s11%x(i,1,1,e)*s11%x(i,1,1,e) + s22%x(i,1,1,e)*s22%x(i,1,1,e) + &
188 s33%x(i,1,1,e)*s33%x(i,1,1,e) + 2.0_rp * (s12%x(i,1,1,e)*s12%x(i,1,1,e) + &
189 s13%x(i,1,1,e)*s13%x(i,1,1,e) + s23%x(i,1,1,e)*s23%x(i,1,1,e))
192 op_wale = sdij_sdij**(3.0_rp / 2.0_rp) / &
193 max((sij_sij**(5.0_rp / 2.0_rp) + sdij_sdij**(5.0_rp / 4.0_rp)),
neko_eps)
196 nut%x(i,1,1,e) = c_w**2 * delta%x(i,1,1,e)**2 * op_wale * coef%mult(i,1,1,e)