Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
les_model_fctry.f90
Go to the documentation of this file.
1! Copyright (c) 2021-2025, 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(les_model) les_model_fctry
34 use vreman, only : vreman_t
35 use smagorinsky, only : smagorinsky_t
37 use sigma, only : sigma_t
39 use wale, only : wale_t
41 implicit none
42
43 ! List of all possible types created by the factory routine
44 character(len=20) :: LES_KNOWN_TYPES(5) = [character(len=20) :: &
45 "vreman", &
46 "smagorinsky", &
47 "dymamic_smagorinsky", &
48 "sigma", &
49 "wale"]
50
51contains
59 module subroutine les_model_factory(object, type_name, fluid, json)
60 class(les_model_t), allocatable, intent(inout) :: object
61 character(len=*), intent(in) :: type_name
62 class(fluid_scheme_base_t), intent(inout) :: fluid
63 type(json_file), intent(inout) :: json
64 character(len=:), allocatable :: type_string
65
66 call les_model_allocator(object, type_name)
67 call object%init(fluid, json)
68 end subroutine les_model_factory
69
73 module subroutine les_model_allocator(object, type_name)
74 class(les_model_t), allocatable, intent(inout) :: object
75 character(len=*), intent(in) :: type_name
76 integer :: i
77
78 if (allocated(object)) deallocate(object)
79
80 select case (trim(type_name))
81 case ('vreman')
82 allocate(vreman_t::object)
83 case ('smagorinsky')
84 allocate(smagorinsky_t::object)
85 case ('dynamic_smagorinsky')
86 allocate(dynamic_smagorinsky_t::object)
87 case ('sigma')
88 allocate(sigma_t::object)
89 case ('wale')
90 allocate(wale_t::object)
91 case default
92 do i = 1, les_model_registry_size
93 if (trim(type_name) == trim(les_model_registry(i)%type_name)) then
94 call les_model_registry(i)%allocator(object)
95 return
96 end if
97 end do
98
99 call neko_type_error("LES model", type_name, les_known_types)
100 end select
101
102 end subroutine les_model_allocator
103
109 module subroutine register_les_model(type_name, allocator)
110 character(len=*), intent(in) :: type_name
111 procedure(les_model_allocate), pointer, intent(in) :: allocator
112 type(allocator_entry), allocatable :: temp(:)
113 integer :: i
114
115 do i = 1, size(les_known_types)
116 if (trim(type_name) .eq. trim(les_known_types(i))) then
117 call neko_type_registration_error("LES model", type_name, .true.)
118 end if
119 end do
120
121 do i = 1, les_model_registry_size
122 if (trim(type_name) .eq. trim(les_model_registry(i)%type_name)) then
123 call neko_type_registration_error("LES model", type_name, .false.)
124 end if
125 end do
126
127 ! Expand registry
128 if (les_model_registry_size == 0) then
129 allocate(les_model_registry(1))
130 else
131 allocate(temp(les_model_registry_size + 1))
132 temp(1:les_model_registry_size) = les_model_registry
133 call move_alloc(temp, les_model_registry)
134 end if
135
136 les_model_registry_size = les_model_registry_size + 1
137 les_model_registry(les_model_registry_size)%type_name = type_name
138 les_model_registry(les_model_registry_size)%allocator => allocator
139 end subroutine register_les_model
140
141end submodule les_model_fctry
Implements dynamic_smagorinsky_t.
Implements les_model_t.
Definition les_model.f90:35
Implements sigma_t.
Definition sigma.f90:35
Implements smagorinsky_t.
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 wale_t.
Definition wale.f90:35
Implements the dynamic Smagorinsky LES model.
Base type of all fluid formulations.
Implements the Sigma LES model.
Definition sigma.f90:52
Implements the smagorinsky LES model.
Implements the Vreman LES model.
Definition vreman.f90:52
Implements the Wale LES model.
Definition wale.f90:52