Neko-TOP
A portable framework for high-order spectral element flow toplogy optimization.
Loading...
Searching...
No Matches
adjoint_fluid_pnpn_bc_fctry.f90
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!
35submodule(adjoint_fluid_pnpn) adjoint_fluid_pnpn_bc_fctry
36 use user_intf, only: user_t
37 use utils, only: neko_type_error
38 use field_dirichlet, only: field_dirichlet_t
39 use inflow, only: inflow_t
40 use blasius, only: blasius_t
41 use dirichlet, only: dirichlet_t
42 use dong_outflow, only: dong_outflow_t
43 use symmetry, only: symmetry_t
44 use non_normal, only: non_normal_t
45 use field_dirichlet_vector, only: field_dirichlet_vector_t
46 implicit none
47
48 ! List of all possible types created by the boundary condition factories
49 character(len=25) :: ADJOINT_FLUID_PNPN_KNOWN_BCS(13) = &
50 [character(len=25) :: &
51 "symmetry", &
52 "velocity_value", &
53 "no_slip", &
54 "outflow", &
55 "normal_outflow", &
56 "outflow+dong", &
57 "normal_outflow+dong", &
58 "shear_stress", &
59 "user_velocity", &
60 "user_pressure", &
61 "blasius_profile", &
62 "user_velocity_pointwise", &
63 "wall_model"]
64
65contains
66
73 module subroutine pressure_bc_factory(object, scheme, json, coef, user)
74 class(bc_t), pointer, intent(inout) :: object
75 type(adjoint_fluid_pnpn_t), intent(in) :: scheme
76 type(json_file), intent(inout) :: json
77 type(coef_t), intent(in) :: coef
78 type(user_t), intent(in) :: user
79 character(len=:), allocatable :: type
80 integer :: i, j, k
81 integer, allocatable :: zone_indices(:)
82
83 call json_get(json, "type", type)
84
85 select case (trim(type))
86 case ("outflow", "normal_outflow")
87 allocate(zero_dirichlet_t::object)
88
89 case ("outflow+dong", "normal_outflow+dong")
90 allocate(dong_outflow_t::object)
91
92 case ("user_pressure")
93 allocate(field_dirichlet_t::object)
94 select type (obj => object)
95 type is (field_dirichlet_t)
96 obj%update => user%dirichlet_conditions
97 call json%add("field_name", scheme%p_adj%name)
98 end select
99
100 case default
101 do i = 1, size(adjoint_fluid_pnpn_known_bcs)
102 if (trim(type) .eq. trim(adjoint_fluid_pnpn_known_bcs(i))) return
103 end do
104 call neko_type_error("adjoint_fluid_pnpn boundary conditions", type, &
105 ADJOINT_FLUID_PNPN_KNOWN_BCS)
106 end select
107
108 call json_get(json, "zone_indices", zone_indices)
109 call object%init(coef, json)
110
111 do i = 1, size(zone_indices)
112 call object%mark_zone(coef%msh%labeled_zones(zone_indices(i)))
113 end do
114 call object%finalize()
115
116 ! All pressure bcs are currently strong, so for all of them we
117 ! mark with value 1 in the mesh
118 do i = 1, size(zone_indices)
119 do j = 1, scheme%msh%nelv
120 do k = 1, 2 * scheme%msh%gdim
121 if (scheme%msh%facet_type(k,j) .eq. -zone_indices(i)) then
122 scheme%msh%facet_type(k, j) = 1
123 end if
124 end do
125 end do
126 end do
127 end subroutine pressure_bc_factory
128
135 module subroutine velocity_bc_factory(object, scheme, json, coef, user)
136 class(bc_t), pointer, intent(inout) :: object
137 type(adjoint_fluid_pnpn_t), intent(in) :: scheme
138 type(json_file), intent(inout) :: json
139 type(coef_t), intent(in) :: coef
140 type(user_t), intent(in) :: user
141 character(len=:), allocatable :: type
142 integer :: i, j, k
143 integer, allocatable :: zone_indices(:)
144
145 call json_get(json, "type", type)
146
147 select case (trim(type))
148 case ("symmetry")
149 allocate(symmetry_t::object)
150 case ("velocity_value")
151 allocate(inflow_t::object)
152 case ("no_slip")
153 allocate(zero_dirichlet_t::object)
154 case ("normal_outflow", "normal_outflow+dong")
155 allocate(non_normal_t::object)
156 case ("blasius_profile")
157 allocate(blasius_t::object)
158 case ("shear_stress")
159 allocate(shear_stress_t::object)
160 case ("wall_model")
161 allocate(wall_model_bc_t::object)
162 ! Kind of hack, but maybe OK? The thing is, we need the nu for
163 ! initing the wall model, and forcing the user duplicate that there
164 ! would be a nightmare.
165 ! call json%add("nu", scheme%mu / scheme%rho)
166
167 case ("user_velocity")
168 allocate(field_dirichlet_vector_t::object)
169 select type (obj => object)
170 type is (field_dirichlet_vector_t)
171 obj%update => user%dirichlet_conditions
172 end select
173
174 ! case ("user_velocity_pointwise")
175 ! allocate(usr_inflow_t::object)
176 ! select type (obj => object)
177 ! type is (usr_inflow_t)
178 ! call obj%set_eval(user%adjoint_user_if)
179 ! call obj%validate()
180 ! end select
181
182 case default
183 do i = 1, size(adjoint_fluid_pnpn_known_bcs)
184 if (trim(type) .eq. trim(adjoint_fluid_pnpn_known_bcs(i))) return
185 end do
186 call neko_type_error("adjoint_fluid_pnpn boundary conditions", type, &
187 ADJOINT_FLUID_PNPN_KNOWN_BCS)
188 end select
189
190 call json_get(json, "zone_indices", zone_indices)
191 call object%init(coef, json)
192 do i = 1, size(zone_indices)
193 call object%mark_zone(coef%msh%labeled_zones(zone_indices(i)))
194 end do
195 call object%finalize()
196
197 ! Exclude these two because they are bcs for the residual, not velocity
198 if (trim(type) .ne. "normal_outflow" .and. &
199 trim(type) .ne. "normal_outflow+dong") then
200 do i = 1, size(zone_indices)
201 do j = 1, scheme%msh%nelv
202 do k = 1, 2 * scheme%msh%gdim
203 if (scheme%msh%facet_type(k,j) .eq. -zone_indices(i)) then
204 scheme%msh%facet_type(k, j) = 2
205 end if
206 end do
207 end do
208 end do
209 end if
210 end subroutine velocity_bc_factory
211
212end submodule adjoint_fluid_pnpn_bc_fctry
Adjoint Pn/Pn formulation.