Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
wall_model_fctry.f90
Go to the documentation of this file.
1! Copyright (c) 2021-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!
33submodule(wall_model) wall_model_fctry
34 use vreman, only : vreman_t
35 use spalding, only : spalding_t
37 use utils, only : neko_type_error
39 implicit none
40
41 ! List of all possible types created by the factory routine
42 character(len=20) :: WALLM_KNOWN_TYPES(2) = [character(len=20) :: &
43 "spalding", &
44 "rough_log_law"]
45
46contains
47
56 module subroutine wall_model_factory(object, scheme_name, coef, msk, facet, &
57 json)
58 class(wall_model_t), allocatable, intent(inout) :: object
59 character(len=*), intent(in) :: scheme_name
60 type(coef_t), intent(in) :: coef
61 integer, intent(in) :: msk(:)
62 integer, intent(in) :: facet(:)
63 type(json_file), intent(inout) :: json
64 character(len=:), allocatable :: type_name
65 character(len=:), allocatable :: type_string
66 integer :: h_index
67
68 call wall_model_allocator(object, type_name)
69
70 call json_get(json, "model", type_name)
71 call json_get(json, "h_index", h_index)
72
73 ! Initialize
74 call object%init(scheme_name, coef, msk, facet, h_index, json)
75
76 end subroutine wall_model_factory
77
81 module subroutine wall_model_allocator(object, type_name)
82 class(wall_model_t), allocatable, intent(inout) :: object
83 character(len=*), intent(in) :: type_name
84 integer :: i
85
86 select case (trim(type_name) )
87 case ("spalding")
88 allocate(spalding_t::object)
89 case ("rough_log_law")
90 allocate(rough_log_law_t::object)
91 case default
92 do i = 1, wall_model_registry_size
93 if (trim(type_name) .eq. trim(wall_model_registry(i)%type_name)) then
94 call wall_model_registry(i)%allocator(object)
95 return
96 end if
97 end do
98 call neko_type_error("wall model", trim(type_name), wallm_known_types)
99 end select
100 end subroutine wall_model_allocator
101
107 module subroutine register_wall_model(type_name, allocator)
108 character(len=*), intent(in) :: type_name
109 procedure(wall_model_allocate), pointer, intent(in) :: allocator
110 type(allocator_entry), allocatable :: temp(:)
111 integer :: i
112
113 do i = 1, size(wallm_known_types)
114 if (trim(type_name) .eq. trim(wallm_known_types(i))) then
115 call neko_type_registration_error("wall model", type_name, .true.)
116 end if
117 end do
118
119 do i = 1, wall_model_registry_size
120 if (trim(type_name) .eq. trim(wall_model_registry(i)%type_name)) then
121 call neko_type_registration_error("wall model", type_name, .false.)
122 end if
123 end do
124
125 ! Expand registry
126 if (wall_model_registry_size .eq. 0) then
127 allocate(wall_model_registry(1))
128 else
129 allocate(temp(wall_model_registry_size + 1))
130 temp(1:wall_model_registry_size) = wall_model_registry
131 call move_alloc(temp, wall_model_registry)
132 end if
133
134 wall_model_registry_size = wall_model_registry_size + 1
135 wall_model_registry(wall_model_registry_size)%type_name = type_name
136 wall_model_registry(wall_model_registry_size)%allocator => allocator
137 end subroutine register_wall_model
138
139end submodule wall_model_fctry
Implements rough_log_law_t.
Implements spalding_t.
Definition spalding.f90:35
Utilities.
Definition utils.f90:35
subroutine, public neko_type_registration_error(base_type, wrong_type, known)
Definition utils.f90:266
subroutine, public neko_type_error(base_type, wrong_type, known_types)
Reports an error allocating a type for a particular base pointer class.
Definition utils.f90:251
Implements vreman_t.
Definition vreman.f90:35
Implements wall_model_t.
Wall model based on the log-law for a rough wall. The formula defining the law is ....
Wall model based on Spalding's law of the wall. Reference: http://dx.doi.org/10.1115/1....
Definition spalding.f90:57
Implements the Vreman LES model.
Definition vreman.f90:52