Neko  0.9.0
A portable framework for high-order spectral element flow simulations
les_model.f90
Go to the documentation of this file.
1 ! Copyright (c) 2023-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 !
33 !
35 module les_model
36  use num_types, only : rp
37  use field, only : field_t, field_ptr_t
38  use json_module, only : json_file
40  use dofmap, only : dofmap_t
41  use coefs, only : coef_t
42  use gs_ops, only : gs_op_add
43  use neko_config, only : neko_bcknd_device
45  use math, only : col2
46  use device_math, only : device_col2
47  implicit none
48  private
49 
51  type, abstract, public :: les_model_t
53  type(field_t), pointer :: nut => null()
55  character(len=:), allocatable :: delta_type
57  type(field_t), pointer :: delta => null()
59  type(coef_t), pointer :: coef => null()
60  contains
62  procedure, pass(this) :: init_base => les_model_init_base
64  procedure, pass(this) :: free_base => les_model_free_base
66  procedure, pass(this) :: compute_delta => les_model_compute_delta
68  procedure(les_model_init), pass(this), deferred :: init
70  procedure(les_model_free), pass(this), deferred :: free
72  procedure(les_model_compute), pass(this), deferred :: compute
73  end type les_model_t
74 
75  abstract interface
76 
79  subroutine les_model_compute(this, t, tstep)
80  import les_model_t, rp
81  class(les_model_t), intent(inout) :: this
82  real(kind=rp), intent(in) :: t
83  integer, intent(in) :: tstep
84  end subroutine les_model_compute
85  end interface
86 
87  abstract interface
88 
92  subroutine les_model_init(this, dofmap, coef, json)
93  import les_model_t, json_file, dofmap_t, coef_t
94  class(les_model_t), intent(inout) :: this
95  type(coef_t), intent(in) :: coef
96  type(dofmap_t), intent(in) :: dofmap
97  type(json_file), intent(inout) :: json
98  end subroutine les_model_init
99  end interface
100 
101  abstract interface
102 
103  subroutine les_model_free(this)
104  import les_model_t
105  class(les_model_t), intent(inout) :: this
106  end subroutine les_model_free
107  end interface
108 
109  interface
110 
116  module subroutine les_model_factory(object, type_name, dofmap, coef, json)
117  class(les_model_t), allocatable, target, intent(inout) :: object
118  character(len=*), intent(in) :: type_name
119  type(dofmap_t), intent(in) :: dofmap
120  type(coef_t), intent(in) :: coef
121  type(json_file), intent(inout) :: json
122  end subroutine les_model_factory
123  end interface
124 
125  public :: les_model_factory
126 
127 contains
132  subroutine les_model_init_base(this, dofmap, coef, nut_name, delta_type)
133  class(les_model_t), intent(inout) :: this
134  type(dofmap_t), intent(in) :: dofmap
135  type(coef_t), target, intent(in) :: coef
136  character(len=*), intent(in) :: nut_name
137  character(len=*), intent(in) :: delta_type
138 
139  if (.not. neko_field_registry%field_exists(trim(nut_name))) then
140  call neko_field_registry%add_field(dofmap, trim(nut_name))
141  end if
142  if (.not. neko_field_registry%field_exists("les_delta")) then
143  call neko_field_registry%add_field(dofmap, "les_delta")
144  end if
145  this%nut => neko_field_registry%get_field(trim(nut_name))
146  this%delta => neko_field_registry%get_field("les_delta")
147  this%coef => coef
148  this%delta_type = delta_type
149 
150  call this%compute_delta()
151  end subroutine les_model_init_base
152 
154  subroutine les_model_free_base(this)
155  class(les_model_t), intent(inout) :: this
156 
157  nullify(this%nut)
158  nullify(this%delta)
159  nullify(this%coef)
160  end subroutine les_model_free_base
161 
168  subroutine les_model_compute_delta(this)
169  class(les_model_t), intent(inout) :: this
170  integer :: e, i, j, k
171  integer :: im, ip, jm, jp, km, kp
172  real(kind=rp) :: di, dj, dk, ndim_inv
173  integer :: lx_half, ly_half, lz_half
174 
175  lx_half = this%coef%Xh%lx / 2
176  ly_half = this%coef%Xh%ly / 2
177  lz_half = this%coef%Xh%lz / 2
178 
179  if (this%delta_type .eq. "elementwise") then
180  ! use a same length scale throughout an entire element
181  ! the length scale is based on maximum GLL spacing
182  do e = 1, this%coef%msh%nelv
183  di = (this%coef%dof%x(lx_half, 1, 1, e) &
184  - this%coef%dof%x(lx_half + 1, 1, 1, e))**2 &
185  + (this%coef%dof%y(lx_half, 1, 1, e) &
186  - this%coef%dof%y(lx_half + 1, 1, 1, e))**2 &
187  + (this%coef%dof%z(lx_half, 1, 1, e) &
188  - this%coef%dof%z(lx_half + 1, 1, 1, e))**2
189 
190  dj = (this%coef%dof%x(1, ly_half, 1, e) &
191  - this%coef%dof%x(1, ly_half + 1, 1, e))**2 &
192  + (this%coef%dof%y(1, ly_half, 1, e) &
193  - this%coef%dof%y(1, ly_half + 1, 1, e))**2 &
194  + (this%coef%dof%z(1, ly_half, 1, e) &
195  - this%coef%dof%z(1, ly_half + 1, 1, e))**2
196 
197  dk = (this%coef%dof%x(1, 1, lz_half, e) &
198  - this%coef%dof%x(1, 1, lz_half + 1, e))**2 &
199  + (this%coef%dof%y(1, 1, lz_half, e) &
200  - this%coef%dof%y(1, 1, lz_half + 1, e))**2 &
201  + (this%coef%dof%z(1, 1, lz_half, e) &
202  - this%coef%dof%z(1, 1, lz_half + 1, e))**2
203  di = sqrt(di)
204  dj = sqrt(dj)
205  dk = sqrt(dk)
206  this%delta%x(:,:,:,e) = (di * dj * dk)**(1.0_rp / 3.0_rp)
207  end do
208 
209  else if (this%delta_type .eq. "pointwise") then
210  do e = 1, this%coef%msh%nelv
211  do k = 1, this%coef%Xh%lz
212  km = max(1, k-1)
213  kp = min(this%coef%Xh%lz, k+1)
214 
215  do j = 1, this%coef%Xh%ly
216  jm = max(1, j-1)
217  jp = min(this%coef%Xh%ly, j+1)
218 
219  do i = 1, this%coef%Xh%lx
220  im = max(1, i-1)
221  ip = min(this%coef%Xh%lx, i+1)
222 
223  di = (this%coef%dof%x(ip, j, k, e) - &
224  this%coef%dof%x(im, j, k, e))**2 &
225  + (this%coef%dof%y(ip, j, k, e) - &
226  this%coef%dof%y(im, j, k, e))**2 &
227  + (this%coef%dof%z(ip, j, k, e) - &
228  this%coef%dof%z(im, j, k, e))**2
229 
230  dj = (this%coef%dof%x(i, jp, k, e) - &
231  this%coef%dof%x(i, jm, k, e))**2 &
232  + (this%coef%dof%y(i, jp, k, e) - &
233  this%coef%dof%y(i, jm, k, e))**2 &
234  + (this%coef%dof%z(i, jp, k, e) - &
235  this%coef%dof%z(i, jm, k, e))**2
236 
237  dk = (this%coef%dof%x(i, j, kp, e) - &
238  this%coef%dof%x(i, j, km, e))**2 &
239  + (this%coef%dof%y(i, j, kp, e) - &
240  this%coef%dof%y(i, j, km, e))**2 &
241  + (this%coef%dof%z(i, j, kp, e) - &
242  this%coef%dof%z(i, j, km, e))**2
243 
244  di = sqrt(di) / (ip - im)
245  dj = sqrt(dj) / (jp - jm)
246  dk = sqrt(dk) / (kp - km)
247  this%delta%x(i,j,k,e) = (di * dj * dk)**(1.0_rp / 3.0_rp)
248 
249  end do
250  end do
251  end do
252  end do
253  end if
254 
255  if (neko_bcknd_device .eq. 1) then
256  call device_memcpy(this%delta%x, this%delta%x_d, this%delta%dof%size(),&
257  host_to_device, sync = .false.)
258  call this%coef%gs_h%op(this%delta%x, this%delta%dof%size(), gs_op_add)
259  call device_col2(this%delta%x_d, this%coef%mult_d, this%delta%dof%size())
260  else
261  call this%coef%gs_h%op(this%delta%x, this%delta%dof%size(), gs_op_add)
262  call col2(this%delta%x, this%coef%mult, this%delta%dof%size())
263  end if
264 
265  end subroutine les_model_compute_delta
266 
267 end module les_model
Copy data between host and device (or device and device)
Definition: device.F90:51
Compute eddy viscosity.
Definition: les_model.f90:79
Common constructor.
Definition: les_model.f90:92
Coefficients.
Definition: coef.f90:34
subroutine, public device_col2(a_d, b_d, n)
Vector multiplication .
Device abstraction, common interface for various accelerators.
Definition: device.F90:34
integer, parameter, public host_to_device
Definition: device.F90:47
Defines a mapping of the degrees of freedom.
Definition: dofmap.f90:35
Defines a registry for storing solution fields.
type(field_registry_t), target, public neko_field_registry
Global field registry.
Defines a field.
Definition: field.f90:34
Defines Gather-scatter operations.
Definition: gs_ops.f90:34
integer, parameter, public gs_op_add
Definition: gs_ops.f90:36
Implements les_model_t.
Definition: les_model.f90:35
subroutine les_model_init_base(this, dofmap, coef, nut_name, delta_type)
LES model factory. Both constructs and initializes the object.
Definition: les_model.f90:133
subroutine les_model_compute_delta(this)
Compute the LES lengthscale. For each GLL point, we take the distance between its neighbours in all 3...
Definition: les_model.f90:169
subroutine les_model_free_base(this)
Destructor for the les_model_t (base) class.
Definition: les_model.f90:155
Definition: math.f90:60
subroutine, public col2(a, b, n)
Vector multiplication .
Definition: math.f90:729
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
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
Definition: coef.f90:55
field_ptr_t, To easily obtain a pointer to a field
Definition: field.f90:80
Base abstract type for LES models based on the Boussinesq approximation.
Definition: les_model.f90:51
#define max(a, b)
Definition: tensor.cu:40