Neko  0.9.99
A portable framework for high-order spectral element flow simulations
symmetry.f90
Go to the documentation of this file.
1 ! Copyright (c) 2020-2021, 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 !
34 module symmetry
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, intrinsic :: iso_c_binding, only : c_ptr
42  implicit none
43  private
44 
46  type, public, extends(bc_t) :: symmetry_t
47  type(dirichlet_t) :: bc_x
48  type(dirichlet_t) :: bc_y
49  type(dirichlet_t) :: bc_z
50  contains
51  procedure, pass(this) :: init => symmetry_init
52  procedure, pass(this) :: apply_scalar => symmetry_apply_scalar
53  procedure, pass(this) :: apply_vector => symmetry_apply_vector
54  procedure, pass(this) :: apply_scalar_dev => symmetry_apply_scalar_dev
55  procedure, pass(this) :: apply_vector_dev => symmetry_apply_vector_dev
57  procedure, pass(this) :: free => symmetry_free
58  end type symmetry_t
59 
60 contains
61 
63  subroutine symmetry_init(this, coef)
64  class(symmetry_t), intent(inout) :: this
65  type(coef_t), target, intent(in) :: coef
66  integer :: i, j, l
67  type(tuple_i4_t), pointer :: bfp(:)
68  real(kind=rp) :: sx, sy, sz
69  real(kind=rp), parameter :: tol = 1d-3
70  type(tuple_i4_t) :: bc_facet
71  integer :: facet, el
72 
73  call this%bc_x%free()
74  call this%bc_y%free()
75  call this%bc_z%free()
76  call this%bc_x%init_base(this%coef)
77  call this%bc_y%init_base(this%coef)
78  call this%bc_z%init_base(this%coef)
79 
80  associate(c => this%coef, nx => this%coef%nx, ny => this%coef%ny, &
81  nz => this%coef%nz)
82  bfp => this%marked_facet%array()
83  do i = 1, this%marked_facet%size()
84  bc_facet = bfp(i)
85  facet = bc_facet%x(1)
86  el = bc_facet%x(2)
87  sx = 0d0
88  sy = 0d0
89  sz = 0d0
90  select case (facet)
91  case (1, 2)
92  do l = 2, c%Xh%lx - 1
93  do j = 2, c%Xh%lx -1
94  sx = sx + abs(abs(nx(l, j, facet, el)) - 1d0)
95  sy = sy + abs(abs(ny(l, j, facet, el)) - 1d0)
96  sz = sz + abs(abs(nz(l, j, facet, el)) - 1d0)
97  end do
98  end do
99  case (3, 4)
100  do l = 2, c%Xh%lx - 1
101  do j = 2, c%Xh%lx - 1
102  sx = sx + abs(abs(nx(l, j, facet, el)) - 1d0)
103  sy = sy + abs(abs(ny(l, j, facet, el)) - 1d0)
104  sz = sz + abs(abs(nz(l, j, facet, el)) - 1d0)
105  end do
106  end do
107  case (5, 6)
108  do l = 2, c%Xh%lx - 1
109  do j = 2, c%Xh%lx - 1
110  sx = sx + abs(abs(nx(l, j, facet, el)) - 1d0)
111  sy = sy + abs(abs(ny(l, j, facet, el)) - 1d0)
112  sz = sz + abs(abs(nz(l, j, facet, el)) - 1d0)
113  end do
114  end do
115  end select
116  sx = sx / (c%Xh%lx - 2)**2
117  sy = sy / (c%Xh%lx - 2)**2
118  sz = sz / (c%Xh%lx - 2)**2
119 
120  if (sx .lt. tol) then
121  call this%bc_x%mark_facet(facet, el)
122  end if
123 
124  if (sy .lt. tol) then
125  call this%bc_y%mark_facet(facet, el)
126  end if
127 
128  if (sz .lt. tol) then
129  call this%bc_z%mark_facet(facet, el)
130  end if
131  end do
132  end associate
133  call this%bc_x%finalize()
134  call this%bc_x%set_g(0.0_rp)
135  call this%bc_y%finalize()
136  call this%bc_y%set_g(0.0_rp)
137  call this%bc_z%finalize()
138  call this%bc_z%set_g(0.0_rp)
139 
140  end subroutine symmetry_init
141 
143  subroutine symmetry_apply_scalar(this, x, n, t, tstep)
144  class(symmetry_t), intent(inout) :: this
145  integer, intent(in) :: n
146  real(kind=rp), intent(inout), dimension(n) :: x
147  real(kind=rp), intent(in), optional :: t
148  integer, intent(in), optional :: tstep
149  end subroutine symmetry_apply_scalar
150 
152  subroutine symmetry_apply_vector(this, x, y, z, n, t, tstep)
153  class(symmetry_t), intent(inout) :: this
154  integer, intent(in) :: n
155  real(kind=rp), intent(inout), dimension(n) :: x
156  real(kind=rp), intent(inout), dimension(n) :: y
157  real(kind=rp), intent(inout), dimension(n) :: z
158  real(kind=rp), intent(in), optional :: t
159  integer, intent(in), optional :: tstep
160 
161  call this%bc_x%apply_scalar(x, n)
162  call this%bc_y%apply_scalar(y, n)
163  call this%bc_z%apply_scalar(z, n)
164 
165  end subroutine symmetry_apply_vector
166 
168  subroutine symmetry_apply_scalar_dev(this, x_d, t, tstep)
169  class(symmetry_t), intent(inout), target :: this
170  type(c_ptr) :: x_d
171  real(kind=rp), intent(in), optional :: t
172  integer, intent(in), optional :: tstep
173  end subroutine symmetry_apply_scalar_dev
174 
176  subroutine symmetry_apply_vector_dev(this, x_d, y_d, z_d, t, tstep)
177  class(symmetry_t), intent(inout), target :: this
178  type(c_ptr) :: x_d
179  type(c_ptr) :: y_d
180  type(c_ptr) :: z_d
181  real(kind=rp), intent(in), optional :: t
182  integer, intent(in), optional :: tstep
183 
184  call device_symmetry_apply_vector(this%bc_x%msk_d, this%bc_y%msk_d, &
185  this%bc_z%msk_d, x_d, y_d, z_d, &
186  this%bc_x%msk(0), &
187  this%bc_y%msk(0), &
188  this%bc_z%msk(0))
189 
190  end subroutine symmetry_apply_vector_dev
191 
193  subroutine symmetry_free(this)
194  class(symmetry_t), target, intent(inout) :: this
195 
196  call this%free_base()
197  call this%bc_x%free()
198  call this%bc_y%free()
199  call this%bc_z%free()
200 
201  end subroutine symmetry_free
202 end 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)
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_init(this, coef)
Initialize symmetry mask for each axis.
Definition: symmetry.f90:64
subroutine symmetry_apply_vector(this, x, y, z, n, t, tstep)
Apply symmetry conditions (axis aligned)
Definition: symmetry.f90:153
subroutine symmetry_apply_scalar(this, x, n, t, tstep)
No-op scalar apply.
Definition: symmetry.f90:144
subroutine symmetry_apply_scalar_dev(this, x_d, t, tstep)
No-op scalar apply (device version)
Definition: symmetry.f90:169
subroutine symmetry_free(this)
Destructor.
Definition: symmetry.f90:194
subroutine symmetry_apply_vector_dev(this, x_d, y_d, z_d, t, tstep)
Apply symmetry conditions (axis aligned) (device version)
Definition: symmetry.f90:177
Implements a n-tuple.
Definition: tuple.f90:34
Base type for a boundary condition.
Definition: bc.f90:51
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:44
Mixed Dirichlet-Neumann symmetry plane condition.
Definition: symmetry.f90:46
Integer based 2-tuple.
Definition: tuple.f90:51