Neko  0.9.0
A portable framework for high-order spectral element flow simulations
non_normal.f90
Go to the documentation of this file.
1 ! Copyright (c) 2020-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 !
34 module non_normal
35  use symmetry, only : symmetry_t
36  use num_types, only : rp
37  use tuple, only : tuple_i4_t
38  use coefs, only : coef_t
39  use, intrinsic :: iso_c_binding
40  implicit none
41  private
42 
44  type, public, extends(symmetry_t) :: non_normal_t
45  contains
47  procedure, pass(this) :: init => non_normal_init
49  procedure, pass(this) :: free => non_normal_free
50  end type non_normal_t
51 
52 contains
53 
55  subroutine non_normal_init(this, coef)
56  class(non_normal_t), intent(inout) :: this
57  type(coef_t), target, intent(in) :: coef
58  integer :: i, j, l
59  type(tuple_i4_t), pointer :: bfp(:)
60  real(kind=rp) :: sx, sy, sz
61  real(kind=rp), parameter :: tol = 1d-3
62  type(tuple_i4_t) :: bc_facet
63  integer :: facet, el
64 
65  call this%bc_x%free()
66  call this%bc_y%free()
67  call this%bc_z%free()
68  call this%bc_x%init_base(this%coef)
69  call this%bc_y%init_base(this%coef)
70  call this%bc_z%init_base(this%coef)
71 
72  associate(c => this%coef, nx => this%coef%nx, ny => this%coef%ny, &
73  nz => this%coef%nz)
74  bfp => this%marked_facet%array()
75  do i = 1, this%marked_facet%size()
76  bc_facet = bfp(i)
77  facet = bc_facet%x(1)
78  el = bc_facet%x(2)
79  sx = 0d0
80  sy = 0d0
81  sz = 0d0
82  select case (facet)
83  case (1, 2)
84  do l = 2, c%Xh%lx - 1
85  do j = 2, c%Xh%lx -1
86  sx = sx + abs(abs(nx(l, j, facet, el)) - 1d0)
87  sy = sy + abs(abs(ny(l, j, facet, el)) - 1d0)
88  sz = sz + abs(abs(nz(l, j, facet, el)) - 1d0)
89  end do
90  end do
91  case (3, 4)
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 (5, 6)
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  end select
108  sx = sx / (c%Xh%lx - 2)**2
109  sy = sy / (c%Xh%lx - 2)**2
110  sz = sz / (c%Xh%lx - 2)**2
111 
112  if (sx .lt. tol) then
113  call this%bc_y%mark_facet(facet, el)
114  call this%bc_z%mark_facet(facet, el)
115  end if
116 
117  if (sy .lt. tol) then
118  call this%bc_x%mark_facet(facet, el)
119  call this%bc_z%mark_facet(facet, el)
120  end if
121 
122  if (sz .lt. tol) then
123  call this%bc_y%mark_facet(facet, el)
124  call this%bc_x%mark_facet(facet, el)
125  end if
126  end do
127  end associate
128  call this%bc_x%finalize()
129  call this%bc_x%set_g(0.0_rp)
130  call this%bc_y%finalize()
131  call this%bc_y%set_g(0.0_rp)
132  call this%bc_z%finalize()
133  call this%bc_z%set_g(0.0_rp)
134  end subroutine non_normal_init
135 
137  subroutine non_normal_free(this)
138  class(non_normal_t), target, intent(inout) :: this
139 
140  call this%symmetry_t%free()
141  end subroutine non_normal_free
142 end module non_normal
Coefficients.
Definition: coef.f90:34
Dirichlet condition on axis aligned plane in the non normal direction.
Definition: non_normal.f90:34
subroutine non_normal_init(this, coef)
Constructor.
Definition: non_normal.f90:56
subroutine non_normal_free(this)
Destructor.
Definition: non_normal.f90:138
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
Implements a n-tuple.
Definition: tuple.f90:34
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
Definition: coef.f90:55
Dirichlet condition in non normal direction of a plane.
Definition: non_normal.f90:44
Mixed Dirichlet-Neumann symmetry plane condition.
Definition: symmetry.f90:46
Integer based 2-tuple.
Definition: tuple.f90:51