Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
neumann.f90
Go to the documentation of this file.
1! Copyright (c) 2024-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!
34module neumann
35 use num_types, only : rp
36 use bc, only : bc_t
37 use, intrinsic :: iso_c_binding, only : c_ptr, c_null_ptr
39 use coefs, only : coef_t
40 use json_module, only : json_file
41 use json_utils, only : json_get
42 use math, only : cfill, copy, abscmp
43 use vector, only: vector_t
48 use time_state, only : time_state_t
49 implicit none
50 private
51
56 type, public, extends(bc_t) :: neumann_t
57 type(vector_t) :: flux_
58 real(kind=rp), private :: init_flux_
59 logical :: uniform_0 = .false.
60 contains
61 procedure, pass(this) :: apply_scalar => neumann_apply_scalar
62 procedure, pass(this) :: apply_vector => neumann_apply_vector
63 procedure, pass(this) :: apply_scalar_dev => neumann_apply_scalar_dev
64 procedure, pass(this) :: apply_vector_dev => neumann_apply_vector_dev
66 procedure, pass(this) :: init => neumann_init
68 procedure, pass(this) :: init_from_components => &
70 procedure, pass(this) :: flux => neumann_flux
72 procedure, pass(this) :: set_flux_scalar => neumann_set_flux_scalar
74 procedure, pass(this) :: set_flux_array => neumann_set_flux_array
76 generic :: set_flux => set_flux_scalar, set_flux_array
78 procedure, pass(this) :: free => neumann_free
80 procedure, pass(this) :: finalize => neumann_finalize
81 end type neumann_t
82
83contains
84
88 subroutine neumann_init(this, coef, json)
89 class(neumann_t), intent(inout), target :: this
90 type(coef_t), target, intent(in) :: coef
91 type(json_file), intent(inout) :: json
92 real(kind=rp) :: flux
93
94 call this%init_base(coef)
95 this%strong = .false.
96
97 call json_get(json, "flux", flux)
98 this%init_flux_ = flux
99 end subroutine neumann_init
100
104 subroutine neumann_init_from_components(this, coef, flux)
105 class(neumann_t), intent(inout), target :: this
106 type(coef_t), intent(in) :: coef
107 real(kind=rp), intent(in) :: flux
108
109 call this%init_base(coef)
110 this%init_flux_ = flux
111 end subroutine neumann_init_from_components
112
115 subroutine neumann_apply_scalar(this, x, n, time, strong)
116 class(neumann_t), intent(inout) :: this
117 integer, intent(in) :: n
118 real(kind=rp), intent(inout), dimension(n) :: x
119 type(time_state_t), intent(in), optional :: time
120 logical, intent(in), optional :: strong
121 integer :: i, m, k, facet
122 ! Store non-linear index
123 integer :: idx(4)
124 logical :: strong_
125
126 if (present(strong)) then
127 strong_ = strong
128 else
129 strong_ = .true.
130 end if
131
132 m = this%msk(0)
133 if (.not. strong_) then
134 do i = 1, m
135 k = this%msk(i)
136 facet = this%facet(i)
137 idx = nonlinear_index(k, this%coef%Xh%lx, this%coef%Xh%lx,&
138 this%coef%Xh%lx)
139 select case (facet)
140 case (1,2)
141 x(k) = x(k) + &
142 this%flux_%x(i)*this%coef%area(idx(2), idx(3), facet, idx(4))
143 case (3,4)
144 x(k) = x(k) + &
145 this%flux_%x(i)*this%coef%area(idx(1), idx(3), facet, idx(4))
146 case (5,6)
147 x(k) = x(k) + &
148 this%flux_%x(i)*this%coef%area(idx(1), idx(2), facet, idx(4))
149 end select
150 end do
151 end if
152 end subroutine neumann_apply_scalar
153
156 subroutine neumann_apply_vector(this, x, y, z, n, time, strong)
157 class(neumann_t), intent(inout) :: this
158 integer, intent(in) :: n
159 real(kind=rp), intent(inout), dimension(n) :: x
160 real(kind=rp), intent(inout), dimension(n) :: y
161 real(kind=rp), intent(inout), dimension(n) :: z
162 type(time_state_t), intent(in), optional :: time
163 logical, intent(in), optional :: strong
164
165 if (.not. this%uniform_0 .and. this%msk(0) .gt. 0) then
166 call neko_error("Neumann bc not implemented for vectors")
167 end if
168 end subroutine neumann_apply_vector
169
172 subroutine neumann_apply_scalar_dev(this, x_d, time, strong, strm)
173 class(neumann_t), intent(inout), target :: this
174 type(c_ptr), intent(inout) :: x_d
175 type(time_state_t), intent(in), optional :: time
176 logical, intent(in), optional :: strong
177 type(c_ptr), intent(inout) :: strm
178 logical :: strong_
179
180 if (present(strong)) then
181 strong_ = strong
182 else
183 strong_ = .true.
184 end if
185
186 if (.not. this%uniform_0 .and. this%msk(0) .gt. 0 .and. &
187 .not. strong_) then
188 call device_neumann_apply_scalar(this%msk_d, this%facet_d, x_d, &
189 this%flux_%x_d, this%coef%area_d, this%coef%Xh%lx, &
190 size(this%msk), strm)
191 end if
192 end subroutine neumann_apply_scalar_dev
193
196 subroutine neumann_apply_vector_dev(this, x_d, y_d, z_d, &
197 time, strong, strm)
198 class(neumann_t), intent(inout), target :: this
199 type(c_ptr), intent(inout) :: x_d
200 type(c_ptr), intent(inout) :: y_d
201 type(c_ptr), intent(inout) :: z_d
202 type(time_state_t), intent(in), optional :: time
203 logical, intent(in), optional :: strong
204 type(c_ptr), intent(inout) :: strm
205
206 if (.not. this%uniform_0 .and. this%msk(0) .gt. 0) then
207 call neko_error("Neumann bc not implemented for vectors.")
208 end if
209
210 end subroutine neumann_apply_vector_dev
211
213 subroutine neumann_free(this)
214 class(neumann_t), target, intent(inout) :: this
215
216 call this%free_base()
217
218 end subroutine neumann_free
219
221 subroutine neumann_finalize(this, only_facets)
222 class(neumann_t), target, intent(inout) :: this
223 logical, optional, intent(in) :: only_facets
224 integer :: i
225
226 if (present(only_facets)) then
227 if (only_facets .eqv. .false.) then
228 call neko_error("For neumann_t, only_facets has to be true.")
229 end if
230 end if
231
232 call this%finalize_base(.true.)
233
234 call this%flux_%init(this%msk(0))
235
236 this%flux_ = this%init_flux_
237 this%uniform_0 = .true.
238
239 do i = 1,this%msk(0)
240 this%uniform_0 = abscmp(this%init_flux_, 0.0_rp) .and. this%uniform_0
241 end do
242 end subroutine neumann_finalize
243
245 !This looks really sketchy IMO /Martin
246 pure function neumann_flux(this) result(flux)
247 class(neumann_t), intent(in) :: this
248 real(kind=rp) :: flux(this%msk(0))
249
250 flux = this%flux_%x
251 end function neumann_flux
252
255 subroutine neumann_set_flux_scalar(this, flux)
256 class(neumann_t), intent(inout) :: this
257 real(kind=rp), intent(in) :: flux
258
259 this%flux_ = flux
260 this%uniform_0 = abscmp(flux, 0.0_rp)
261
262 end subroutine neumann_set_flux_scalar
263
266 subroutine neumann_set_flux_array(this, flux)
267 class(neumann_t), intent(inout) :: this
268 type(vector_t), intent(in) :: flux
269 integer :: i
270
271 this%flux_ = flux
272
273 ! We assume that passing a homogeneous array is so rare that it is not
274 ! worth doing a check, which requires a loop over the flux values.
275 ! Note that this routine can be called at each timestep, e.g. by wall
276 ! models.
277 this%uniform_0 = .false.
278
279 end subroutine neumann_set_flux_array
280
281end module neumann
__inline__ __device__ void nonlinear_index(const int idx, const int lx, int *index)
Definition bc_utils.h:44
Copy data between host and device (or device and device)
Definition device.F90:66
Retrieves a parameter by name or throws an error.
Defines a boundary condition.
Definition bc.f90:34
Coefficients.
Definition coef.f90:34
subroutine, public device_copy(a_d, b_d, n, strm)
Copy a vector .
subroutine, public device_cfill(a_d, c, n, strm)
Set all elements to a constant c .
subroutine, public device_neumann_apply_scalar(msk, facet, x, flux, area, lx, m, strm)
Device abstraction, common interface for various accelerators.
Definition device.F90:34
integer, parameter, public device_to_host
Definition device.F90:47
Utilities for retrieving parameters from the case files.
Definition math.f90:60
subroutine, public cfill(a, c, n)
Set all elements to a constant c .
Definition math.f90:493
subroutine, public copy(a, b, n)
Copy a vector .
Definition math.f90:255
Build configurations.
integer, parameter neko_bcknd_device
Defines a Neumann boundary condition.
Definition neumann.f90:34
subroutine neumann_set_flux_array(this, flux)
Set the flux using an array of values.
Definition neumann.f90:267
subroutine neumann_apply_vector_dev(this, x_d, y_d, z_d, time, strong, strm)
Boundary condition apply for a generic Neumann condition to vectors x, y and z (device version)
Definition neumann.f90:198
subroutine neumann_init_from_components(this, coef, flux)
Constructor from components.
Definition neumann.f90:105
subroutine neumann_finalize(this, only_facets)
Finalize by setting the flux.
Definition neumann.f90:222
subroutine neumann_init(this, coef, json)
Constructor.
Definition neumann.f90:89
subroutine neumann_apply_vector(this, x, y, z, n, time, strong)
Boundary condition apply for a generic Neumann condition to vectors x, y and z.
Definition neumann.f90:157
subroutine neumann_apply_scalar(this, x, n, time, strong)
Boundary condition apply for a generic Neumann condition to a vector x.
Definition neumann.f90:116
subroutine neumann_set_flux_scalar(this, flux)
Set the flux using a scalar.
Definition neumann.f90:256
pure real(kind=rp) function, dimension(this%msk(0)) neumann_flux(this)
Get the flux.
Definition neumann.f90:247
subroutine neumann_apply_scalar_dev(this, x_d, time, strong, strm)
Boundary condition apply for a generic Neumann condition to a vector x (device version)
Definition neumann.f90:173
subroutine neumann_free(this)
Destructor.
Definition neumann.f90:214
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
Module with things related to the simulation time.
Utilities.
Definition utils.f90:35
Defines a vector.
Definition vector.f90:34
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 Neumann boundary condition for scalar fields. Sets the flux of the field to the chosen value.
Definition neumann.f90:56
A struct that contains all info about the time, expand as needed.