Neko  0.9.99
A portable framework for high-order spectral element flow simulations
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 !
33 !
36  use num_types, only : rp
37  use field, only : field_t
38  use les_model, only : les_model_t
39  use dofmap , only : dofmap_t
41  use json_module, only : json_file
42  use utils, only : neko_error
43  use neko_config, only : neko_bcknd_device
45  use coefs, only : coef_t
46  implicit none
47  private
48 
51  type, public, extends(les_model_t) :: smagorinsky_t
53  real(kind=rp) :: c_s
54  contains
56  procedure, pass(this) :: init => smagorinsky_init
58  procedure, pass(this) :: init_from_components => &
61  procedure, pass(this) :: free => smagorinsky_free
63  procedure, pass(this) :: compute => smagorinsky_compute
64  end type smagorinsky_t
65 
66 contains
71  subroutine smagorinsky_init(this, dofmap, coef, json)
72  class(smagorinsky_t), intent(inout) :: this
73  type(dofmap_t), intent(in) :: dofmap
74  type(coef_t), intent(in) :: coef
75  type(json_file), intent(inout) :: json
76  character(len=:), allocatable :: nut_name
77  real(kind=rp) :: c_s
78  character(len=:), allocatable :: delta_type
79 
80  call json_get_or_default(json, "nut_field", nut_name, "nut")
81  call json_get_or_default(json, "delta_type", delta_type, "pointwise")
82  call json_get_or_default(json, "c_s", c_s, 0.17_rp)
83 
84  call smagorinsky_init_from_components(this, dofmap, coef, c_s, nut_name, &
85  delta_type)
86  end subroutine smagorinsky_init
87 
93  subroutine smagorinsky_init_from_components(this, dofmap, coef, c_s, &
94  nut_name, delta_type)
95  class(smagorinsky_t), intent(inout) :: this
96  type(dofmap_t), intent(in) :: dofmap
97  type(coef_t), intent(in) :: coef
98  real(kind=rp) :: c_s
99  character(len=*), intent(in) :: nut_name
100  character(len=*), intent(in) :: delta_type
101 
102  call this%free()
103 
104  call this%init_base(dofmap, coef, nut_name, delta_type)
105  this%c_s = c_s
106 
107  end subroutine smagorinsky_init_from_components
108 
110  subroutine smagorinsky_free(this)
111  class(smagorinsky_t), intent(inout) :: this
112 
113  call this%free_base()
114  end subroutine smagorinsky_free
115 
119  subroutine smagorinsky_compute(this, t, tstep)
120  class(smagorinsky_t), intent(inout) :: this
121  real(kind=rp), intent(in) :: t
122  integer, intent(in) :: tstep
123 
124  if (neko_bcknd_device .eq. 1) then
125  call neko_error("Smagorinsky model not implemented on accelarators.")
126  else
127  call smagorinsky_compute_cpu(t, tstep, this%coef, this%nut, this%delta,&
128  this%c_s)
129  end if
130 
131  end subroutine smagorinsky_compute
132 
133 end module 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
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
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
Implements the CPU kernel for the smagorinsky_t type.
subroutine, public smagorinsky_compute_cpu(t, tstep, coef, nut, delta, c_s)
Compute eddy viscosity on the CPU.
Implements smagorinsky_t.
Definition: smagorinsky.f90:35
subroutine smagorinsky_init_from_components(this, dofmap, coef, c_s, nut_name, delta_type)
Constructor from components.
Definition: smagorinsky.f90:95
subroutine smagorinsky_free(this)
Destructor for the les_model_t (base) class.
subroutine smagorinsky_compute(this, t, tstep)
Compute eddy viscosity.
subroutine smagorinsky_init(this, dofmap, coef, json)
Constructor.
Definition: smagorinsky.f90:72
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
Base abstract type for LES models based on the Boussinesq approximation.
Definition: les_model.f90:51
Implements the smagorinsky LES model.
Definition: smagorinsky.f90:51