ADflow  v1.0
ADflow is a finite volume RANS solver tailored for gradient-based aerodynamic design optimization.
kw.F90
Go to the documentation of this file.
1 module kw
2 
3 contains
4 
5  subroutine kw_block(resOnly)
6  !
7  ! kw solves the transport equations for the standard and
8  ! modified k-omega turbulence models in a decoupled manner
9  ! using a diagonal dominant adi-scheme.
10  !
11  use constants
12  use blockpointers, only: il, jl, kl
13  use inputtimespectral
14  use iteration
15  use utils, only: setpointers
16  use turbutils, only: kweddyviscosity
17  use turbbcroutines, only: bcturbtreatment, applyallturbbcthisblock
18  implicit none
19  !
20  ! Subroutine argument.
21  !
22  logical, intent(in) :: resOnly
23 
24  ! Set the arrays for the boundary condition treatment.
25 
26  call bcturbtreatment
27 
28  ! Solve the transport equations for k and omega.
29 
30  call kwsolve(resonly)
31 
32  ! The eddy viscosity and the boundary conditions are only
33  ! applied if an actual update has been computed in kwSolve.
34 
35  if (.not. resonly) then
36 
37  ! Compute the corresponding eddy viscosity.
38 
39  call kweddyviscosity(2, il, 2, jl, 2, kl)
40 
41  ! Set the halo values for the turbulent variables.
42  ! We are on the finest mesh, so the second layer of halo
43  ! cells must be computed as well.
44 
45  call applyallturbbcthisblock(.true.)
46 
47  end if
48 
49  end subroutine kw_block
50 
51  subroutine kwsolve(resOnly)
52  !
53  ! kwSolve solves the k-omega transport equations of both
54  ! the original and modified k-omega models
55  ! in a coupled manner using a diagonal dominant ADI-scheme.
56  !
57  use blockpointers
58  use constants
59  use flowvarrefstate
60  use inputiteration
61  use inputphysics
62  use paramturb
63  use turbmod, only: prod, dvt, sct, kwcd, sig1, sig2, vort
64  use turbutils, only: prodsmag2, prodwmag2, prodkatolaunder, &
65  turbadvection, unsteadyturbterm, tdia3, kwcdterm
66  use turbcurvefits, only: curvetupyp
67  implicit none
68  !
69  ! Subroutine arguments.
70  !
71  logical, intent(in) :: resOnly
72  !
73  ! Local variables.
74  !
75  integer(kind=intType) :: i, j, k, nn
76 
77  real(kind=realtype) :: rkwgam1
78  real(kind=realtype) :: rhoi, ss, spk, sdk
79  real(kind=realtype) :: voli, volmi, volpi
80  real(kind=realtype) :: xm, ym, zm, xp, yp, zp, xa, ya, za
81  real(kind=realtype) :: ttm, ttp, mulm, mulp, muem, muep
82  real(kind=realtype) :: c1m, c1p, c10, c2m, c2p, c20
83  real(kind=realtype) :: b1, b2, c1, c2, d1, d2
84  real(kind=realtype) :: qs, uu, um, up, factor, utau, rblank
85 
86  real(kind=realtype), dimension(itu1:itu2) :: tup
87 
88  real(kind=realtype), dimension(2:il, 2:jl, 2:kl, 2, 2) :: qq
89  real(kind=realtype), dimension(2, 2:max(il, jl, kl)) :: bb, dd, ff
90  real(kind=realtype), dimension(2, 2, 2:max(il, jl, kl)) :: cc
91 
92  real(kind=realtype), dimension(:, :, :), pointer :: ddw, ww, ddvt
93  real(kind=realtype), dimension(:, :), pointer :: rrlv
94  real(kind=realtype), dimension(:, :), pointer :: dd2wall
95 
96  logical, dimension(2:jl, 2:kl), target :: flagI2, flagIl
97  logical, dimension(2:il, 2:kl), target :: flagJ2, flagJl
98  logical, dimension(2:il, 2:jl), target :: flagK2, flagKl
99 
100  logical, dimension(:, :), pointer :: flag
101 
102 
103  ! Set model constants
104 
105  rkwgam1 = rkwbeta1 / rkwbetas - rkwsigw1 * rkwk * rkwk / sqrt(rkwbetas)
106  sig1 = rkwsigk1
107  sig2 = rkwsigw1
108 
109  ! Set the pointer for dvt in dw, such that the code is more
110  ! readable. Also set the pointers for the production term,
111  ! vorticity and the cross diffusion term.
112 
113  dvt => scratch(1:, 1:, 1:, idvt:)
114  prod => scratch(1:, 1:, 1:, iprod)
115  vort => prod
116  kwcd => scratch(1:, 1:, 1:, icd)
117  !
118  ! Production term.
119  !
120  select case (turbprod)
121  case (strain)
122  call prodsmag2
123 
124  case (vorticity)
125  call prodwmag2
126 
127  case (katolaunder)
128  call prodkatolaunder
129 
130  end select
131  !
132  ! Source terms.
133  ! Determine the source term and its derivative w.r.t. k and
134  ! omega for all internal cells of the block.
135  !
136  do k = 2, kl
137  do j = 2, jl
138  do i = 2, il
139 
140  ! Compute the source terms for both the k and the omega
141  ! equation. Note that prod contains the unscaled production
142  ! term. Furthermore the production term of k is limited to
143  ! a certain times the destruction term.
144 
145  rhoi = one / w(i, j, k, irho)
146  ss = prod(i, j, k)
147  spk = rev(i, j, k) * ss * rhoi
148  sdk = rkwbetas * w(i, j, k, itu1) * w(i, j, k, itu2)
149  spk = min(spk, pklim * sdk)
150 
151  dvt(i, j, k, 1) = spk - sdk
152  dvt(i, j, k, 2) = rkwgam1 * ss - rkwbeta1 * w(i, j, k, itu2)**2
153 
154  ! Compute the source term jacobian. Note that only the
155  ! destruction terms are linearized to increase the diagonal
156  ! dominance of the matrix. Furthermore minus the source
157  ! term jacobian is stored.
158 
159  qq(i, j, k, 1, 1) = rkwbetas * w(i, j, k, itu2)
160  ! qq(i,j,k,1,2) = rkwBetas*w(i,j,k,itu1)
161  qq(i, j, k, 1, 2) = zero
162  qq(i, j, k, 2, 1) = zero
163  qq(i, j, k, 2, 2) = two * rkwbeta1 * w(i, j, k, itu2)
164 
165  end do
166  end do
167  end do
168  !
169  ! Compute the cross-diffusion term for the modified version of
170  ! the k-omega model. It should cure the free-stream dependency
171  ! of the original model.
172  !
173  if (turbmodel == komegamodified) then
174 
175  ! Compute the cross diffusion term.
176 
177  call kwcdterm
178 
179  ! Multiply it with the correct constant and take it only into
180  ! account if it is positive.
181 
182  do k = 2, kl
183  do j = 2, jl
184  do i = 2, il
185  dvt(i, j, k, 2) = dvt(i, j, k, 2) &
186  + rkwsigd1 * max(zero, kwcd(i, j, k))
187  end do
188  end do
189  end do
190 
191  end if
192  !
193  ! Advection and unsteady terms.
194  !
195  nn = itu1 - 1
196  call turbadvection(2_inttype, 2_inttype, nn, qq)
197 
198  call unsteadyturbterm(2_inttype, 2_inttype, nn, qq)
199  !
200  ! Viscous terms in k-direction.
201  !
202  do k = 2, kl
203  do j = 2, jl
204  do i = 2, il
205 
206  ! Compute the metrics in zeta-direction, i.e. along the
207  ! line k = constant.
208 
209  voli = one / vol(i, j, k)
210  volmi = two / (vol(i, j, k) + vol(i, j, k - 1))
211  volpi = two / (vol(i, j, k) + vol(i, j, k + 1))
212 
213  xm = sk(i, j, k - 1, 1) * volmi
214  ym = sk(i, j, k - 1, 2) * volmi
215  zm = sk(i, j, k - 1, 3) * volmi
216  xp = sk(i, j, k, 1) * volpi
217  yp = sk(i, j, k, 2) * volpi
218  zp = sk(i, j, k, 3) * volpi
219 
220  xa = half * (sk(i, j, k, 1) + sk(i, j, k - 1, 1)) * voli
221  ya = half * (sk(i, j, k, 2) + sk(i, j, k - 1, 2)) * voli
222  za = half * (sk(i, j, k, 3) + sk(i, j, k - 1, 3)) * voli
223  ttm = xm * xa + ym * ya + zm * za
224  ttp = xp * xa + yp * ya + zp * za
225 
226  ! Computation of the viscous terms in zeta-direction; note
227  ! that cross-derivatives are neglected, i.e. the mesh is
228  ! assumed to be orthogonal.
229  ! The second derivative in zeta-direction is constructed as
230  ! the central difference of the first order derivatives, i.e.
231  ! d^2/dzeta^2 = d/dzeta (d/dzeta k+1/2 - d/dzeta k-1/2).
232  ! In this way the metric as well as the varying viscosity
233  ! can be taken into account; the latter appears inside the
234  ! d/dzeta derivative. The whole term is divided by rho to
235  ! obtain the diffusion term for k and omega.
236 
237  rhoi = one / w(i, j, k, irho)
238  mulm = half * (rlv(i, j, k - 1) + rlv(i, j, k))
239  mulp = half * (rlv(i, j, k + 1) + rlv(i, j, k))
240  muem = half * (rev(i, j, k - 1) + rev(i, j, k))
241  muep = half * (rev(i, j, k + 1) + rev(i, j, k))
242 
243  c1m = ttm * (mulm + sig1 * muem) * rhoi
244  c1p = ttp * (mulp + sig1 * muep) * rhoi
245  c10 = c1m + c1p
246 
247  c2m = ttm * (mulm + sig2 * muem) * rhoi
248  c2p = ttp * (mulp + sig2 * muep) * rhoi
249  c20 = c2m + c2p
250 
251  ! Update the residual for this cell and store the possible
252  ! coefficients for the matrix in b1, b2, c1, c2, d1 and d2.
253 
254  dvt(i, j, k, 1) = dvt(i, j, k, 1) + c1m * w(i, j, k - 1, itu1) &
255  - c10 * w(i, j, k, itu1) + c1p * w(i, j, k + 1, itu1)
256  dvt(i, j, k, 2) = dvt(i, j, k, 2) + c2m * w(i, j, k - 1, itu2) &
257  - c20 * w(i, j, k, itu2) + c2p * w(i, j, k + 1, itu2)
258 
259  b1 = -c1m
260  c1 = c10
261  d1 = -c1p
262 
263  b2 = -c2m
264  c2 = c20
265  d2 = -c2p
266 
267  ! Update the central jacobian. For nonboundary cells this
268  ! is simply c1 and c2. For boundary cells this is slightly
269  ! more complicated, because the boundary conditions are
270  ! treated implicitly and the off-diagonal terms b1, b2 and
271  ! d1, d2 must be taken into account.
272  ! The boundary conditions are only treated implicitly if
273  ! the diagonal dominance of the matrix is increased.
274 
275  if (k == 2) then
276  qq(i, j, k, 1, 1) = qq(i, j, k, 1, 1) + c1 &
277  - b1 * max(bmtk1(i, j, itu1, itu1), zero)
278  qq(i, j, k, 1, 2) = qq(i, j, k, 1, 2) - b1 * bmtk1(i, j, itu1, itu2)
279  qq(i, j, k, 2, 1) = qq(i, j, k, 2, 1) - b2 * bmtk1(i, j, itu2, itu1)
280  qq(i, j, k, 2, 2) = qq(i, j, k, 2, 2) + c2 &
281  - b2 * max(bmtk1(i, j, itu2, itu2), zero)
282  else if (k == kl) then
283  qq(i, j, k, 1, 1) = qq(i, j, k, 1, 1) + c1 &
284  - d1 * max(bmtk2(i, j, itu1, itu1), zero)
285  qq(i, j, k, 1, 2) = qq(i, j, k, 1, 2) - d1 * bmtk2(i, j, itu1, itu2)
286  qq(i, j, k, 2, 1) = qq(i, j, k, 2, 1) - d2 * bmtk2(i, j, itu2, itu1)
287  qq(i, j, k, 2, 2) = qq(i, j, k, 2, 2) + c2 &
288  - d2 * max(bmtk2(i, j, itu2, itu2), zero)
289  else
290  qq(i, j, k, 1, 1) = qq(i, j, k, 1, 1) + c1
291  qq(i, j, k, 2, 2) = qq(i, j, k, 2, 2) + c2
292  end if
293 
294  end do
295  end do
296  end do
297  !
298  ! Viscous terms in j-direction.
299  !
300  do k = 2, kl
301  do j = 2, jl
302  do i = 2, il
303 
304  ! Compute the metrics in eta-direction, i.e. along the
305  ! line j = constant.
306 
307  voli = one / vol(i, j, k)
308  volmi = two / (vol(i, j, k) + vol(i, j - 1, k))
309  volpi = two / (vol(i, j, k) + vol(i, j + 1, k))
310 
311  xm = sj(i, j - 1, k, 1) * volmi
312  ym = sj(i, j - 1, k, 2) * volmi
313  zm = sj(i, j - 1, k, 3) * volmi
314  xp = sj(i, j, k, 1) * volpi
315  yp = sj(i, j, k, 2) * volpi
316  zp = sj(i, j, k, 3) * volpi
317 
318  xa = half * (sj(i, j, k, 1) + sj(i, j - 1, k, 1)) * voli
319  ya = half * (sj(i, j, k, 2) + sj(i, j - 1, k, 2)) * voli
320  za = half * (sj(i, j, k, 3) + sj(i, j - 1, k, 3)) * voli
321  ttm = xm * xa + ym * ya + zm * za
322  ttp = xp * xa + yp * ya + zp * za
323 
324  ! Computation of the viscous terms in eta-direction; note
325  ! that cross-derivatives are neglected, i.e. the mesh is
326  ! assumed to be orthogonal.
327  ! The second derivative in eta-direction is constructed as
328  ! the central difference of the first order derivatives, i.e.
329  ! d^2/deta^2 = d/deta (d/deta j+1/2 - d/deta j-1/2).
330  ! In this way the metric as well as the varying viscosity
331  ! can be taken into account; the latter appears inside the
332  ! d/deta derivative. The whole term is divided by rho to
333  ! obtain the diffusion term for k and omega.
334 
335  rhoi = one / w(i, j, k, irho)
336  mulm = half * (rlv(i, j - 1, k) + rlv(i, j, k))
337  mulp = half * (rlv(i, j + 1, k) + rlv(i, j, k))
338  muem = half * (rev(i, j - 1, k) + rev(i, j, k))
339  muep = half * (rev(i, j + 1, k) + rev(i, j, k))
340 
341  c1m = ttm * (mulm + sig1 * muem) * rhoi
342  c1p = ttp * (mulp + sig1 * muep) * rhoi
343  c10 = c1m + c1p
344 
345  c2m = ttm * (mulm + sig2 * muem) * rhoi
346  c2p = ttp * (mulp + sig2 * muep) * rhoi
347  c20 = c2m + c2p
348 
349  ! Update the residual for this cell and store the possible
350  ! coefficients for the matrix in b1, b2, c1, c2, d1 and d2.
351 
352  dvt(i, j, k, 1) = dvt(i, j, k, 1) + c1m * w(i, j - 1, k, itu1) &
353  - c10 * w(i, j, k, itu1) + c1p * w(i, j + 1, k, itu1)
354  dvt(i, j, k, 2) = dvt(i, j, k, 2) + c2m * w(i, j - 1, k, itu2) &
355  - c20 * w(i, j, k, itu2) + c2p * w(i, j + 1, k, itu2)
356 
357  b1 = -c1m
358  c1 = c10
359  d1 = -c1p
360 
361  b2 = -c2m
362  c2 = c20
363  d2 = -c2p
364 
365  ! Update the central jacobian. For nonboundary cells this
366  ! is simply c1 and c2. For boundary cells this is slightly
367  ! more complicated, because the boundary conditions are
368  ! treated implicitly and the off-diagonal terms b1, b2 and
369  ! d1, d2 must be taken into account.
370  ! The boundary conditions are only treated implicitly if
371  ! the diagonal dominance of the matrix is increased.
372 
373  if (j == 2) then
374  qq(i, j, k, 1, 1) = qq(i, j, k, 1, 1) + c1 &
375  - b1 * max(bmtj1(i, k, itu1, itu1), zero)
376  qq(i, j, k, 1, 2) = qq(i, j, k, 1, 2) - b1 * bmtj1(i, k, itu1, itu2)
377  qq(i, j, k, 2, 1) = qq(i, j, k, 2, 1) - b2 * bmtj1(i, k, itu2, itu1)
378  qq(i, j, k, 2, 2) = qq(i, j, k, 2, 2) + c2 &
379  - b2 * max(bmtj1(i, k, itu2, itu2), zero)
380  else if (j == jl) then
381  qq(i, j, k, 1, 1) = qq(i, j, k, 1, 1) + c1 &
382  - d1 * max(bmtj2(i, k, itu1, itu1), zero)
383  qq(i, j, k, 1, 2) = qq(i, j, k, 1, 2) - d1 * bmtj2(i, k, itu1, itu2)
384  qq(i, j, k, 2, 1) = qq(i, j, k, 2, 1) - d2 * bmtj2(i, k, itu2, itu1)
385  qq(i, j, k, 2, 2) = qq(i, j, k, 2, 2) + c2 &
386  - d2 * max(bmtj2(i, k, itu2, itu2), zero)
387  else
388  qq(i, j, k, 1, 1) = qq(i, j, k, 1, 1) + c1
389  qq(i, j, k, 2, 2) = qq(i, j, k, 2, 2) + c2
390  end if
391 
392  end do
393  end do
394  end do
395  !
396  ! Viscous terms in i-direction.
397  !
398  do k = 2, kl
399  do j = 2, jl
400  do i = 2, il
401 
402  ! Compute the metrics in xi-direction, i.e. along the
403  ! line i = constant.
404 
405  voli = one / vol(i, j, k)
406  volmi = two / (vol(i, j, k) + vol(i - 1, j, k))
407  volpi = two / (vol(i, j, k) + vol(i + 1, j, k))
408 
409  xm = si(i - 1, j, k, 1) * volmi
410  ym = si(i - 1, j, k, 2) * volmi
411  zm = si(i - 1, j, k, 3) * volmi
412  xp = si(i, j, k, 1) * volpi
413  yp = si(i, j, k, 2) * volpi
414  zp = si(i, j, k, 3) * volpi
415 
416  xa = half * (si(i, j, k, 1) + si(i - 1, j, k, 1)) * voli
417  ya = half * (si(i, j, k, 2) + si(i - 1, j, k, 2)) * voli
418  za = half * (si(i, j, k, 3) + si(i - 1, j, k, 3)) * voli
419  ttm = xm * xa + ym * ya + zm * za
420  ttp = xp * xa + yp * ya + zp * za
421 
422  ! Computation of the viscous terms in xi-direction; note
423  ! that cross-derivatives are neglected, i.e. the mesh is
424  ! assumed to be orthogonal.
425  ! The second derivative in xi-direction is constructed as
426  ! the central difference of the first order derivatives, i.e.
427  ! d^2/dxi^2 = d/dxi (d/dxi i+1/2 - d/dxi i-1/2).
428  ! In this way the metric as well as the varying viscosity
429  ! can be taken into account; the latter appears inside the
430  ! d/dxi derivative. The whole term is divided by rho to
431  ! obtain the diffusion term for k and omega.
432 
433  rhoi = one / w(i, j, k, irho)
434  mulm = half * (rlv(i - 1, j, k) + rlv(i, j, k))
435  mulp = half * (rlv(i + 1, j, k) + rlv(i, j, k))
436  muem = half * (rev(i - 1, j, k) + rev(i, j, k))
437  muep = half * (rev(i + 1, j, k) + rev(i, j, k))
438 
439  c1m = ttm * (mulm + sig1 * muem) * rhoi
440  c1p = ttp * (mulp + sig1 * muep) * rhoi
441  c10 = c1m + c1p
442 
443  c2m = ttm * (mulm + sig2 * muem) * rhoi
444  c2p = ttp * (mulp + sig2 * muep) * rhoi
445  c20 = c2m + c2p
446 
447  ! Update the residual for this cell and store the possible
448  ! coefficients for the matrix in b1, b2, c1, c2, d1 and d2.
449 
450  dvt(i, j, k, 1) = dvt(i, j, k, 1) + c1m * w(i - 1, j, k, itu1) &
451  - c10 * w(i, j, k, itu1) + c1p * w(i + 1, j, k, itu1)
452  dvt(i, j, k, 2) = dvt(i, j, k, 2) + c2m * w(i - 1, j, k, itu2) &
453  - c20 * w(i, j, k, itu2) + c2p * w(i + 1, j, k, itu2)
454 
455  b1 = -c1m
456  c1 = c10
457  d1 = -c1p
458 
459  b2 = -c2m
460  c2 = c20
461  d2 = -c2p
462 
463  ! Update the central jacobian. For nonboundary cells this
464  ! is simply c1 and c2. For boundary cells this is slightly
465  ! more complicated, because the boundary conditions are
466  ! treated implicitly and the off-diagonal terms b1, b2 and
467  ! d1, d2 must be taken into account.
468  ! The boundary conditions are only treated implicitly if
469  ! the diagonal dominance of the matrix is increased.
470 
471  if (i == 2) then
472  qq(i, j, k, 1, 1) = qq(i, j, k, 1, 1) + c1 &
473  - b1 * max(bmti1(j, k, itu1, itu1), zero)
474  qq(i, j, k, 1, 2) = qq(i, j, k, 1, 2) - b1 * bmti1(j, k, itu1, itu2)
475  qq(i, j, k, 2, 1) = qq(i, j, k, 2, 1) - b2 * bmti1(j, k, itu2, itu1)
476  qq(i, j, k, 2, 2) = qq(i, j, k, 2, 2) + c2 &
477  - b2 * max(bmti1(j, k, itu2, itu2), zero)
478  else if (i == il) then
479  qq(i, j, k, 1, 1) = qq(i, j, k, 1, 1) + c1 &
480  - d1 * max(bmti2(j, k, itu1, itu1), zero)
481  qq(i, j, k, 1, 2) = qq(i, j, k, 1, 2) - d1 * bmti2(j, k, itu1, itu2)
482  qq(i, j, k, 2, 1) = qq(i, j, k, 2, 1) - d2 * bmti2(j, k, itu2, itu1)
483  qq(i, j, k, 2, 2) = qq(i, j, k, 2, 2) + c2 &
484  - d2 * max(bmti2(j, k, itu2, itu2), zero)
485  else
486  qq(i, j, k, 1, 1) = qq(i, j, k, 1, 1) + c1
487  qq(i, j, k, 2, 2) = qq(i, j, k, 2, 2) + c2
488  end if
489 
490  end do
491  end do
492  end do
493 
494  ! Multiply the residual by the volume and store this in dw; this
495  ! is done for monitoring reasons only. The multiplication with the
496  ! volume is present to be consistent with the flow residuals; also
497  ! the negative value is taken, again to be consistent with the
498  ! flow equations. Also multiply by iblank so that no updates occur
499  ! in holes or the overset boundary.
500 
501  do k = 2, kl
502  do j = 2, jl
503  do i = 2, il
504  rblank = real(iblank(i, j, k), realtype)
505  dw(i, j, k, itu1) = -volref(i, j, k) * dvt(i, j, k, 1) * rblank
506  dw(i, j, k, itu2) = -volref(i, j, k) * dvt(i, j, k, 2) * rblank
507  end do
508  end do
509  end do
510 
511  ! Initialize the wall function flags to .false.
512 
513  flagi2 = .false.
514  flagil = .false.
515  flagj2 = .false.
516  flagjl = .false.
517  flagk2 = .false.
518  flagkl = .false.
519 
520  ! Modify the rhs of the 1st internal cell, if wall functions
521  ! are used; their value is determined by the table.
522 
523  testwallfunctions: if (wallfunctions) then
524 
525  bocos: do nn = 1, nviscbocos
526 
527  ! Determine the block face on which the subface is located
528  ! and set some variables. As flag points to the entire array
529  ! flagI2, etc., its starting indices are the starting indices
530  ! of its target and not 1.
531 
532  select case (bcfaceid(nn))
533  case (imin)
534  flag => flagi2
535  ddw => dw(2, 1:, 1:, 1:); ddvt => dvt(2, 1:, 1:, 1:)
536  ww => w(2, 1:, 1:, 1:); rrlv => rlv(2, 1:, 1:)
537  dd2wall => d2wall(2, :, :)
538 
539  case (imax)
540  flag => flagil
541  ddw => dw(il, 1:, 1:, 1:); ddvt => dvt(il, 1:, 1:, 1:)
542  ww => w(il, 1:, 1:, 1:); rrlv => rlv(il, 1:, 1:)
543  dd2wall => d2wall(il, :, :)
544 
545  case (jmin)
546  flag => flagj2
547  ddw => dw(1:, 2, 1:, 1:); ddvt => dvt(1:, 2, 1:, 1:)
548  ww => w(1:, 2, 1:, 1:); rrlv => rlv(1:, 2, 1:)
549  dd2wall => d2wall(:, 2, :)
550 
551  case (jmax)
552  flag => flagjl
553  ddw => dw(1:, jl, 1:, 1:); ddvt => dvt(1:, jl, 1:, 1:)
554  ww => w(1:, jl, 1:, 1:); rrlv => rlv(1:, jl, 1:)
555  dd2wall => d2wall(:, jl, :)
556 
557  case (kmin)
558  flag => flagk2
559  ddw => dw(1:, 1:, 2, 1:); ddvt => dvt(1:, 1:, 2, 1:)
560  ww => w(1:, 1:, 2, 1:); rrlv => rlv(1:, 1:, 2)
561  dd2wall => d2wall(:, :, 2)
562 
563  case (kmax)
564  flag => flagkl
565  ddw => dw(1:, 1:, kl, 1:); ddvt => dvt(1:, 1:, kl, 1:)
566  ww => w(1:, 1:, kl, 1:); rrlv => rlv(1:, 1:, kl)
567  dd2wall => d2wall(:, :, kl)
568 
569  end select
570 
571  ! Loop over the owned faces of this subface. Therefore the
572  ! nodal range of BCData must be used. The offset of +1 is
573  ! present, because the starting index of the cell range is
574  ! 1 larger than the starting index of the nodal range.
575 
576  do j = (bcdata(nn)%jnBeg + 1), bcdata(nn)%jnEnd
577  do i = (bcdata(nn)%inBeg + 1), bcdata(nn)%inEnd
578 
579  ! Set ddw to zero for proper monitoring of the
580  ! convergence.
581 
582  ddw(i, j, itu1) = zero
583  ddw(i, j, itu2) = zero
584 
585  ! Enforce k and omega in the 1st internal cell from
586  ! the wall function table. There is an offset of -1 in
587  ! the wall distance. Note that the offset compared to
588  ! the current value must be stored. Also note that the
589  ! curve fits contain the non-dimensional values.
590 
591  utau = viscsubface(nn)%utau(i, j)
592  yp = ww(i, j, irho) * dd2wall(i - 1, j - 1) * utau / rrlv(i, j)
593 
594  call curvetupyp(tup, yp, itu1, itu2)
595 
596  tup(itu1) = tup(itu1) * utau**2
597  tup(itu2) = tup(itu2) * utau**2 / rrlv(i, j) * ww(i, j, irho)
598 
599  ddvt(i, j, 1) = tup(itu1) - ww(i, j, itu1)
600  ddvt(i, j, 2) = tup(itu2) - ww(i, j, itu2)
601 
602  ! Set the wall flag to .true.
603 
604  flag(i, j) = .true.
605 
606  end do
607  end do
608 
609  end do bocos
610  end if testwallfunctions
611 
612  ! Return if only the residual must be computed.
613 
614  if (resonly) return
615 
616  ! For implicit relaxation take the local time step into account,
617  ! where dt is the inverse of the central jacobian times the cfl
618  ! number. The following system is solved:
619  ! (i/dt + cc + bb + dd)*dw = rhs, in which i/dt = cc/cfl. As in
620  ! the rest of the algorithm only the modified central jacobian is
621  ! used, stored it now.
622 
623  ! Compute the factor multiplying the central jacobian, which
624  ! is 1 + 1/cfl (implicit relaxation only).
625 
626  factor = one
627  if (turbrelax == turbrelaximplicit) &
628  factor = one + (one - alfaturb) / alfaturb
629 
630  do k = 2, kl
631  do j = 2, jl
632  do i = 2, il
633  qq(i, j, k, 1, 1) = factor * qq(i, j, k, 1, 1)
634  qq(i, j, k, 1, 2) = factor * qq(i, j, k, 1, 2)
635  qq(i, j, k, 2, 1) = factor * qq(i, j, k, 2, 1)
636  qq(i, j, k, 2, 2) = factor * qq(i, j, k, 2, 2)
637 
638  ! Set qq to 1 if the value is determined by the table.
639 
640  if ((i == 2 .and. flagi2(j, k)) .or. &
641  (i == il .and. flagil(j, k)) .or. &
642  (j == 2 .and. flagj2(i, k)) .or. &
643  (j == jl .and. flagjl(i, k)) .or. &
644  (k == 2 .and. flagk2(i, j)) .or. &
645  (k == kl .and. flagkl(i, j))) then
646  qq(i, j, k, 1, 1) = one
647  qq(i, j, k, 1, 2) = zero
648  qq(i, j, k, 2, 1) = zero
649  qq(i, j, k, 2, 2) = one
650  end if
651  end do
652  end do
653  end do
654 
655  ! Initialize the grid velocity to zero. This value will be used
656  ! if the block is not moving.
657 
658  qs = zero
659  !
660  ! dd-ADI step in j-direction. There is no particular reason to
661  ! start in j-direction, it just happened to be so. As we solve
662  ! in j-direction, the j-loop is the innermost loop.
663  !
664  do k = 2, kl
665  do i = 2, il
666  do j = 2, jl
667 
668  ! More or less the same code is executed here as above when
669  ! the residual was built. However, now the off-diagonal
670  ! terms for the dd-ADI must be built and stored. This could
671  ! have been done earlier, but then all the coefficients had
672  ! to be stored. To save memory, they are recomputed.
673  ! Consequently, see the j-loop to build the residual for
674  ! the comments.
675 
676  voli = one / vol(i, j, k)
677  volmi = two / (vol(i, j, k) + vol(i, j - 1, k))
678  volpi = two / (vol(i, j, k) + vol(i, j + 1, k))
679 
680  xm = sj(i, j - 1, k, 1) * volmi
681  ym = sj(i, j - 1, k, 2) * volmi
682  zm = sj(i, j - 1, k, 3) * volmi
683  xp = sj(i, j, k, 1) * volpi
684  yp = sj(i, j, k, 2) * volpi
685  zp = sj(i, j, k, 3) * volpi
686 
687  xa = half * (sj(i, j, k, 1) + sj(i, j - 1, k, 1)) * voli
688  ya = half * (sj(i, j, k, 2) + sj(i, j - 1, k, 2)) * voli
689  za = half * (sj(i, j, k, 3) + sj(i, j - 1, k, 3)) * voli
690  ttm = xm * xa + ym * ya + zm * za
691  ttp = xp * xa + yp * ya + zp * za
692 
693  ! Off-diagonal terms due to the diffusion terms
694  ! in j-direction.
695 
696  rhoi = one / w(i, j, k, irho)
697  mulm = half * (rlv(i, j - 1, k) + rlv(i, j, k))
698  mulp = half * (rlv(i, j + 1, k) + rlv(i, j, k))
699  muem = half * (rev(i, j - 1, k) + rev(i, j, k))
700  muep = half * (rev(i, j + 1, k) + rev(i, j, k))
701 
702  c1m = ttm * (mulm + sig1 * muem) * rhoi
703  c1p = ttp * (mulp + sig1 * muep) * rhoi
704 
705  c2m = ttm * (mulm + sig2 * muem) * rhoi
706  c2p = ttp * (mulp + sig2 * muep) * rhoi
707 
708  bb(1, j) = -c1m
709  dd(1, j) = -c1p
710  bb(2, j) = -c2m
711  dd(2, j) = -c2p
712 
713  ! Compute the grid velocity if present.
714  ! It is taken as the average of j and j-1,
715 
716  if (addgridvelocities) &
717  qs = half * (sfacej(i, j, k) + sfacej(i, j - 1, k)) * voli
718 
719  ! Off-diagonal terms due to the advection term in
720  ! j-direction. First order approximation.
721 
722  uu = xa * w(i, j, k, ivx) + ya * w(i, j, k, ivy) + za * w(i, j, k, ivz) - qs
723  um = zero
724  up = zero
725  if (uu < zero) um = uu
726  if (uu > zero) up = uu
727 
728  bb(1, j) = bb(1, j) - up
729  dd(1, j) = dd(1, j) + um
730  bb(2, j) = bb(2, j) - up
731  dd(2, j) = dd(2, j) + um
732 
733  ! Store the central jacobian and rhs in cc and ff.
734  ! Multiply the off-diagonal terms and rhs by the iblank
735  ! value so the update determined for iblank = 0 is zero.
736 
737  rblank = real(iblank(i, j, k), realtype)
738 
739  cc(1, 1, j) = qq(i, j, k, 1, 1)
740  cc(1, 2, j) = qq(i, j, k, 1, 2) * rblank
741  cc(2, 1, j) = qq(i, j, k, 2, 1) * rblank
742  cc(2, 2, j) = qq(i, j, k, 2, 2)
743 
744  ff(1, j) = dvt(i, j, k, 1) * rblank
745  ff(2, j) = dvt(i, j, k, 2) * rblank
746 
747  bb(:, j) = bb(:, j) * rblank
748  dd(:, j) = dd(:, j) * rblank
749 
750  ! Set off diagonal terms to zero if wall function are used.
751 
752  if ((i == 2 .and. flagi2(j, k)) .or. &
753  (i == il .and. flagil(j, k)) .or. &
754  (j == 2 .and. flagj2(i, k)) .or. &
755  (j == jl .and. flagjl(i, k)) .or. &
756  (k == 2 .and. flagk2(i, j)) .or. &
757  (k == kl .and. flagkl(i, j))) then
758  bb(1, j) = zero
759  dd(1, j) = zero
760  bb(2, j) = zero
761  dd(2, j) = zero
762  end if
763 
764  end do
765 
766  ! Solve the tri-diagonal system in j-direction.
767 
768  call tdia3(2_inttype, jl, bb, cc, dd, ff)
769 
770  ! Determine the new rhs for the next direction.
771 
772  do j = 2, jl
773  dvt(i, j, k, 1) = qq(i, j, k, 1, 1) * ff(1, j) + qq(i, j, k, 1, 2) * ff(2, j)
774  dvt(i, j, k, 2) = qq(i, j, k, 2, 1) * ff(1, j) + qq(i, j, k, 2, 2) * ff(2, j)
775  end do
776 
777  end do
778  end do
779  !
780  ! dd-ADI step in i-direction. As we solve in i-direction, the
781  ! i-loop is the innermost loop.
782  !
783  do k = 2, kl
784  do j = 2, jl
785  do i = 2, il
786 
787  ! More or less the same code is executed here as above when
788  ! the residual was built. However, now the off-diagonal
789  ! terms for the dd-ADI must be built and stored. This could
790  ! have been done earlier, but then all the coefficients had
791  ! to be stored. To save memory, they are recomputed.
792  ! Consequently, see the i-loop to build the residual for
793  ! the comments.
794 
795  voli = one / vol(i, j, k)
796  volmi = two / (vol(i, j, k) + vol(i - 1, j, k))
797  volpi = two / (vol(i, j, k) + vol(i + 1, j, k))
798 
799  xm = si(i - 1, j, k, 1) * volmi
800  ym = si(i - 1, j, k, 2) * volmi
801  zm = si(i - 1, j, k, 3) * volmi
802  xp = si(i, j, k, 1) * volpi
803  yp = si(i, j, k, 2) * volpi
804  zp = si(i, j, k, 3) * volpi
805 
806  xa = half * (si(i, j, k, 1) + si(i - 1, j, k, 1)) * voli
807  ya = half * (si(i, j, k, 2) + si(i - 1, j, k, 2)) * voli
808  za = half * (si(i, j, k, 3) + si(i - 1, j, k, 3)) * voli
809  ttm = xm * xa + ym * ya + zm * za
810  ttp = xp * xa + yp * ya + zp * za
811 
812  ! Off-diagonal terms due to the diffusion terms
813  ! in i-direction.
814 
815  rhoi = one / w(i, j, k, irho)
816  mulm = half * (rlv(i - 1, j, k) + rlv(i, j, k))
817  mulp = half * (rlv(i + 1, j, k) + rlv(i, j, k))
818  muem = half * (rev(i - 1, j, k) + rev(i, j, k))
819  muep = half * (rev(i + 1, j, k) + rev(i, j, k))
820 
821  c1m = ttm * (mulm + sig1 * muem) * rhoi
822  c1p = ttp * (mulp + sig1 * muep) * rhoi
823 
824  c2m = ttm * (mulm + sig2 * muem) * rhoi
825  c2p = ttp * (mulp + sig2 * muep) * rhoi
826 
827  bb(1, i) = -c1m
828  dd(1, i) = -c1p
829  bb(2, i) = -c2m
830  dd(2, i) = -c2p
831 
832  ! Compute the grid velocity if present.
833  ! It is taken as the average of i and i-1,
834 
835  if (addgridvelocities) &
836  qs = half * (sfacei(i, j, k) + sfacei(i - 1, j, k)) * voli
837 
838  ! Off-diagonal terms due to the advection term in
839  ! i-direction. First order approximation.
840 
841  uu = xa * w(i, j, k, ivx) + ya * w(i, j, k, ivy) + za * w(i, j, k, ivz) - qs
842  um = zero
843  up = zero
844  if (uu < zero) um = uu
845  if (uu > zero) up = uu
846 
847  bb(1, i) = bb(1, i) - up
848  dd(1, i) = dd(1, i) + um
849  bb(2, i) = bb(2, i) - up
850  dd(2, i) = dd(2, i) + um
851 
852  ! Store the central jacobian and rhs in cc and ff.
853  ! Multiply the off-diagonal terms and rhs by the iblank
854  ! value so the update determined for iblank = 0 is zero.
855 
856  rblank = real(iblank(i, j, k), realtype)
857 
858  cc(1, 1, i) = qq(i, j, k, 1, 1)
859  cc(1, 2, i) = qq(i, j, k, 1, 2) * rblank
860  cc(2, 1, i) = qq(i, j, k, 2, 1) * rblank
861  cc(2, 2, i) = qq(i, j, k, 2, 2)
862 
863  ff(1, i) = dvt(i, j, k, 1) * rblank
864  ff(2, i) = dvt(i, j, k, 2) * rblank
865 
866  bb(:, i) = bb(:, i) * rblank
867  dd(:, i) = dd(:, i) * rblank
868 
869  ! Set off diagonal terms to zero if wall function are used.
870 
871  if ((i == 2 .and. flagi2(j, k)) .or. &
872  (i == il .and. flagil(j, k)) .or. &
873  (j == 2 .and. flagj2(i, k)) .or. &
874  (j == jl .and. flagjl(i, k)) .or. &
875  (k == 2 .and. flagk2(i, j)) .or. &
876  (k == kl .and. flagkl(i, j))) then
877  bb(1, i) = zero
878  dd(1, i) = zero
879  bb(2, i) = zero
880  dd(2, i) = zero
881  end if
882 
883  end do
884 
885  ! Solve the tri-diagonal system in i-direction.
886 
887  call tdia3(2_inttype, il, bb, cc, dd, ff)
888 
889  ! Determine the new rhs for the next direction.
890 
891  do i = 2, il
892  dvt(i, j, k, 1) = qq(i, j, k, 1, 1) * ff(1, i) + qq(i, j, k, 1, 2) * ff(2, i)
893  dvt(i, j, k, 2) = qq(i, j, k, 2, 1) * ff(1, i) + qq(i, j, k, 2, 2) * ff(2, i)
894  end do
895 
896  end do
897  end do
898  !
899  ! dd-adi step in k-direction. As we solve in k-direction, the
900  ! k-loop is the innermost loop.
901  !
902  do j = 2, jl
903  do i = 2, il
904  do k = 2, kl
905 
906  ! More or less the same code is executed here as above when
907  ! the residual was built. However, now the off-diagonal
908  ! terms for the dd-ADI must be built and stored. This could
909  ! have been done earlier, but then all the coefficients had
910  ! to be stored. To save memory, they are recomputed.
911  ! Consequently, see the k-loop to build the residual for
912  ! the comments.
913 
914  voli = one / vol(i, j, k)
915  volmi = two / (vol(i, j, k) + vol(i, j, k - 1))
916  volpi = two / (vol(i, j, k) + vol(i, j, k + 1))
917 
918  xm = sk(i, j, k - 1, 1) * volmi
919  ym = sk(i, j, k - 1, 2) * volmi
920  zm = sk(i, j, k - 1, 3) * volmi
921  xp = sk(i, j, k, 1) * volpi
922  yp = sk(i, j, k, 2) * volpi
923  zp = sk(i, j, k, 3) * volpi
924 
925  xa = half * (sk(i, j, k, 1) + sk(i, j, k - 1, 1)) * voli
926  ya = half * (sk(i, j, k, 2) + sk(i, j, k - 1, 2)) * voli
927  za = half * (sk(i, j, k, 3) + sk(i, j, k - 1, 3)) * voli
928  ttm = xm * xa + ym * ya + zm * za
929  ttp = xp * xa + yp * ya + zp * za
930 
931  ! Off-diagonal terms due to the diffusion terms
932  ! in k-direction.
933 
934  rhoi = one / w(i, j, k, irho)
935  mulm = half * (rlv(i, j, k - 1) + rlv(i, j, k))
936  mulp = half * (rlv(i, j, k + 1) + rlv(i, j, k))
937  muem = half * (rev(i, j, k - 1) + rev(i, j, k))
938  muep = half * (rev(i, j, k + 1) + rev(i, j, k))
939 
940  c1m = ttm * (mulm + sig1 * muem) * rhoi
941  c1p = ttp * (mulp + sig1 * muep) * rhoi
942 
943  c2m = ttm * (mulm + sig2 * muem) * rhoi
944  c2p = ttp * (mulp + sig2 * muep) * rhoi
945 
946  bb(1, k) = -c1m
947  dd(1, k) = -c1p
948  bb(2, k) = -c2m
949  dd(2, k) = -c2p
950 
951  ! Compute the grid velocity if present.
952  ! It is taken as the average of k and k-1,
953 
954  if (addgridvelocities) &
955  qs = half * (sfacek(i, j, k) + sfacek(i, j, k - 1)) * voli
956 
957  ! Off-diagonal terms due to the advection term in
958  ! k-direction. First order approximation.
959 
960  uu = xa * w(i, j, k, ivx) + ya * w(i, j, k, ivy) + za * w(i, j, k, ivz) - qs
961  um = zero
962  up = zero
963  if (uu < zero) um = uu
964  if (uu > zero) up = uu
965 
966  bb(1, k) = bb(1, k) - up
967  dd(1, k) = dd(1, k) + um
968  bb(2, k) = bb(2, k) - up
969  dd(2, k) = dd(2, k) + um
970 
971  ! Store the central jacobian and rhs in cc and ff.
972  ! Multiply the off-diagonal terms and rhs by the iblank
973  ! value so the update determined for iblank = 0 is zero.
974 
975  rblank = real(iblank(i, j, k), realtype)
976 
977  cc(1, 1, k) = qq(i, j, k, 1, 1)
978  cc(1, 2, k) = qq(i, j, k, 1, 2) * rblank
979  cc(2, 1, k) = qq(i, j, k, 2, 1) * rblank
980  cc(2, 2, k) = qq(i, j, k, 2, 2)
981 
982  ff(1, k) = dvt(i, j, k, 1) * rblank
983  ff(2, k) = dvt(i, j, k, 2) * rblank
984 
985  bb(:, k) = bb(:, k) * rblank
986  dd(:, k) = dd(:, k) * rblank
987 
988  ! Set off diagonal terms to zero if wall function are used.
989 
990  if ((i == 2 .and. flagi2(j, k)) .or. &
991  (i == il .and. flagil(j, k)) .or. &
992  (j == 2 .and. flagj2(i, k)) .or. &
993  (j == jl .and. flagjl(i, k)) .or. &
994  (k == 2 .and. flagk2(i, j)) .or. &
995  (k == kl .and. flagkl(i, j))) then
996  bb(1, k) = zero
997  dd(1, k) = zero
998  bb(2, k) = zero
999  dd(2, k) = zero
1000  end if
1001 
1002  end do
1003 
1004  ! Solve the tri-diagonal system in k-direction.
1005 
1006  call tdia3(2_inttype, kl, bb, cc, dd, ff)
1007 
1008  ! Store the update in dvt.
1009 
1010  do k = 2, kl
1011  dvt(i, j, k, 1) = ff(1, k)
1012  dvt(i, j, k, 2) = ff(2, k)
1013  end do
1014 
1015  end do
1016  end do
1017  !
1018  ! Update the turbulent variables. For explicit relaxation the
1019  ! update must be relaxed; for implicit relaxation this has been
1020  ! done via the time step.
1021  !
1022  factor = one
1023  if (turbrelax == turbrelaxexplicit) factor = alfaturb
1024 
1025  do k = 2, kl
1026  do j = 2, jl
1027  do i = 2, il
1028  w(i, j, k, itu1) = w(i, j, k, itu1) + factor * dvt(i, j, k, 1)
1029  w(i, j, k, itu1) = max(w(i, j, k, itu1), zero)
1030 
1031  w(i, j, k, itu2) = w(i, j, k, itu2) + factor * dvt(i, j, k, 2)
1032  w(i, j, k, itu2) = max(w(i, j, k, itu2), 1.e-5_realtype * winf(itu2))
1033  end do
1034  end do
1035  end do
1036 
1037  end subroutine kwsolve
1038 end module kw
bcdata
Definition: BCData.F90:1
blockpointers::bmtk2
real(kind=realtype), dimension(:, :, :, :), pointer bmtk2
Definition: blockPointers.F90:162
blockpointers::sfacek
real(kind=realtype), dimension(:, :, :), pointer sfacek
Definition: blockPointers.F90:119
blockpointers::addgridvelocities
logical addgridvelocities
Definition: blockPointers.F90:117
blockpointers::jl
integer(kind=inttype) jl
Definition: blockPointers.F90:38
blockpointers::bmti1
real(kind=realtype), dimension(:, :, :, :), pointer bmti1
Definition: blockPointers.F90:157
blockpointers::nviscbocos
integer(kind=inttype) nviscbocos
Definition: blockPointers.F90:47
blockpointers::bmtj1
real(kind=realtype), dimension(:, :, :, :), pointer bmtj1
Definition: blockPointers.F90:159
blockpointers::bmti2
real(kind=realtype), dimension(:, :, :, :), pointer bmti2
Definition: blockPointers.F90:158
blockpointers::w
real(kind=realtype), dimension(:, :, :, :), pointer w
Definition: blockPointers.F90:120
blockpointers::scratch
real(kind=realtype), dimension(:, :, :, :), pointer scratch
Definition: blockPointers.F90:129
blockpointers::sfacei
real(kind=realtype), dimension(:, :, :), pointer sfacei
Definition: blockPointers.F90:119
blockpointers::viscsubface
type(viscsubfacetype), dimension(:), pointer viscsubface
Definition: blockPointers.F90:85
blockpointers::d2wall
real(kind=realtype), dimension(:, :, :), pointer d2wall
Definition: blockPointers.F90:154
blockpointers::iblank
integer(kind=inttype), dimension(:, :, :), pointer iblank
Definition: blockPointers.F90:71
blockpointers::rlv
real(kind=realtype), dimension(:, :, :), pointer rlv
Definition: blockPointers.F90:125
blockpointers::bcfaceid
integer(kind=inttype), dimension(:), pointer bcfaceid
Definition: blockPointers.F90:50
blockpointers::si
real(kind=realtype), dimension(:, :, :, :), pointer si
Definition: blockPointers.F90:96
blockpointers::volref
real(kind=realtype), dimension(:, :, :), pointer volref
Definition: blockPointers.F90:98
blockpointers::sj
real(kind=realtype), dimension(:, :, :, :), pointer sj
Definition: blockPointers.F90:96
blockpointers::rev
real(kind=realtype), dimension(:, :, :), pointer rev
Definition: blockPointers.F90:125
blockpointers::bmtj2
real(kind=realtype), dimension(:, :, :, :), pointer bmtj2
Definition: blockPointers.F90:160
blockpointers::dw
real(kind=realtype), dimension(:, :, :, :), pointer dw
Definition: blockPointers.F90:128
blockpointers::sk
real(kind=realtype), dimension(:, :, :, :), pointer sk
Definition: blockPointers.F90:96
blockpointers::vol
real(kind=realtype), dimension(:, :, :), pointer vol
Definition: blockPointers.F90:97
blockpointers::sfacej
real(kind=realtype), dimension(:, :, :), pointer sfacej
Definition: blockPointers.F90:119
blockpointers::bmtk1
real(kind=realtype), dimension(:, :, :, :), pointer bmtk1
Definition: blockPointers.F90:161
blockpointers::kl
integer(kind=inttype) kl
Definition: blockPointers.F90:38
blockpointers::il
integer(kind=inttype) il
Definition: blockPointers.F90:38
constants::strain
integer(kind=inttype), parameter strain
Definition: constants.F90:137
constants::zero
real(kind=realtype), parameter zero
Definition: constants.F90:71
constants::imax
integer(kind=inttype), parameter imax
Definition: constants.F90:293
constants::idvt
integer, parameter idvt
Definition: constants.F90:51
constants::kmin
integer(kind=inttype), parameter kmin
Definition: constants.F90:296
constants::jmax
integer(kind=inttype), parameter jmax
Definition: constants.F90:295
constants::itu2
integer, parameter itu2
Definition: constants.F90:42
constants::icd
integer, parameter icd
Definition: constants.F90:56
constants::itu1
integer, parameter itu1
Definition: constants.F90:40
constants::irho
integer, parameter irho
Definition: constants.F90:34
constants::vorticity
integer(kind=inttype), parameter vorticity
Definition: constants.F90:137
constants::komegamodified
integer(kind=inttype), parameter komegamodified
Definition: constants.F90:128
constants::ivx
integer, parameter ivx
Definition: constants.F90:35
constants::iprod
integer, parameter iprod
Definition: constants.F90:55
constants::one
real(kind=realtype), parameter one
Definition: constants.F90:72
constants::half
real(kind=realtype), parameter half
Definition: constants.F90:80
constants::turbrelaximplicit
integer(kind=inttype), parameter turbrelaximplicit
Definition: constants.F90:223
constants::turbrelaxexplicit
integer(kind=inttype), parameter turbrelaxexplicit
Definition: constants.F90:223
constants::katolaunder
integer(kind=inttype), parameter katolaunder
Definition: constants.F90:137
constants::imin
integer(kind=inttype), parameter imin
Definition: constants.F90:292
constants::ivz
integer, parameter ivz
Definition: constants.F90:37
constants::two
real(kind=realtype), parameter two
Definition: constants.F90:73
constants::kmax
integer(kind=inttype), parameter kmax
Definition: constants.F90:297
constants::ivy
integer, parameter ivy
Definition: constants.F90:36
constants::jmin
integer(kind=inttype), parameter jmin
Definition: constants.F90:294
flowvarrefstate::winf
real(kind=realtype), dimension(:), allocatable winf
Definition: flowVarRefState.F90:74
inputiteration::alfaturb
real(kind=realtype) alfaturb
Definition: inputParam.F90:276
inputiteration::turbrelax
integer(kind=inttype) turbrelax
Definition: inputParam.F90:269
inputphysics::turbprod
integer(kind=inttype) turbprod
Definition: inputParam.F90:586
inputphysics::turbmodel
integer(kind=inttype) turbmodel
Definition: inputParam.F90:586
inputphysics::wallfunctions
logical wallfunctions
Definition: inputParam.F90:591
inputphysics::pklim
real(kind=realtype) pklim
Definition: inputParam.F90:598
kw
Definition: kw.F90:1
kw::kw_block
subroutine kw_block(resOnly)
Definition: kw.F90:6
kw::kwsolve
subroutine kwsolve(resOnly)
Definition: kw.F90:52
paramturb::rkwsigd1
real(kind=realtype), parameter rkwsigd1
Definition: paramTurb.F90:27
paramturb::rkwbetas
real(kind=realtype), parameter rkwbetas
Definition: paramTurb.F90:29
paramturb::rkwsigk1
real(kind=realtype), parameter rkwsigk1
Definition: paramTurb.F90:25
paramturb::rkwk
real(kind=realtype), parameter rkwk
Definition: paramTurb.F90:24
paramturb::rkwsigw1
real(kind=realtype), parameter rkwsigw1
Definition: paramTurb.F90:26
paramturb::rkwbeta1
real(kind=realtype), parameter rkwbeta1
Definition: paramTurb.F90:28
turbbcroutines::bcturbtreatment
subroutine bcturbtreatment
Definition: turbBCRoutines.F90:663
turbbcroutines::applyallturbbcthisblock
subroutine applyallturbbcthisblock(secondHalo)
Definition: turbBCRoutines.F90:50
turbcurvefits::curvetupyp
subroutine curvetupyp(tup, yp, ntu1, ntu2)
Definition: turbCurveFits.F90:105
turbmod::vort
real(kind=realtype), dimension(:, :, :), pointer vort
Definition: turbMod.F90:34
turbmod::sig2
real(kind=realtype) sig2
Definition: turbMod.F90:16
turbmod::kwcd
real(kind=realtype), dimension(:, :, :), pointer kwcd
Definition: turbMod.F90:37
turbmod::prod
real(kind=realtype), dimension(:, :, :), pointer prod
Definition: turbMod.F90:35
turbmod::sig1
real(kind=realtype) sig1
Definition: turbMod.F90:16
turbmod::dvt
real(kind=realtype), dimension(:, :, :, :), pointer dvt
Definition: turbMod.F90:33
turbmod::sct
real(kind=realtype), dimension(:, :, :), pointer sct
Definition: turbMod.F90:39
turbutils::unsteadyturbterm
subroutine unsteadyturbterm(mAdv, nAdv, offset, qq)
Definition: turbUtils.F90:412
turbutils::kwcdterm
subroutine kwcdterm
Definition: turbUtils.F90:2075
turbutils::prodwmag2
subroutine prodwmag2
Definition: turbUtils.F90:240
turbutils::tdia3
subroutine tdia3(nb, ne, l, c, u, r)
Definition: turbUtils.F90:1564
turbutils::turbadvection
subroutine turbadvection(mAdv, nAdv, offset, qq)
Definition: turbUtils.F90:829
turbutils::prodsmag2
subroutine prodsmag2
Definition: turbUtils.F90:130
turbutils::prodkatolaunder
subroutine prodkatolaunder
Definition: turbUtils.F90:6
turbutils::kweddyviscosity
subroutine kweddyviscosity(iBeg, iEnd, jBeg, jEnd, kBeg, kEnd)
Definition: turbUtils.F90:715
utils
Definition: utils.F90:1
utils::setpointers
subroutine setpointers(nn, mm, ll)
Definition: utils.F90:3237
Generated by doxygen 1.9.1