Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
zero_dirichlet.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 num_types, only : rp
37 use bc, only : bc_t
38 use, intrinsic :: iso_c_binding, only : c_ptr
39 use coefs, only : coef_t
40 use json_module, only : json_file
41 use time_state, only : time_state_t
42 implicit none
43 private
44
48 type, public, extends(bc_t) :: zero_dirichlet_t
49 contains
50 procedure, pass(this) :: apply_scalar => zero_dirichlet_apply_scalar
51 procedure, pass(this) :: apply_vector => zero_dirichlet_apply_vector
52 procedure, pass(this) :: apply_scalar_dev => &
54 procedure, pass(this) :: apply_vector_dev => &
57 procedure, pass(this) :: init => zero_dirichlet_init
59 procedure, pass(this) :: init_from_components => &
62 procedure, pass(this) :: free => zero_dirichlet_free
64 procedure, pass(this) :: finalize => zero_dirichlet_finalize
65 end type zero_dirichlet_t
66
67contains
68
72 subroutine zero_dirichlet_init(this, coef, json)
73 class(zero_dirichlet_t), intent(inout), target :: this
74 type(coef_t), target, intent(in) :: coef
75 type(json_file), intent(inout) :: json
76
77 call this%init_from_components(coef)
78 end subroutine zero_dirichlet_init
79
83 class(zero_dirichlet_t), intent(inout), target :: this
84 type(coef_t), target, intent(in) :: coef
85
86 call this%init_base(coef)
88
91 subroutine zero_dirichlet_apply_scalar(this, x, n, time, strong)
92 class(zero_dirichlet_t), intent(inout) :: this
93 integer, intent(in) :: n
94 real(kind=rp), intent(inout), dimension(n) :: x
95 type(time_state_t), intent(in), optional :: time
96 logical, intent(in), optional :: strong
97 integer :: i, m, k
98 logical :: strong_
99
100 if (present(strong)) then
101 strong_ = strong
102 else
103 strong_ = .true.
104 end if
105
106 m = this%msk(0)
107
108 if (strong_) then
109 do i = 1, m
110 k = this%msk(i)
111 x(k) = 0d0
112 end do
113 end if
114
115 end subroutine zero_dirichlet_apply_scalar
116
118 subroutine zero_dirichlet_apply_vector(this, x, y, z, n, time, strong)
119 class(zero_dirichlet_t), intent(inout) :: this
120 integer, intent(in) :: n
121 real(kind=rp), intent(inout), dimension(n) :: x
122 real(kind=rp), intent(inout), dimension(n) :: y
123 real(kind=rp), intent(inout), dimension(n) :: z
124 type(time_state_t), intent(in), optional :: time
125 logical, intent(in), optional :: strong
126 integer :: i, m, k
127 logical :: strong_
128
129 if (present(strong)) then
130 strong_ = strong
131 else
132 strong_ = .true.
133 end if
134
135 if (strong_) then
136 m = this%msk(0)
137 do i = 1, m
138 k = this%msk(i)
139 x(k) = 0d0
140 y(k) = 0d0
141 z(k) = 0d0
142 end do
143 end if
144
145 end subroutine zero_dirichlet_apply_vector
146
148 subroutine zero_dirichlet_apply_scalar_dev(this, x_d, time, strong, strm)
149 class(zero_dirichlet_t), intent(inout), target :: this
150 type(c_ptr), intent(inout) :: x_d
151 type(time_state_t), intent(in), optional :: time
152 logical, intent(in), optional :: strong
153 type(c_ptr), intent(inout) :: strm
154 logical :: strong_
155
156 if (present(strong)) then
157 strong_ = strong
158 else
159 strong_ = .true.
160 end if
161
162 if (strong_ .and. (this%msk(0) .gt. 0)) then
163 call device_zero_dirichlet_apply_scalar(this%msk_d, x_d, &
164 size(this%msk), strm)
165 end if
166
168
170 subroutine zero_dirichlet_apply_vector_dev(this, x_d, y_d, z_d, time, &
171 strong, strm)
172 class(zero_dirichlet_t), intent(inout), target :: this
173 type(c_ptr), intent(inout) :: x_d
174 type(c_ptr), intent(inout) :: y_d
175 type(c_ptr), intent(inout) :: z_d
176 type(time_state_t), intent(in), optional :: time
177 logical, intent(in), optional :: strong
178 type(c_ptr), intent(inout) :: strm
179 logical :: strong_
180
181 if (present(strong)) then
182 strong_ = strong
183 else
184 strong_ = .true.
185 end if
186
187 if (strong_ .and. (this%msk(0) .gt. 0)) then
188 call device_zero_dirichlet_apply_vector(this%msk_d, x_d, y_d, z_d, &
189 size(this%msk), strm)
190 end if
191
193
195 subroutine zero_dirichlet_free(this)
196 class(zero_dirichlet_t), target, intent(inout) :: this
197
198 call this%free_base()
199
200 end subroutine zero_dirichlet_free
201
203 subroutine zero_dirichlet_finalize(this, only_facets)
204 class(zero_dirichlet_t), target, intent(inout) :: this
205 logical, optional, intent(in) :: only_facets
206 logical :: only_facets_
207
208 if (present(only_facets)) then
209 only_facets_ = only_facets
210 else
211 only_facets_ = .false.
212 end if
213
214 call this%finalize_base(only_facets_)
215 end subroutine zero_dirichlet_finalize
216
217end module zero_dirichlet
Defines a boundary condition.
Definition bc.f90:34
Coefficients.
Definition coef.f90:34
subroutine, public device_zero_dirichlet_apply_vector(msk, x, y, z, m, strm)
subroutine, public device_zero_dirichlet_apply_scalar(msk, x, m, strm)
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
Module with things related to the simulation time.
Defines a zero-valued Dirichlet boundary condition.
subroutine zero_dirichlet_apply_scalar(this, x, n, time, strong)
Apply boundary condition to a scalar field. to a vector x.
subroutine zero_dirichlet_apply_vector_dev(this, x_d, y_d, z_d, time, strong, strm)
Apply boundary condition to a vector field, device version.
subroutine zero_dirichlet_init(this, coef, json)
Constructor.
subroutine zero_dirichlet_apply_scalar_dev(this, x_d, time, strong, strm)
Apply boundary condition to a scalar field, device version.
subroutine zero_dirichlet_free(this)
Destructor.
subroutine zero_dirichlet_init_from_components(this, coef)
Constructor.
subroutine zero_dirichlet_finalize(this, only_facets)
Finalize.
subroutine zero_dirichlet_apply_vector(this, x, y, z, n, time, strong)
Apply boundary condition to a vector field.
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
A struct that contains all info about the time, expand as needed.
Zero-valued Dirichlet boundary condition. Used for no-slip walls, but also for various auxillary cond...