Loading [MathJax]/extensions/tex2jax.js
Neko 0.9.99
A portable framework for high-order spectral element flow simulations
All Classes Namespaces Files Functions Variables Typedefs Enumerator Macros Pages
vreman.f90
Go to the documentation of this file.
1! Copyright (c) 2023, 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!
35module vreman
36 use num_types, only : rp
37 use les_model, only : les_model_t
38 use dofmap , only : dofmap_t
40 use json_module, only : json_file
44 use coefs, only : coef_t
45 use logger, only : log_size, neko_log
46 implicit none
47 private
48
51 type, public, extends(les_model_t) :: vreman_t
53 real(kind=rp) :: c
54 contains
56 procedure, pass(this) :: init => vreman_init
58 procedure, pass(this) :: init_from_components => &
61 procedure, pass(this) :: free => vreman_free
63 procedure, pass(this) :: compute => vreman_compute
64 end type vreman_t
65
66contains
71 subroutine vreman_init(this, dofmap, coef, json)
72 class(vreman_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
78 character(len=:), allocatable :: delta_type
79 character(len=LOG_SIZE) :: log_buf
80
81 call json_get_or_default(json, "nut_field", nut_name, "nut")
82 call json_get_or_default(json, "delta_type", delta_type, "pointwise")
83 ! Based on the Smagorinsky Cs = 0.17.
84 call json_get_or_default(json, "c", c, 0.07_rp)
85
86 call neko_log%section('LES model')
87 write(log_buf, '(A)') 'Model : Vreman'
88 call neko_log%message(log_buf)
89 write(log_buf, '(A, A)') 'Delta evaluation : ', delta_type
90 call neko_log%message(log_buf)
91 write(log_buf, '(A, E15.7)') 'c : ', c
92 call neko_log%message(log_buf)
93 call neko_log%end_section()
94
95 call vreman_init_from_components(this, dofmap, coef, c, nut_name, &
96 delta_type)
97 end subroutine vreman_init
98
104 subroutine vreman_init_from_components(this, dofmap, coef, c, nut_name, &
105 delta_type)
106 class(vreman_t), intent(inout) :: this
107 type(dofmap_t), intent(in) :: dofmap
108 type(coef_t), intent(in) :: coef
109 real(kind=rp) :: c
110 character(len=*), intent(in) :: nut_name
111 character(len=*), intent(in) :: delta_type
112
113 call this%free()
114
115 call this%init_base(dofmap, coef, nut_name, delta_type)
116 this%c = c
117
118 end subroutine vreman_init_from_components
119
121 subroutine vreman_free(this)
122 class(vreman_t), intent(inout) :: this
123
124 call this%free_base()
125 end subroutine vreman_free
126
130 subroutine vreman_compute(this, t, tstep)
131 class(vreman_t), intent(inout) :: this
132 real(kind=rp), intent(in) :: t
133 integer, intent(in) :: tstep
134
135 if (neko_bcknd_device .eq. 1) then
136 call vreman_compute_device(t, tstep, this%coef, this%nut, this%delta,&
137 this%c)
138 else
139 call vreman_compute_cpu(t, tstep, this%coef, this%nut, this%delta,&
140 this%c)
141 end if
142
143 end subroutine vreman_compute
144
145end module vreman
Retrieves a parameter by name or assigns a provided default value. In the latter case also adds the m...
Coefficients.
Definition coef.f90:34
Defines a mapping of the degrees of freedom.
Definition dofmap.f90:35
Utilities for retrieving parameters from the case files.
Implements les_model_t.
Definition les_model.f90:35
Logging routines.
Definition log.f90:34
type(log_t), public neko_log
Global log stream.
Definition log.f90:65
integer, parameter, public log_size
Definition log.f90:42
Build configurations.
integer, parameter neko_bcknd_device
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
Implements the CPU kernel for the vreman_t type.
subroutine, public vreman_compute_cpu(t, tstep, coef, nut, delta, c)
Compute eddy viscosity on the CPU.
Implements the device kernel for the vreman_t type.
subroutine, public vreman_compute_device(t, tstep, coef, nut, delta, c)
Compute eddy viscosity on the device.
Implements vreman_t.
Definition vreman.f90:35
subroutine vreman_init(this, dofmap, coef, json)
Constructor.
Definition vreman.f90:72
subroutine vreman_free(this)
Destructor for the les_model_t (base) class.
Definition vreman.f90:122
subroutine vreman_compute(this, t, tstep)
Compute eddy viscosity.
Definition vreman.f90:131
subroutine vreman_init_from_components(this, dofmap, coef, c, nut_name, delta_type)
Constructor from components.
Definition vreman.f90:106
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 Vreman LES model.
Definition vreman.f90:51