Neko  0.9.0
A portable framework for high-order spectral element flow simulations
dynamic_smagorinsky.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 math
37  use field_list, only : field_list_t
38  use field, only : field_t
39  use les_model, only : les_model_t
40  use dofmap , only : dofmap_t
42  use json_module, only : json_file
43  use utils, only : neko_error
44  use neko_config, only : neko_bcknd_device
45  use coefs, only : coef_t
48  implicit none
49  private
50 
53  type, public, extends(les_model_t) :: dynamic_smagorinsky_t
55  type(field_t) :: c_dyn
57  character(len=:), allocatable :: test_filter_type
58  type(elementwise_filter_t) :: test_filter
62  type(field_t) :: mij(6)
64  type(field_t) :: lij(6)
66  type(field_t) :: num
68  type(field_t) :: den
69  contains
71  procedure, pass(this) :: init => dynamic_smagorinsky_init
73  procedure, pass(this) :: free => dynamic_smagorinsky_free
75  procedure, pass(this) :: compute => dynamic_smagorinsky_compute
76 
77  end type dynamic_smagorinsky_t
78 
79 contains
84  subroutine dynamic_smagorinsky_init(this, dofmap, coef, json)
85  class(dynamic_smagorinsky_t), intent(inout) :: this
86  type(dofmap_t), intent(in) :: dofmap
87  type(coef_t), intent(in) :: coef
88  type(json_file), intent(inout) :: json
89  character(len=:), allocatable :: nut_name !! The name of the SGS viscosity field.
90  integer :: i
91  character(len=:), allocatable :: delta_type
92 
93  call json_get_or_default(json, "nut_field", nut_name, "nut")
94  call json_get_or_default(json, "delta_type", delta_type, "pointwise")
95 
96  call this%free()
97  call this%init_base(dofmap, coef, nut_name, delta_type)
98  this%test_filter_type = "nonBoyd"
99  ! Filter assumes lx = ly = lz
100  call this%test_filter%init(dofmap%xh%lx, this%test_filter_type)
101  call set_ds_filt(this%test_filter)
102 
103  call this%c_dyn%init(dofmap, "ds_c_dyn")
104  call this%num%init(dofmap, "ds_num")
105  call this%den%init(dofmap, "ds_den")
106 
107  do i = 1, 6
108  call this%mij(i)%init(dofmap)
109  call this%lij(i)%init(dofmap)
110  end do
111 
112  end subroutine dynamic_smagorinsky_init
113 
115  subroutine dynamic_smagorinsky_free(this)
116  class(dynamic_smagorinsky_t), intent(inout) :: this
117  integer :: i
118 
119  call this%c_dyn%free()
120  do i = 1, 6
121  call this%mij(i)%free()
122  call this%lij(i)%free()
123  end do
124  call this%num%free()
125  call this%den%free()
126  call this%test_filter%free()
127  call this%free_base()
128 
129  end subroutine dynamic_smagorinsky_free
130 
134  subroutine dynamic_smagorinsky_compute(this, t, tstep)
135  class(dynamic_smagorinsky_t), intent(inout) :: this
136  real(kind=rp), intent(in) :: t
137  integer, intent(in) :: tstep
138 
139  if (neko_bcknd_device .eq. 1) then
140  call neko_error("Dynamic Smagorinsky model not implemented on &
141  &accelarators.")
142  else
143  call dynamic_smagorinsky_compute_cpu(t, tstep, this%coef, this%nut, &
144  this%delta, this%c_dyn, this%test_filter, &
145  this%mij, this%lij, this%num, this%den)
146  end if
147 
148  end subroutine dynamic_smagorinsky_compute
149 
151  subroutine set_ds_filt(filter_1d)
152  type(elementwise_filter_t), intent(inout) :: filter_1d
153  integer :: i
154 
155  if (filter_1d%nx .le. 2) then
156  call neko_error("Dynamic Smagorinsky model error: test filter is not &
157  &defined for the current polynomial order")
158  end if
159  if (mod(filter_1d%nx,2) .eq. 0) then ! number of grid spacing is odd
160  ! cutoff at polynomial order int((filter_1d%nx)/2)
161  filter_1d%trnsfr(int((filter_1d%nx)/2)) = 0.95_rp
162  filter_1d%trnsfr(int((filter_1d%nx)/2) + 1) = 0.50_rp
163  filter_1d%trnsfr(int((filter_1d%nx)/2) + 2) = 0.05_rp
164  if ((int((filter_1d%nx)/2) + 2) .lt. filter_1d%nx) then
165  do i = int((filter_1d%nx)/2) + 3, filter_1d%nx
166  filter_1d%trnsfr(i) = 0.0_rp
167  end do
168  end if
169  ! make delta_ratio = (nx-1)/(nt-1) as close to 2
170  filter_1d%nt = int(filter_1d%nx/2) + 1
171  else ! number of grid spacing is even
172  ! cutoff at polynomial order int((filter_1d%nx-1)/2)
173  filter_1d%trnsfr(int((filter_1d%nx-1)/2)) = 0.95_rp
174  filter_1d%trnsfr(int((filter_1d%nx-1)/2) + 1) = 0.50_rp
175  filter_1d%trnsfr(int((filter_1d%nx-1)/2) + 2) = 0.05_rp
176  if ((int((filter_1d%nx-1)/2) + 2) .lt. filter_1d%nx) then
177  do i = int((filter_1d%nx-1)/2) + 3, filter_1d%nx
178  filter_1d%trnsfr(i) = 0.0_rp
179  end do
180  end if
181  ! make delta_ratio = (nx-1)/(nt-1) = 2
182  filter_1d%nt = int((filter_1d%nx-1)/2) + 1
183  end if
184 
185  call filter_1d%build_1d()
186 
187  end subroutine set_ds_filt
188 
189 end module dynamic_smagorinsky
Retrieves a parameter by name or assigns a provided default value. In the latter case also adds the m...
Definition: json_utils.f90:54
Retrieves a parameter by name or throws an error.
Definition: json_utils.f90:45
Coefficients.
Definition: coef.f90:34
Defines a mapping of the degrees of freedom.
Definition: dofmap.f90:35
Implements the CPU kernel for the smagorinsky_t type.
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 dynamic_smagorinsky_t.
subroutine set_ds_filt(filter_1d)
Set up the test filter.
subroutine dynamic_smagorinsky_init(this, dofmap, coef, json)
Constructor.
subroutine dynamic_smagorinsky_free(this)
Destructor for the les_model_t (base) class.
subroutine dynamic_smagorinsky_compute(this, t, tstep)
Compute eddy viscosity.
Implements explicit_filter_t.
Defines a field.
Definition: field.f90:34
Utilities for retrieving parameters from the case files.
Definition: json_utils.f90:34
Implements les_model_t.
Definition: les_model.f90:35
Definition: math.f90:60
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
Utilities.
Definition: utils.f90:35
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
Definition: coef.f90:55
Implements the dynamic Smagorinsky LES model.
Implements the explicit filter for SEM.
field_list_t, To be able to group fields together
Definition: field_list.f90:13
Base abstract type for LES models based on the Boussinesq approximation.
Definition: les_model.f90:51