Neko  0.9.0
A portable framework for high-order spectral element flow simulations
fluid_stats.f90
Go to the documentation of this file.
1 ! Copyright (c) 2022-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 !
37  use mean_field, only : mean_field_t
40  use num_types, only : rp
41  use math, only : invers2, col2, addcol3, col3, copy, subcol3
42  use operators, only : opgrad
43  use coefs, only : coef_t
44  use field, only : field_t
45  use field_list, only : field_list_t
46  use stats_quant, only : stats_quant_t
47  use device
48  use neko_config, only : neko_bcknd_device
49  use utils, only : neko_warning
50  implicit none
51  private
52 
53  type, public, extends(stats_quant_t) :: fluid_stats_t
55  type(field_t) :: stats_u
56  type(field_t) :: stats_v
57  type(field_t) :: stats_w
58  type(field_t) :: stats_p
59  type(field_t) :: stats_work
60 
62  type(field_t), pointer :: u
63  type(field_t), pointer :: v
64  type(field_t), pointer :: w
65  type(field_t), pointer :: p
66 
67  type(mean_field_t) :: u_mean
68  type(mean_field_t) :: v_mean
69  type(mean_field_t) :: w_mean
70  type(mean_field_t) :: p_mean
72  type(mean_field_t) :: uu
73  type(mean_field_t) :: vv
74  type(mean_field_t) :: ww
75  type(mean_field_t) :: uv
76  type(mean_field_t) :: uw
77  type(mean_field_t) :: vw
79  type(mean_field_t) :: uuu
80  type(mean_field_t) :: vvv
81  type(mean_field_t) :: www
82  type(mean_field_t) :: uuv
83  type(mean_field_t) :: uuw
84  type(mean_field_t) :: uvv
85  type(mean_field_t) :: uvw
86  type(mean_field_t) :: vvw
87  type(mean_field_t) :: uww
88  type(mean_field_t) :: vww
90  type(mean_field_t) :: uuuu
91  type(mean_field_t) :: vvvv
92  type(mean_field_t) :: wwww
94  type(mean_field_t) :: pp
95  type(mean_field_t) :: ppp
96  type(mean_field_t) :: pppp
98  type(mean_field_t) :: pu
99  type(mean_field_t) :: pv
100  type(mean_field_t) :: pw
101 
103  type(mean_field_t) :: pdudx
104  type(mean_field_t) :: pdudy
105  type(mean_field_t) :: pdudz
106  type(mean_field_t) :: pdvdx
107  type(mean_field_t) :: pdvdy
108  type(mean_field_t) :: pdvdz
109  type(mean_field_t) :: pdwdx
110  type(mean_field_t) :: pdwdy
111  type(mean_field_t) :: pdwdz
112 
114  type(mean_field_t) :: e11
115  type(mean_field_t) :: e22
116  type(mean_field_t) :: e33
117  type(mean_field_t) :: e12
118  type(mean_field_t) :: e13
119  type(mean_field_t) :: e23
121  type(field_t) :: dudx
122  type(field_t) :: dudy
123  type(field_t) :: dudz
124  type(field_t) :: dvdx
125  type(field_t) :: dvdy
126  type(field_t) :: dvdz
127  type(field_t) :: dwdx
128  type(field_t) :: dwdy
129  type(field_t) :: dwdz
130 
132  type(coef_t), pointer :: coef
134  integer :: n_stats = 44
137  character(5) :: stat_set
140  type(field_list_t) :: stat_fields
141  contains
143  procedure, pass(this) :: init => fluid_stats_init
145  procedure, pass(this) :: free => fluid_stats_free
147  procedure, pass(this) :: update => fluid_stats_update
149  procedure, pass(this) :: reset => fluid_stats_reset
150  ! Convert computed weak gradients to strong.
151  procedure, pass(this) :: make_strong_grad => fluid_stats_make_strong_grad
154  procedure, pass(this) :: post_process => fluid_stats_post_process
155  end type fluid_stats_t
156 
157 contains
158 
167  subroutine fluid_stats_init(this, coef, u, v, w, p, set)
168  class(fluid_stats_t), intent(inout), target:: this
169  type(coef_t), target, optional :: coef
170  type(field_t), target, intent(inout) :: u, v, w, p
171  character(*), intent(in), optional :: set
172 
173  call this%free()
174  this%coef => coef
175 
176  this%u => u
177  this%v => v
178  this%w => w
179  this%p => p
180 
181  if (present(set)) then
182  this%stat_set = trim(set)
183  if (this%stat_set .eq. 'basic') then
184  this%n_stats = 11
185  end if
186  else
187  this%stat_set = 'full'
188  this%n_stats = 44
189  end if
190 
191  call this%stats_work%init(this%u%dof, 'stats')
192  call this%stats_u%init(this%u%dof, 'u temp')
193  call this%stats_v%init(this%u%dof, 'v temp')
194  call this%stats_w%init(this%u%dof, 'w temp')
195  call this%stats_p%init(this%u%dof, 'p temp')
196  call this%u_mean%init(this%u)
197  call this%v_mean%init(this%v)
198  call this%w_mean%init(this%w)
199  call this%p_mean%init(this%p)
200  call this%uu%init(this%stats_u, 'uu')
201  call this%vv%init(this%stats_v, 'vv')
202  call this%ww%init(this%stats_w, 'ww')
203  call this%uv%init(this%stats_work, 'uv')
204  call this%uw%init(this%stats_work, 'uw')
205  call this%vw%init(this%stats_work, 'vw')
206  call this%pp%init(this%stats_p, 'pp')
207 
208  if (this%n_stats .eq. 44) then
209  call this%dudx%init(this%u%dof, 'dudx')
210  call this%dudy%init(this%u%dof, 'dudy')
211  call this%dudz%init(this%u%dof, 'dudz')
212  call this%dvdx%init(this%u%dof, 'dvdx')
213  call this%dvdy%init(this%u%dof, 'dvdy')
214  call this%dvdz%init(this%u%dof, 'dvdz')
215  call this%dwdx%init(this%u%dof, 'dwdx')
216  call this%dwdy%init(this%u%dof, 'dwdy')
217  call this%dwdz%init(this%u%dof, 'dwdz')
218 
219  call this%uuu%init(this%stats_work, 'uuu')
220  call this%vvv%init(this%stats_work, 'vvv')
221  call this%www%init(this%stats_work, 'www')
222  call this%uuv%init(this%stats_work, 'uuv')
223  call this%uuw%init(this%stats_work, 'uuw')
224  call this%uvv%init(this%stats_work, 'uvv')
225  call this%uvw%init(this%stats_work, 'uvw')
226  call this%vvw%init(this%stats_work, 'vvw')
227  call this%uww%init(this%stats_work, 'uww')
228  call this%vww%init(this%stats_work, 'vww')
229  call this%uuuu%init(this%stats_work, 'uuuu')
230  call this%vvvv%init(this%stats_work, 'vvvv')
231  call this%wwww%init(this%stats_work, 'wwww')
233  call this%ppp%init(this%stats_work, 'ppp')
234  call this%pppp%init(this%stats_work, 'pppp')
236  call this%pu%init(this%stats_work, 'pu')
237  call this%pv%init(this%stats_work, 'pv')
238  call this%pw%init(this%stats_work, 'pw')
239 
240  call this%pdudx%init(this%stats_work, 'pdudx')
241  call this%pdudy%init(this%stats_work, 'pdudy')
242  call this%pdudz%init(this%stats_work, 'pdudz')
243  call this%pdvdx%init(this%stats_work, 'pdvdx')
244  call this%pdvdy%init(this%stats_work, 'pdvdy')
245  call this%pdvdz%init(this%stats_work, 'pdvdz')
246  call this%pdwdx%init(this%stats_work, 'pdwdx')
247  call this%pdwdy%init(this%stats_work, 'pdwdy')
248  call this%pdwdz%init(this%stats_work, 'pdwdz')
249 
250  call this%e11%init(this%stats_work, 'e11')
251  call this%e22%init(this%stats_work, 'e22')
252  call this%e33%init(this%stats_work, 'e33')
253  call this%e12%init(this%stats_work, 'e12')
254  call this%e13%init(this%stats_work, 'e13')
255  call this%e23%init(this%stats_work, 'e23')
256  end if
257 
258  allocate(this%stat_fields%items(this%n_stats))
259 
260  call this%stat_fields%assign_to_field(1, this%p_mean%mf)
261  call this%stat_fields%assign_to_field(2, this%u_mean%mf)
262  call this%stat_fields%assign_to_field(3, this%v_mean%mf)
263  call this%stat_fields%assign_to_field(4, this%w_mean%mf)
264  call this%stat_fields%assign_to_field(5, this%pp%mf)
265  call this%stat_fields%assign_to_field(6, this%uu%mf)
266  call this%stat_fields%assign_to_field(7, this%vv%mf)
267  call this%stat_fields%assign_to_field(8, this%ww%mf)
268  call this%stat_fields%assign_to_field(9, this%uv%mf)
269  call this%stat_fields%assign_to_field(10, this%uw%mf)
270  call this%stat_fields%assign_to_field(11, this%vw%mf)
271 
272  if (this%n_stats .eq. 44) then
273  call this%stat_fields%assign_to_field(12, this%uuu%mf)
274  call this%stat_fields%assign_to_field(13, this%vvv%mf)
275  call this%stat_fields%assign_to_field(14, this%www%mf)
276  call this%stat_fields%assign_to_field(15, this%uuv%mf)
277  call this%stat_fields%assign_to_field(16, this%uuw%mf)
278  call this%stat_fields%assign_to_field(17, this%uvv%mf)
279  call this%stat_fields%assign_to_field(18, this%uvw%mf)
280  call this%stat_fields%assign_to_field(19, this%vvw%mf)
281  call this%stat_fields%assign_to_field(20, this%uww%mf)
282  call this%stat_fields%assign_to_field(21, this%vww%mf)
283  call this%stat_fields%assign_to_field(22, this%uuuu%mf)
284  call this%stat_fields%assign_to_field(23, this%vvvv%mf)
285  call this%stat_fields%assign_to_field(24, this%wwww%mf)
286  call this%stat_fields%assign_to_field(25, this%ppp%mf)
287  call this%stat_fields%assign_to_field(26, this%pppp%mf)
288  call this%stat_fields%assign_to_field(27, this%pu%mf)
289  call this%stat_fields%assign_to_field(28, this%pv%mf)
290  call this%stat_fields%assign_to_field(29, this%pw%mf)
291 
292  call this%stat_fields%assign_to_field(30, this%pdudx%mf)
293  call this%stat_fields%assign_to_field(31, this%pdudy%mf)
294  call this%stat_fields%assign_to_field(32, this%pdudz%mf)
295  call this%stat_fields%assign_to_field(33, this%pdvdx%mf)
296  call this%stat_fields%assign_to_field(34, this%pdvdy%mf)
297  call this%stat_fields%assign_to_field(35, this%pdvdz%mf)
298  call this%stat_fields%assign_to_field(36, this%pdwdx%mf)
299  call this%stat_fields%assign_to_field(37, this%pdwdy%mf)
300  call this%stat_fields%assign_to_field(38, this%pdwdz%mf)
301  call this%stat_fields%assign_to_field(39, this%e11%mf)
302  call this%stat_fields%assign_to_field(40, this%e22%mf)
303  call this%stat_fields%assign_to_field(41, this%e33%mf)
304  call this%stat_fields%assign_to_field(42, this%e12%mf)
305  call this%stat_fields%assign_to_field(43, this%e13%mf)
306  call this%stat_fields%assign_to_field(44, this%e23%mf)
307  end if
308 
309  end subroutine fluid_stats_init
310 
313  subroutine fluid_stats_update(this, k)
314  class(fluid_stats_t), intent(inout) :: this
315  real(kind=rp), intent(in) :: k
316  integer :: n
317 
318  associate(stats_work => this%stats_work, stats_u => this%stats_u, &
319  stats_v => this%stats_v, stats_w => this%stats_w, &
320  stats_p => this%stats_p)
321  n = stats_work%dof%size()
322 
324  if (neko_bcknd_device .eq. 1) then
325 
326  call this%u_mean%update(k)
327  call this%v_mean%update(k)
328  call this%w_mean%update(k)
329  call this%p_mean%update(k)
330 
331  call device_col3(stats_u%x_d, this%u%x_d, this%u%x_d, n)
332  call device_col3(stats_v%x_d, this%v%x_d, this%v%x_d, n)
333  call device_col3(stats_w%x_d, this%w%x_d, this%w%x_d, n)
334  call device_col3(stats_p%x_d, this%p%x_d, this%p%x_d, n)
335 
336  call this%uu%update(k)
337  call this%vv%update(k)
338  call this%ww%update(k)
339  call this%pp%update(k)
340 
341  call device_col3(stats_work%x_d, this%u%x_d, this%v%x_d, n)
342  call this%uv%update(k)
343  call device_col3(stats_work%x_d, this%u%x_d, this%w%x_d, n)
344  call this%uw%update(k)
345  call device_col3(stats_work%x_d, this%v%x_d, this%w%x_d, n)
346  call this%vw%update(k)
347  if (this%n_stats .eq. 11) return
348  call device_col2(stats_work%x_d, this%u%x_d, n)
349  call this%uvw%update(k)
350  call device_col3(stats_work%x_d, this%stats_u%x_d, this%u%x_d, n)
351  call this%uuu%update(k)
352  call device_col3(stats_work%x_d, this%stats_v%x_d, this%v%x_d, n)
353  call this%vvv%update(k)
354  call device_col3(stats_work%x_d, this%stats_w%x_d, this%w%x_d, n)
355  call this%www%update(k)
356  call device_col3(stats_work%x_d, this%stats_u%x_d, this%v%x_d, n)
357  call this%uuv%update(k)
358  call device_col3(stats_work%x_d, this%stats_u%x_d, this%w%x_d, n)
359  call this%uuw%update(k)
360  call device_col3(stats_work%x_d, this%stats_v%x_d, this%u%x_d, n)
361  call this%uvv%update(k)
362  call device_col3(stats_work%x_d, this%stats_v%x_d, this%w%x_d, n)
363  call this%vvw%update(k)
364  call device_col3(stats_work%x_d, this%stats_w%x_d, this%u%x_d, n)
365  call this%uww%update(k)
366  call device_col3(stats_work%x_d, this%stats_w%x_d, this%v%x_d, n)
367  call this%vww%update(k)
368 
369  call device_col3(stats_work%x_d, this%stats_u%x_d, this%stats_u%x_d, n)
370  call this%uuuu%update(k)
371  call device_col3(stats_work%x_d, this%stats_v%x_d, this%stats_v%x_d, n)
372  call this%vvvv%update(k)
373  call device_col3(stats_work%x_d, this%stats_w%x_d, this%stats_w%x_d, n)
374  call this%wwww%update(k)
375 
376  call device_col3(stats_work%x_d, this%stats_p%x_d, this%p%x_d, n)
377  call this%ppp%update(k)
378  call device_col3(stats_work%x_d, this%stats_p%x_d, this%stats_p%x_d, n)
379  call this%pppp%update(k)
380 
381  call device_col3(stats_work%x_d, this%p%x_d, this%u%x_d, n)
382  call this%pu%update(k)
383  call device_col3(stats_work%x_d, this%p%x_d, this%v%x_d, n)
384  call this%pv%update(k)
385  call device_col3(stats_work%x_d, this%p%x_d, this%w%x_d, n)
386  call this%pw%update(k)
387 
388  else
389 
390  call this%u_mean%update(k)
391  call this%v_mean%update(k)
392  call this%w_mean%update(k)
393  call this%p_mean%update(k)
394  call col3(stats_u%x, this%u%x, this%u%x, n)
395  call col3(stats_v%x, this%v%x, this%v%x, n)
396  call col3(stats_w%x, this%w%x, this%w%x, n)
397  call col3(stats_p%x, this%p%x, this%p%x, n)
398 
399  call this%uu%update(k)
400  call this%vv%update(k)
401  call this%ww%update(k)
402  call this%pp%update(k)
403 
404  call col3(stats_work%x, this%u%x, this%v%x, n)
405  call this%uv%update(k)
406  call col3(stats_work%x, this%u%x, this%w%x, n)
407  call this%uw%update(k)
408  call col3(stats_work%x, this%v%x, this%w%x, n)
409  call this%vw%update(k)
410 
411  if (this%n_stats .eq. 11) return
412 
413  call col2(stats_work%x, this%u%x, n)
414  call this%uvw%update(k)
415  call col3(stats_work%x, this%stats_u%x, this%u%x, n)
416  call this%uuu%update(k)
417  call col3(stats_work%x, this%stats_v%x, this%v%x, n)
418  call this%vvv%update(k)
419  call col3(stats_work%x, this%stats_w%x, this%w%x, n)
420  call this%www%update(k)
421  call col3(stats_work%x, this%stats_u%x, this%v%x, n)
422  call this%uuv%update(k)
423  call col3(stats_work%x, this%stats_u%x, this%w%x, n)
424  call this%uuw%update(k)
425  call col3(stats_work%x, this%stats_v%x, this%u%x, n)
426  call this%uvv%update(k)
427  call col3(stats_work%x, this%stats_v%x, this%w%x, n)
428  call this%vvw%update(k)
429  call col3(stats_work%x, this%stats_w%x, this%u%x, n)
430  call this%uww%update(k)
431  call col3(stats_work%x, this%stats_w%x, this%v%x, n)
432  call this%vww%update(k)
433 
434  call col3(stats_work%x, this%stats_u%x, this%stats_u%x, n)
435  call this%uuuu%update(k)
436  call col3(stats_work%x, this%stats_v%x, this%stats_v%x, n)
437  call this%vvvv%update(k)
438  call col3(stats_work%x, this%stats_w%x, this%stats_w%x, n)
439  call this%wwww%update(k)
440 
441  call col3(stats_work%x, this%stats_p%x, this%p%x, n)
442  call this%ppp%update(k)
443  call col3(stats_work%x, this%stats_p%x, this%stats_p%x, n)
444  call this%pppp%update(k)
445 
446  call col3(stats_work%x, this%p%x, this%u%x,n)
447  call this%pu%update(k)
448  call col3(stats_work%x, this%p%x, this%v%x,n)
449  call this%pv%update(k)
450  call col3(stats_work%x, this%p%x, this%w%x,n)
451  call this%pw%update(k)
452 
453 
454  end if
455  call opgrad(this%dudx%x, this%dudy%x, this%dudz%x, this%u%x, this%coef)
456  call opgrad(this%dvdx%x, this%dvdy%x, this%dvdz%x, this%v%x, this%coef)
457  call opgrad(this%dwdx%x, this%dwdy%x, this%dwdz%x, this%w%x, this%coef)
458 
459  if (neko_bcknd_device .eq. 1) then
460  call device_col3(stats_work%x_d, this%dudx%x_d, this%p%x_d, n)
461  call this%pdudx%update(k)
462  call device_col3(stats_work%x_d, this%dudy%x_d, this%p%x_d, n)
463  call this%pdudy%update(k)
464  call device_col3(stats_work%x_d, this%dudz%x_d, this%p%x_d, n)
465  call this%pdudz%update(k)
466 
467  call device_col3(stats_work%x_d, this%dvdx%x_d, this%p%x_d, n)
468  call this%pdvdx%update(k)
469  call device_col3(stats_work%x_d, this%dvdy%x_d, this%p%x_d, n)
470  call this%pdvdy%update(k)
471  call device_col3(stats_work%x_d, this%dvdz%x_d, this%p%x_d, n)
472  call this%pdvdz%update(k)
473 
474  call device_col3(stats_work%x_d, this%dwdx%x_d, this%p%x_d, n)
475  call this%pdwdx%update(k)
476  call device_col3(stats_work%x_d, this%dwdy%x_d, this%p%x_d, n)
477  call this%pdwdy%update(k)
478  call device_col3(stats_work%x_d, this%dwdz%x_d, this%p%x_d, n)
479  call this%pdwdz%update(k)
480 
481  call device_col3(this%stats_work%x_d, this%dudx%x_d, this%dudx%x_d, n)
482  call device_addcol3(this%stats_work%x_d, this%dudy%x_d, &
483  this%dudy%x_d, n)
484  call device_addcol3(this%stats_work%x_d, this%dudz%x_d, &
485  this%dudz%x_d, n)
486  call this%e11%update(k)
487  call device_col3(this%stats_work%x_d, this%dvdx%x_d, this%dvdx%x_d, n)
488  call device_addcol3(this%stats_work%x_d, this%dvdy%x_d, &
489  this%dvdy%x_d, n)
490  call device_addcol3(this%stats_work%x_d, this%dvdz%x_d, &
491  this%dvdz%x_d, n)
492  call this%e22%update(k)
493  call device_col3(this%stats_work%x_d, this%dwdx%x_d, this%dwdx%x_d, n)
494  call device_addcol3(this%stats_work%x_d, this%dwdy%x_d, &
495  this%dwdy%x_d, n)
496  call device_addcol3(this%stats_work%x_d, this%dwdz%x_d, &
497  this%dwdz%x_d, n)
498  call this%e33%update(k)
499  call device_col3(this%stats_work%x_d, this%dudx%x_d, &
500  this%dvdx%x_d, n)
501  call device_addcol3(this%stats_work%x_d, this%dudy%x_d, &
502  this%dvdy%x_d, n)
503  call device_addcol3(this%stats_work%x_d, this%dudz%x_d, &
504  this%dvdz%x_d, n)
505  call this%e12%update(k)
506  call device_col3(this%stats_work%x_d, this%dudx%x_d, this%dwdx%x_d, n)
507  call device_addcol3(this%stats_work%x_d, this%dudy%x_d, &
508  this%dwdy%x_d, n)
509  call device_addcol3(this%stats_work%x_d, this%dudz%x_d, &
510  this%dwdz%x_d, n)
511  call this%e13%update(k)
512  call device_col3(this%stats_work%x_d, this%dvdx%x_d, this%dwdx%x_d, n)
513  call device_addcol3(this%stats_work%x_d, this%dvdy%x_d, &
514  this%dwdy%x_d, n)
515  call device_addcol3(this%stats_work%x_d, this%dvdz%x_d, &
516  this%dwdz%x_d, n)
517  call this%e23%update(k)
518  else
519  call col3(stats_work%x, this%dudx%x, this%p%x, n)
520  call this%pdudx%update(k)
521  call col3(stats_work%x, this%dudy%x, this%p%x, n)
522  call this%pdudy%update(k)
523  call col3(stats_work%x, this%dudz%x, this%p%x, n)
524  call this%pdudz%update(k)
525 
526  call col3(stats_work%x, this%dvdx%x, this%p%x, n)
527  call this%pdvdx%update(k)
528  call col3(stats_work%x, this%dvdy%x, this%p%x, n)
529  call this%pdvdy%update(k)
530  call col3(stats_work%x, this%dvdz%x, this%p%x, n)
531  call this%pdvdz%update(k)
532 
533  call col3(stats_work%x, this%dwdx%x, this%p%x, n)
534  call this%pdwdx%update(k)
535  call col3(stats_work%x, this%dwdy%x, this%p%x, n)
536  call this%pdwdy%update(k)
537  call col3(stats_work%x, this%dwdz%x, this%p%x, n)
538  call this%pdwdz%update(k)
539 
540  call col3(this%stats_work%x, this%dudx%x, this%dudx%x, n)
541  call addcol3(this%stats_work%x, this%dudy%x, this%dudy%x, n)
542  call addcol3(this%stats_work%x, this%dudz%x, this%dudz%x, n)
543  call this%e11%update(k)
544  call col3(this%stats_work%x, this%dvdx%x, this%dvdx%x, n)
545  call addcol3(this%stats_work%x, this%dvdy%x, this%dvdy%x, n)
546  call addcol3(this%stats_work%x, this%dvdz%x, this%dvdz%x, n)
547  call this%e22%update(k)
548  call col3(this%stats_work%x, this%dwdx%x, this%dwdx%x, n)
549  call addcol3(this%stats_work%x, this%dwdy%x, this%dwdy%x, n)
550  call addcol3(this%stats_work%x, this%dwdz%x, this%dwdz%x, n)
551  call this%e33%update(k)
552  call col3(this%stats_work%x, this%dudx%x, this%dvdx%x, n)
553  call addcol3(this%stats_work%x, this%dudy%x, this%dvdy%x, n)
554  call addcol3(this%stats_work%x, this%dudz%x, this%dvdz%x, n)
555  call this%e12%update(k)
556  call col3(this%stats_work%x,this%dudx%x, this%dwdx%x,n)
557  call addcol3(this%stats_work%x,this%dudy%x, this%dwdy%x,n)
558  call addcol3(this%stats_work%x,this%dudz%x, this%dwdz%x,n)
559  call this%e13%update(k)
560  call col3(this%stats_work%x, this%dvdx%x, this%dwdx%x, n)
561  call addcol3(this%stats_work%x, this%dvdy%x, this%dwdy%x, n)
562  call addcol3(this%stats_work%x, this%dvdz%x, this%dwdz%x, n)
563  call this%e23%update(k)
564 
565  end if
566  end associate
567 
568  end subroutine fluid_stats_update
569 
570 
572  subroutine fluid_stats_free(this)
573  class(fluid_stats_t), intent(inout) :: this
574 
575  call this%stats_work%free()
576  call this%stats_u%free()
577  call this%stats_v%free()
578  call this%stats_w%free()
579 
580  call this%u_mean%free()
581  call this%v_mean%free()
582  call this%w_mean%free()
583  call this%p_mean%free()
584 
585  call this%uu%free()
586  call this%vv%free()
587  call this%ww%free()
588  call this%uv%free()
589  call this%uw%free()
590  call this%vw%free()
591  call this%pp%free()
592 
593  call this%dUdx%free()
594  call this%dUdy%free()
595  call this%dUdz%free()
596  call this%dVdx%free()
597  call this%dVdy%free()
598  call this%dVdz%free()
599  call this%dWdx%free()
600  call this%dWdy%free()
601  call this%dWdz%free()
602 
603  end subroutine fluid_stats_free
604 
606  subroutine fluid_stats_reset(this)
607  class(fluid_stats_t), intent(inout), target:: this
608 
609  call this%p_mean%reset()
610  call this%u_mean%reset()
611  call this%v_mean%reset()
612  call this%w_mean%reset()
613 
614  call this%uu%reset()
615  call this%vv%reset()
616  call this%ww%reset()
617  call this%uv%reset()
618  call this%uw%reset()
619  call this%vw%reset()
620  call this%pp%reset()
621  if (this%n_stats .eq. 44) then
622  call this%uuu%reset()
623  call this%vvv%reset()
624  call this%www%reset()
625  call this%uuv%reset()
626  call this%uuw%reset()
627  call this%uvv%reset()
628  call this%uvw%reset()
629  call this%vvw%reset()
630  call this%uww%reset()
631  call this%vww%reset()
632  call this%uuuu%reset()
633  call this%vvvv%reset()
634  call this%wwww%reset()
635  call this%ppp%reset()
636  call this%pppp%reset()
637  call this%pu%reset()
638  call this%pv%reset()
639  call this%pw%reset()
640 
641  call this%pdudx%reset()
642  call this%pdudy%reset()
643  call this%pdudz%reset()
644  call this%pdvdx%reset()
645  call this%pdvdy%reset()
646  call this%pdvdz%reset()
647  call this%pdwdx%reset()
648  call this%pdwdy%reset()
649  call this%pdwdz%reset()
650 
651  call this%e11%reset()
652  call this%e22%reset()
653  call this%e33%reset()
654  call this%e12%reset()
655  call this%e13%reset()
656  call this%e23%reset()
657  end if
658 
659  end subroutine fluid_stats_reset
660 
661  ! Convert computed weak gradients to strong.
663  class(fluid_stats_t) :: this
664  integer :: n
665 
666  if (this%n_stats .eq. 11) return
667 
668  n = size(this%coef%B)
669 
670  if (neko_bcknd_device .eq. 1) then
671  call device_cfill(this%stats_work%x_d, 1.0_rp, n)
672  call device_invcol2(this%stats_work%x_d, this%coef%B_d, n)
673  call device_col2(this%pdudx%mf%x_d, this%stats_work%x_d, n)
674  call device_col2(this%pdudy%mf%x_d, this%stats_work%x_d, n)
675  call device_col2(this%pdudz%mf%x_d, this%stats_work%x_d, n)
676  call device_col2(this%pdvdx%mf%x_d, this%stats_work%x_d, n)
677  call device_col2(this%pdvdy%mf%x_d, this%stats_work%x_d, n)
678  call device_col2(this%pdvdz%mf%x_d, this%stats_work%x_d, n)
679  call device_col2(this%pdwdx%mf%x_d, this%stats_work%x_d, n)
680  call device_col2(this%pdwdy%mf%x_d, this%stats_work%x_d, n)
681  call device_col2(this%pdwdz%mf%x_d, this%stats_work%x_d, n)
682 
683  call device_col2(this%stats_work%x_d, this%stats_work%x_d, n)
684  call device_col2(this%e11%mf%x_d, this%stats_work%x_d, n)
685  call device_col2(this%e22%mf%x_d, this%stats_work%x_d, n)
686  call device_col2(this%e33%mf%x_d, this%stats_work%x_d, n)
687  call device_col2(this%e12%mf%x_d, this%stats_work%x_d, n)
688  call device_col2(this%e13%mf%x_d, this%stats_work%x_d, n)
689  call device_col2(this%e23%mf%x_d, this%stats_work%x_d, n)
690 
691 
692  else
693  call invers2(this%stats_work%x, this%coef%B, n)
694  call col2(this%pdudx%mf%x, this%stats_work%x, n)
695  call col2(this%pdudy%mf%x, this%stats_work%x, n)
696  call col2(this%pdudz%mf%x, this%stats_work%x, n)
697  call col2(this%pdvdx%mf%x, this%stats_work%x, n)
698  call col2(this%pdvdy%mf%x, this%stats_work%x, n)
699  call col2(this%pdvdz%mf%x, this%stats_work%x, n)
700  call col2(this%pdwdx%mf%x, this%stats_work%x, n)
701  call col2(this%pdwdy%mf%x, this%stats_work%x, n)
702  call col2(this%pdwdz%mf%x, this%stats_work%x, n)
703 
704  call col2(this%stats_work%x, this%stats_work%x, n)
705  call col2(this%e11%mf%x, this%stats_work%x, n)
706  call col2(this%e22%mf%x, this%stats_work%x, n)
707  call col2(this%e33%mf%x, this%stats_work%x, n)
708  call col2(this%e12%mf%x, this%stats_work%x, n)
709  call col2(this%e13%mf%x, this%stats_work%x, n)
710  call col2(this%e23%mf%x, this%stats_work%x, n)
711 
712  end if
713 
714  end subroutine fluid_stats_make_strong_grad
715 
718  subroutine fluid_stats_post_process(this, mean, reynolds, pressure_flatness,&
719  pressure_skewness, skewness_tensor, mean_vel_grad, dissipation_tensor)
720  class(fluid_stats_t) :: this
721  type(field_list_t), intent(inout), optional :: mean
722  type(field_list_t), intent(inout), optional :: reynolds
723  type(field_list_t), intent(inout), optional :: pressure_skewness
724  type(field_list_t), intent(inout), optional :: pressure_flatness
725  type(field_list_t), intent(inout), optional :: skewness_tensor
726  type(field_list_t), intent(inout), optional :: mean_vel_grad
727  type(field_list_t), intent(inout), optional :: dissipation_tensor
728  integer :: n
729 
730  if (present(mean)) then
731  n = mean%item_size(1)
732  call copy(mean%items(1)%ptr%x, this%u_mean%mf%x, n)
733  call copy(mean%items(2)%ptr%x, this%v_mean%mf%x, n)
734  call copy(mean%items(3)%ptr%x, this%w_mean%mf%x, n)
735  call copy(mean%items(4)%ptr%x, this%p_mean%mf%x, n)
736  end if
737 
738  if (present(reynolds)) then
739  n = reynolds%item_size(1)
740  call copy(reynolds%items(1)%ptr%x, this%pp%mf%x, n)
741  call subcol3(reynolds%items(1)%ptr%x, this%p_mean%mf%x, &
742  this%p_mean%mf%x, n)
743 
744  call copy(reynolds%items(2)%ptr%x, this%uu%mf%x, n)
745  call subcol3(reynolds%items(2)%ptr%x, this%u_mean%mf%x, &
746  this%u_mean%mf%x, n)
747 
748  call copy(reynolds%items(3)%ptr%x, this%vv%mf%x, n)
749  call subcol3(reynolds%items(3)%ptr%x, this%v_mean%mf%x, &
750  this%v_mean%mf%x,n)
751 
752  call copy(reynolds%items(4)%ptr%x, this%ww%mf%x, n)
753  call subcol3(reynolds%items(4)%ptr%x, this%w_mean%mf%x, &
754  this%w_mean%mf%x,n)
755 
756  call copy(reynolds%items(5)%ptr%x, this%uv%mf%x, n)
757  call subcol3(reynolds%items(5)%ptr%x, this%u_mean%mf%x, &
758  this%v_mean%mf%x, n)
759 
760  call copy(reynolds%items(6)%ptr%x, this%uw%mf%x, n)
761  call subcol3(reynolds%items(6)%ptr%x, this%u_mean%mf%x, &
762  this%w_mean%mf%x, n)
763 
764  call copy(reynolds%items(7)%ptr%x, this%vw%mf%x, n)
765  call subcol3(reynolds%items(7)%ptr%x, this%v_mean%mf%x, &
766  this%w_mean%mf%x, n)
767  end if
768  if (present(pressure_skewness)) then
769 
770  call neko_warning('Presssure skewness stat not implemented'// &
771  ' in fluid_stats, process stats in python instead')
772 
773  end if
774 
775  if (present(pressure_flatness)) then
776  call neko_warning('Presssure flatness stat not implemented'// &
777  ' in fluid_stats, process stats in python instead')
778 
779  end if
780 
781  if (present(skewness_tensor)) then
782  call neko_warning('Skewness tensor stat not implemented'// &
783  ' in fluid_stats, process stats in python instead')
784  end if
785 
786  if (present(mean_vel_grad)) then
787  !Compute gradient of mean flow
788  n = mean_vel_grad%item_size(1)
789  if (neko_bcknd_device .eq. 1) then
790  call device_memcpy(this%u_mean%mf%x, this%u_mean%mf%x_d, n, &
791  host_to_device, sync = .false.)
792  call device_memcpy(this%v_mean%mf%x, this%v_mean%mf%x_d, n, &
793  host_to_device, sync = .false.)
794  call device_memcpy(this%w_mean%mf%x, this%w_mean%mf%x_d, n, &
795  host_to_device, sync = .false.)
796  call opgrad(this%dudx%x, this%dudy%x, this%dudz%x, &
797  this%u_mean%mf%x, this%coef)
798  call opgrad(this%dvdx%x, this%dvdy%x, this%dvdz%x, &
799  this%v_mean%mf%x, this%coef)
800  call opgrad(this%dwdx%x, this%dwdy%x, this%dwdz%x, &
801  this%w_mean%mf%x, this%coef)
802  call device_memcpy(this%dudx%x, this%dudx%x_d, n, &
803  device_to_host, sync = .false.)
804  call device_memcpy(this%dvdx%x, this%dvdx%x_d, n, &
805  device_to_host, sync = .false.)
806  call device_memcpy(this%dwdx%x, this%dwdx%x_d, n, &
807  device_to_host, sync = .false.)
808  call device_memcpy(this%dudy%x, this%dudy%x_d, n, &
809  device_to_host, sync = .false.)
810  call device_memcpy(this%dvdy%x, this%dvdy%x_d, n, &
811  device_to_host, sync = .false.)
812  call device_memcpy(this%dwdy%x, this%dwdy%x_d, n, &
813  device_to_host, sync = .false.)
814  call device_memcpy(this%dudz%x, this%dudz%x_d, n, &
815  device_to_host, sync = .false.)
816  call device_memcpy(this%dvdz%x, this%dvdz%x_d, n, &
817  device_to_host, sync = .false.)
818  call device_memcpy(this%dwdz%x, this%dwdz%x_d, n, &
819  device_to_host, sync = .true.)
820  else
821  call opgrad(this%dudx%x, this%dudy%x, this%dudz%x, &
822  this%u_mean%mf%x, this%coef)
823  call opgrad(this%dvdx%x, this%dvdy%x, this%dvdz%x, &
824  this%v_mean%mf%x, this%coef)
825  call opgrad(this%dwdx%x, this%dwdy%x, this%dwdz%x, &
826  this%w_mean%mf%x, this%coef)
827  end if
828  call invers2(this%stats_work%x, this%coef%B,n)
829  call col3(mean_vel_grad%items(1)%ptr%x, this%dudx%x, &
830  this%stats_work%x, n)
831  call col3(mean_vel_grad%items(2)%ptr%x, this%dudy%x, &
832  this%stats_work%x, n)
833  call col3(mean_vel_grad%items(3)%ptr%x, this%dudz%x, &
834  this%stats_work%x, n)
835  call col3(mean_vel_grad%items(4)%ptr%x, this%dvdx%x, &
836  this%stats_work%x, n)
837  call col3(mean_vel_grad%items(5)%ptr%x, this%dvdy%x, &
838  this%stats_work%x, n)
839  call col3(mean_vel_grad%items(6)%ptr%x, this%dvdz%x, &
840  this%stats_work%x, n)
841  call col3(mean_vel_grad%items(7)%ptr%x, this%dwdx%x, &
842  this%stats_work%x, n)
843  call col3(mean_vel_grad%items(8)%ptr%x, this%dwdy%x, &
844  this%stats_work%x, n)
845  call col3(mean_vel_grad%items(9)%ptr%x, this%dwdz%x, &
846  this%stats_work%x, n)
847 
848  end if
849 
850  if (present(dissipation_tensor)) then
851  call neko_warning('Dissipation tensor stat not implemented'// &
852  ' in fluid_stats, process stats in python instead')
853  end if
854 
855  end subroutine fluid_stats_post_process
856 
857 end module fluid_stats
Coefficients.
Definition: coef.f90:34
subroutine, public device_addcol3(a_d, b_d, c_d, n)
Returns .
subroutine, public device_col2(a_d, b_d, n)
Vector multiplication .
subroutine, public device_col3(a_d, b_d, c_d, n)
Vector multiplication with 3 vectors .
subroutine, public device_cfill(a_d, c, n)
Set all elements to a constant c .
subroutine, public device_invcol2(a_d, b_d, n)
Vector division .
Device abstraction, common interface for various accelerators.
Definition: device.F90:34
Defines a field.
Definition: field.f90:34
Computes various statistics for the fluid fields. We use the Reynolds decomposition for a field u = ...
Definition: fluid_stats.f90:36
subroutine fluid_stats_reset(this)
Resets all the computed means values and sampling times to zero.
subroutine fluid_stats_make_strong_grad(this)
subroutine fluid_stats_update(this, k)
Updates all fields with a new sample.
subroutine fluid_stats_free(this)
Destructor.
subroutine fluid_stats_post_process(this, mean, reynolds, pressure_flatness, pressure_skewness, skewness_tensor, mean_vel_grad, dissipation_tensor)
Compute certain physical statistical quantities based on existing mean fields.
subroutine fluid_stats_init(this, coef, u, v, w, p, set)
Constructor. Initialize the fields associated with fluid_stats.
Definition: math.f90:60
subroutine, public invers2(a, b, n)
Compute inverted vector .
Definition: math.f90:500
subroutine, public subcol3(a, b, c, n)
Returns .
Definition: math.f90:756
subroutine, public addcol3(a, b, c, n)
Returns .
Definition: math.f90:801
subroutine, public col2(a, b, n)
Vector multiplication .
Definition: math.f90:729
subroutine, public copy(a, b, n)
Copy a vector .
Definition: math.f90:239
subroutine, public col3(a, b, c, n)
Vector multiplication with 3 vectors .
Definition: math.f90:742
Implements mean_field_t.
Definition: mean_field.f90:35
Build configurations.
Definition: neko_config.f90:34
integer, parameter neko_bcknd_device
Definition: neko_config.f90:44
integer, parameter, public rp
Global precision used in computations.
Definition: num_types.f90:12
Operators.
Definition: operators.f90:34
subroutine, public opgrad(ux, uy, uz, u, coef, es, ee)
Compute the weak gradient of a scalar field, i.e. the gradient multiplied by the mass matrix.
Definition: operators.f90:171
Defines a statistical quantity.
Definition: stats_quant.f90:34
Utilities.
Definition: utils.f90:35
subroutine, public neko_warning(warning_msg)
Reports a warning to standard output.
Definition: utils.f90:245
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
Definition: coef.f90:55
field_list_t, To be able to group fields together
Definition: field_list.f90:13
Computes the temporal mean of a field.
Definition: mean_field.f90:45
Abstract type defining a statistical quantity.
Definition: stats_quant.f90:40