Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
symmetry.f90
Go to the documentation of this file.
1! Copyright (c) 2020-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!
36 use dirichlet, only : dirichlet_t
37 use num_types, only : rp
38 use bc, only : bc_t
39 use tuple, only : tuple_i4_t
40 use coefs, only : coef_t
41 use json_module, only : json_file
43 use, intrinsic :: iso_c_binding, only : c_ptr
44 use time_state, only : time_state_t
45 implicit none
46 private
47
50 type, public, extends(bc_t) :: symmetry_t
51 type(zero_dirichlet_t) :: bc_x
52 type(zero_dirichlet_t) :: bc_y
53 type(zero_dirichlet_t) :: bc_z
54 contains
55 procedure, pass(this) :: apply_scalar => symmetry_apply_scalar
56 procedure, pass(this) :: apply_vector => symmetry_apply_vector
57 procedure, pass(this) :: apply_scalar_dev => symmetry_apply_scalar_dev
58 procedure, pass(this) :: apply_vector_dev => symmetry_apply_vector_dev
60 procedure, pass(this) :: init => symmetry_init
62 procedure, pass(this) :: init_from_components => &
65 procedure, pass(this) :: free => symmetry_free
67 procedure, pass(this) :: get_normal_axis => symmetry_get_normal_axis
69 procedure, pass(this) :: finalize => symmetry_finalize
70 end type symmetry_t
71
72contains
76 subroutine symmetry_init(this, coef, json)
77 class(symmetry_t), intent(inout), target :: this
78 type(coef_t), target, intent(in) :: coef
79 type(json_file), intent(inout) ::json
80
81 call this%init_from_components(coef)
82
83 end subroutine symmetry_init
84
87 subroutine symmetry_init_from_components(this, coef)
88 class(symmetry_t), intent(inout), target :: this
89 type(coef_t), target, intent(in) :: coef
90
91 call this%free()
92 this%strong = .false.
93
94 call this%init_base(coef)
95 call this%bc_x%init_from_components(this%coef)
96 call this%bc_y%init_from_components(this%coef)
97 call this%bc_z%init_from_components(this%coef)
98
100
104 subroutine symmetry_finalize(this, only_facets)
105 class(symmetry_t), target, intent(inout) :: this
106 logical, optional, intent(in) :: only_facets
107 logical :: only_facets_ = .false.
108 integer :: i, m, j, l
109 type(tuple_i4_t), pointer :: bfp(:)
110 real(kind=rp) :: sx, sy, sz
111 real(kind=rp), parameter :: tol = 1e-3_rp
112 type(tuple_i4_t) :: bc_facet
113 integer :: facet, el
114
115 if (present(only_facets)) then
116 only_facets_ = only_facets
117 else
118 only_facets_ = .false.
119 end if
120
121 call this%finalize_base(only_facets_)
122
123 associate(c => this%coef, nx => this%coef%nx, ny => this%coef%ny, &
124 nz => this%coef%nz)
125 bfp => this%marked_facet%array()
126 do i = 1, this%marked_facet%size()
127 bc_facet = bfp(i)
128 facet = bc_facet%x(1)
129 el = bc_facet%x(2)
130 call this%get_normal_axis(sx, sy, sz, facet, el)
131
132 if (sx .lt. tol) then
133 call this%bc_x%mark_facet(facet, el)
134 end if
135
136 if (sy .lt. tol) then
137 call this%bc_y%mark_facet(facet, el)
138 end if
139
140 if (sz .lt. tol) then
141 call this%bc_z%mark_facet(facet, el)
142 end if
143 end do
144 end associate
145 call this%bc_x%finalize(only_facets_)
146 call this%bc_y%finalize(only_facets_)
147 call this%bc_z%finalize(only_facets_)
148 end subroutine symmetry_finalize
149
153 subroutine symmetry_get_normal_axis(this, sx, sy, sz, facet, el)
154 class(symmetry_t), target, intent(inout) :: this
155 real(kind=rp), intent(out) :: sx, sy, sz
156 integer, intent(in) :: facet
157 integer, intent(in) :: el
158 integer :: i, m, j, l
159 type(tuple_i4_t), pointer :: bfp(:)
160 ! TODO is the tolerance too large?
161 real(kind=rp), parameter :: tol = 1d-3
162 type(tuple_i4_t) :: bc_facet
163
164 associate(c => this%coef, nx => this%coef%nx, ny => this%coef%ny, &
165 nz => this%coef%nz)
166 sx = 0.0_rp
167 sy = 0d0
168 sz = 0d0
169 select case (facet)
170 case (1, 2)
171 do l = 2, c%Xh%lx - 1
172 do j = 2, c%Xh%lx -1
173 sx = sx + abs(abs(nx(l, j, facet, el)) - 1d0)
174 sy = sy + abs(abs(ny(l, j, facet, el)) - 1d0)
175 sz = sz + abs(abs(nz(l, j, facet, el)) - 1d0)
176 end do
177 end do
178 case (3, 4)
179 do l = 2, c%Xh%lx - 1
180 do j = 2, c%Xh%lx - 1
181 sx = sx + abs(abs(nx(l, j, facet, el)) - 1d0)
182 sy = sy + abs(abs(ny(l, j, facet, el)) - 1d0)
183 sz = sz + abs(abs(nz(l, j, facet, el)) - 1d0)
184 end do
185 end do
186 case (5, 6)
187 do l = 2, c%Xh%lx - 1
188 do j = 2, c%Xh%lx - 1
189 sx = sx + abs(abs(nx(l, j, facet, el)) - 1d0)
190 sy = sy + abs(abs(ny(l, j, facet, el)) - 1d0)
191 sz = sz + abs(abs(nz(l, j, facet, el)) - 1d0)
192 end do
193 end do
194 end select
195 sx = sx / (c%Xh%lx - 2)**2
196 sy = sy / (c%Xh%lx - 2)**2
197 sz = sz / (c%Xh%lx - 2)**2
198 end associate
199 end subroutine symmetry_get_normal_axis
200
206 subroutine symmetry_apply_scalar(this, x, n, time, strong)
207 class(symmetry_t), intent(inout) :: this
208 integer, intent(in) :: n
209 real(kind=rp), intent(inout), dimension(n) :: x
210 type(time_state_t), intent(in), optional :: time
211 logical, intent(in), optional :: strong
212 end subroutine symmetry_apply_scalar
213
221 subroutine symmetry_apply_vector(this, x, y, z, n, time, strong)
222 class(symmetry_t), intent(inout) :: this
223 integer, intent(in) :: n
224 real(kind=rp), intent(inout), dimension(n) :: x
225 real(kind=rp), intent(inout), dimension(n) :: y
226 real(kind=rp), intent(inout), dimension(n) :: z
227 type(time_state_t), intent(in), optional :: time
228 logical, intent(in), optional :: strong
229 logical :: strong_
230
231 if (present(strong)) then
232 strong_ = strong
233 else
234 strong_ = .true.
235 end if
236
237 if (strong_) then
238 call this%bc_x%apply_scalar(x, n)
239 call this%bc_y%apply_scalar(y, n)
240 call this%bc_z%apply_scalar(z, n)
241 end if
242
243 end subroutine symmetry_apply_vector
244
249 subroutine symmetry_apply_scalar_dev(this, x_d, time, strong, strm)
250 class(symmetry_t), intent(inout), target :: this
251 type(c_ptr), intent(inout) :: x_d
252 type(time_state_t), intent(in), optional :: time
253 logical, intent(in), optional :: strong
254 type(c_ptr), intent(inout) :: strm
255 end subroutine symmetry_apply_scalar_dev
256
263 subroutine symmetry_apply_vector_dev(this, x_d, y_d, z_d, &
264 time, strong, strm)
265 class(symmetry_t), intent(inout), target :: this
266 type(c_ptr), intent(inout) :: x_d
267 type(c_ptr), intent(inout) :: y_d
268 type(c_ptr), intent(inout) :: z_d
269 type(time_state_t), intent(in), optional :: time
270 logical, intent(in), optional :: strong
271 type(c_ptr), intent(inout) :: strm
272 logical :: strong_
273
274 if (present(strong)) then
275 strong_ = strong
276 else
277 strong_ = .true.
278 end if
279
280 if (strong_ .and. (this%msk(0) .gt. 0)) then
281 call device_symmetry_apply_vector(this%bc_x%msk_d, this%bc_y%msk_d, &
282 this%bc_z%msk_d, x_d, y_d, z_d, &
283 this%bc_x%msk(0), this%bc_y%msk(0), this%bc_z%msk(0), strm)
284 end if
285 end subroutine symmetry_apply_vector_dev
286
288 subroutine symmetry_free(this)
289 class(symmetry_t), target, intent(inout) :: this
290
291 call this%free_base()
292 call this%bc_x%free()
293 call this%bc_y%free()
294 call this%bc_z%free()
295
296 end subroutine symmetry_free
297end module symmetry
Defines a boundary condition.
Definition bc.f90:34
Coefficients.
Definition coef.f90:34
subroutine device_symmetry_apply_vector(xmsk, ymsk, zmsk, x, y, z, m, n, l, strm)
Defines a dirichlet boundary condition.
Definition dirichlet.f90:34
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
Mixed Dirichlet-Neumann axis aligned symmetry plane.
Definition symmetry.f90:34
subroutine symmetry_apply_vector(this, x, y, z, n, time, strong)
Apply symmetry conditions (axis aligned)
Definition symmetry.f90:222
subroutine symmetry_apply_scalar(this, x, n, time, strong)
No-op scalar apply.
Definition symmetry.f90:207
subroutine symmetry_init_from_components(this, coef)
Constructor from components.
Definition symmetry.f90:88
subroutine symmetry_init(this, coef, json)
Constructor.
Definition symmetry.f90:77
subroutine symmetry_apply_scalar_dev(this, x_d, time, strong, strm)
No-op scalar apply (device version)
Definition symmetry.f90:250
subroutine symmetry_apply_vector_dev(this, x_d, y_d, z_d, time, strong, strm)
Apply symmetry conditions (axis aligned) (device version)
Definition symmetry.f90:265
subroutine symmetry_free(this)
Destructor.
Definition symmetry.f90:289
subroutine symmetry_finalize(this, only_facets)
Finalize. Marks the appropriate faces for application of a homogeneous Dirchlet condition based on th...
Definition symmetry.f90:105
subroutine symmetry_get_normal_axis(this, sx, sy, sz, facet, el)
Compute the average normal for a facet of an element.
Definition symmetry.f90:154
Module with things related to the simulation time.
Implements a n-tuple.
Definition tuple.f90:34
Defines a zero-valued Dirichlet boundary condition.
Base type for a boundary condition.
Definition bc.f90:61
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
Definition coef.f90:55
Generic Dirichlet boundary condition on .
Definition dirichlet.f90:48
Mixed Dirichlet-Neumann symmetry plane condition.
Definition symmetry.f90:50
A struct that contains all info about the time, expand as needed.
Integer based 2-tuple.
Definition tuple.f90:51
Zero-valued Dirichlet boundary condition. Used for no-slip walls, but also for various auxillary cond...