Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
smagorinsky_cpu.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!
35 use num_types, only : rp
36 use field_list, only : field_list_t
39 use field, only : field_t
40 use operators, only : strain_rate
41 use coefs, only : coef_t
42 use gs_ops, only : gs_op_add
43 use math, only : col2
44 implicit none
45 private
46
48
49contains
50
59 subroutine smagorinsky_compute_cpu(if_ext, t, tstep, coef, nut, delta, c_s)
60 logical, intent(in) :: if_ext
61 real(kind=rp), intent(in) :: t
62 integer, intent(in) :: tstep
63 type(coef_t), intent(in) :: coef
64 type(field_t), intent(inout) :: nut
65 type(field_t), intent(in) :: delta
66 real(kind=rp), intent(in) :: c_s
67 type(field_t), pointer :: u, v, w
68 ! strain rate tensor
69 type(field_t), pointer :: s11, s22, s33, s12, s13, s23
70 real(kind=rp) :: s_abs
71 integer :: temp_indices(6)
72 integer :: e, i
73
74 if (if_ext .eqv. .true.) then
75 u => neko_field_registry%get_field_by_name("u_e")
76 v => neko_field_registry%get_field_by_name("v_e")
77 w => neko_field_registry%get_field_by_name("w_e")
78 else
79 u => neko_field_registry%get_field_by_name("u")
80 v => neko_field_registry%get_field_by_name("v")
81 w => neko_field_registry%get_field_by_name("w")
82 end if
83
84 call neko_scratch_registry%request_field(s11, temp_indices(1))
85 call neko_scratch_registry%request_field(s22, temp_indices(2))
86 call neko_scratch_registry%request_field(s33, temp_indices(3))
87 call neko_scratch_registry%request_field(s12, temp_indices(4))
88 call neko_scratch_registry%request_field(s13, temp_indices(5))
89 call neko_scratch_registry%request_field(s23, temp_indices(6))
90
91 ! Compute the strain rate tensor
92 call strain_rate(s11%x, s22%x, s33%x, s12%x, s13%x, s23%x, u, v, w, coef)
93
94 call coef%gs_h%op(s11, gs_op_add)
95 call coef%gs_h%op(s22, gs_op_add)
96 call coef%gs_h%op(s33, gs_op_add)
97 call coef%gs_h%op(s12, gs_op_add)
98 call coef%gs_h%op(s13, gs_op_add)
99 call coef%gs_h%op(s23, gs_op_add)
100
101 do concurrent(e = 1:coef%msh%nelv)
102 do concurrent(i = 1:coef%Xh%lxyz)
103 s_abs = sqrt(2.0_rp * (s11%x(i,1,1,e)*s11%x(i,1,1,e) + &
104 s22%x(i,1,1,e)*s22%x(i,1,1,e) + &
105 s33%x(i,1,1,e)*s33%x(i,1,1,e)) + &
106 4.0_rp * (s12%x(i,1,1,e)*s12%x(i,1,1,e) + &
107 s13%x(i,1,1,e)*s13%x(i,1,1,e) + &
108 s23%x(i,1,1,e)*s23%x(i,1,1,e)))
109
110 nut%x(i,1,1,e) = c_s**2 * delta%x(i,1,1,e)**2 * s_abs &
111 * coef%mult(i,1,1,e)
112 end do
113 end do
114
115 call coef%gs_h%op(nut, gs_op_add)
116 call col2(nut%x, coef%mult, nut%dof%size())
117
118 call neko_scratch_registry%relinquish_field(temp_indices)
119 end subroutine smagorinsky_compute_cpu
120
121end module smagorinsky_cpu
122
Coefficients.
Definition coef.f90:34
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
Definition math.f90:60
subroutine, public col2(a, b, n)
Vector multiplication .
Definition math.f90:860
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
Operators.
Definition operators.f90:34
subroutine, public strain_rate(s11, s22, s33, s12, s13, s23, u, v, w, coef)
Compute the strain rate tensor, i.e 0.5 * du_i/dx_j + du_j/dx_i.
Defines a registry for storing and requesting temporary fields This can be used when you have a funct...
type(scratch_registry_t), target, public neko_scratch_registry
Global scratch registry.
Implements the CPU kernel for the smagorinsky_t type.
subroutine, public smagorinsky_compute_cpu(if_ext, t, tstep, coef, nut, delta, c_s)
Compute eddy viscosity on the CPU.
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