Loading [MathJax]/extensions/tex2jax.js
Neko 0.9.99
A portable framework for high-order spectral element flow simulations
All Classes Namespaces Files Functions Variables Typedefs Enumerator Macros Pages
neumann.f90
Go to the documentation of this file.
1! Copyright (c) 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!
34module neumann
35 use num_types, only : rp
36 use bc, only : bc_t
37 use, intrinsic :: iso_c_binding, only : c_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
43 use math, only : cfill, copy
44 implicit none
45 private
46
51 type, public, extends(bc_t) :: neumann_t
52 real(kind=rp), allocatable :: flux_(:)
53 real(kind=rp), private :: init_flux_
54 contains
55 procedure, pass(this) :: apply_scalar => neumann_apply_scalar
56 procedure, pass(this) :: apply_vector => neumann_apply_vector
57 procedure, pass(this) :: apply_scalar_dev => neumann_apply_scalar_dev
58 procedure, pass(this) :: apply_vector_dev => neumann_apply_vector_dev
60 procedure, pass(this) :: init => neumann_init
62 procedure, pass(this) :: init_from_components => &
64 procedure, pass(this) :: flux => neumann_flux
66 procedure, pass(this) :: set_flux_scalar => neumann_set_flux_scalar
68 procedure, pass(this) :: set_flux_array => neumann_set_flux_array
70 generic :: set_flux => set_flux_scalar, set_flux_array
72 procedure, pass(this) :: free => neumann_free
74 procedure, pass(this) :: finalize => neumann_finalize
75 end type neumann_t
76
77contains
78
82 subroutine neumann_init(this, coef, json)
83 class(neumann_t), intent(inout), target :: this
84 type(coef_t), intent(in) :: coef
85 type(json_file), intent(inout) :: json
86 real(kind=rp) :: flux
87
88 call this%init_base(coef)
89 this%strong = .false.
90
91 call json_get(json, "flux", flux)
92 this%init_flux_ = flux
93 end subroutine neumann_init
94
98 subroutine neumann_init_from_components(this, coef, flux)
99 class(neumann_t), intent(inout), target :: this
100 type(coef_t), intent(in) :: coef
101 real(kind=rp), intent(in) :: flux
102
103 call this%init_base(coef)
104 this%init_flux_ = flux
105 end subroutine neumann_init_from_components
106
109 subroutine neumann_apply_scalar(this, x, n, t, tstep, strong)
110 class(neumann_t), intent(inout) :: this
111 integer, intent(in) :: n
112 real(kind=rp), intent(inout), dimension(n) :: x
113 real(kind=rp), intent(in), optional :: t
114 integer, intent(in), optional :: tstep
115 logical, intent(in), optional :: strong
116 integer :: i, m, k, facet
117 ! Store non-linear index
118 integer :: idx(4)
119 logical :: strong_ = .true.
120
121 if (present(strong)) strong_ = strong
122
123 m = this%msk(0)
124 if (.not. strong_) then
125 do i = 1, m
126 k = this%msk(i)
127 facet = this%facet(i)
128 idx = nonlinear_index(k, this%coef%Xh%lx, this%coef%Xh%lx,&
129 this%coef%Xh%lx)
130 select case (facet)
131 case (1,2)
132 x(k) = x(k) + this%flux_(i)*this%coef%area(idx(2), idx(3), facet,&
133 idx(4))
134 case (3,4)
135 x(k) = x(k) + this%flux_(i)*this%coef%area(idx(1), idx(3), facet,&
136 idx(4))
137 case (5,6)
138 x(k) = x(k) + this%flux_(i)*this%coef%area(idx(1), idx(2), facet,&
139 idx(4))
140 end select
141 end do
142 end if
143 end subroutine neumann_apply_scalar
144
147 subroutine neumann_apply_vector(this, x, y, z, n, t, tstep, strong)
148 class(neumann_t), intent(inout) :: this
149 integer, intent(in) :: n
150 real(kind=rp), intent(inout), dimension(n) :: x
151 real(kind=rp), intent(inout), dimension(n) :: y
152 real(kind=rp), intent(inout), dimension(n) :: z
153 real(kind=rp), intent(in), optional :: t
154 integer, intent(in), optional :: tstep
155 logical, intent(in), optional :: strong
156
157 call neko_error("Neumann bc not implemented for vectors")
158
159 end subroutine neumann_apply_vector
160
163 subroutine neumann_apply_scalar_dev(this, x_d, t, tstep, strong)
164 class(neumann_t), intent(inout), target :: this
165 type(c_ptr) :: x_d
166 real(kind=rp), intent(in), optional :: t
167 integer, intent(in), optional :: tstep
168 logical, intent(in), optional :: strong
169
170 call neko_error("Neumann bc not implemented on the device")
171
172 end subroutine neumann_apply_scalar_dev
173
176 subroutine neumann_apply_vector_dev(this, x_d, y_d, z_d, t, tstep, strong)
177 class(neumann_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 logical, intent(in), optional :: strong
184
185 call neko_error("Neumann bc not implemented on the device")
186
187 end subroutine neumann_apply_vector_dev
188
190 subroutine neumann_free(this)
191 class(neumann_t), target, intent(inout) :: this
192
193 call this%free_base()
194
195 end subroutine neumann_free
196
198 subroutine neumann_finalize(this)
199 class(neumann_t), target, intent(inout) :: this
200
201 call this%finalize_base()
202 allocate(this%flux_(this%msk(0)))
203
204 call cfill(this%flux_, this%init_flux_, this%msk(0))
205 end subroutine neumann_finalize
206
208 pure function neumann_flux(this) result(flux)
209 class(neumann_t), intent(in) :: this
210 real(kind=rp) :: flux(this%msk(0))
211
212 flux = this%flux_
213 end function neumann_flux
214
216 subroutine neumann_set_flux_scalar(this, flux)
217 class(neumann_t), intent(inout) :: this
218 real(kind=rp), intent(in) :: flux
219
220 this%flux_ = flux
221 end subroutine neumann_set_flux_scalar
222
225 subroutine neumann_set_flux_array(this, flux)
226 class(neumann_t), intent(inout) :: this
227 real(kind=rp), intent(in) :: flux(this%msk(0))
228
229 call copy(this%flux_, flux, this%msk(0))
230 end subroutine neumann_set_flux_array
231
232end module neumann
__device__ void nonlinear_index(const int idx, const int lx, int *index)
Retrieves a parameter by name or throws an error.
Defines a boundary condition.
Definition bc.f90:34
Coefficients.
Definition coef.f90:34
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:347
subroutine, public copy(a, b, n)
Copy a vector .
Definition math.f90:238
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:226
subroutine neumann_init_from_components(this, coef, flux)
Constructor from components.
Definition neumann.f90:99
subroutine neumann_init(this, coef, json)
Constructor.
Definition neumann.f90:83
subroutine neumann_apply_vector_dev(this, x_d, y_d, z_d, t, tstep, strong)
Boundary condition apply for a generic Neumann condition to vectors x, y and z (device version)
Definition neumann.f90:177
subroutine neumann_apply_scalar(this, x, n, t, tstep, strong)
Boundary condition apply for a generic Neumann condition to a vector x.
Definition neumann.f90:110
subroutine neumann_set_flux_scalar(this, flux)
Set the flux using a scalar.
Definition neumann.f90:217
pure real(kind=rp) function, dimension(this%msk(0)) neumann_flux(this)
Get the flux.
Definition neumann.f90:209
subroutine neumann_apply_scalar_dev(this, x_d, t, tstep, strong)
Boundary condition apply for a generic Neumann condition to a vector x (device version)
Definition neumann.f90:164
subroutine neumann_free(this)
Destructor.
Definition neumann.f90:191
subroutine neumann_finalize(this)
Finalize by setting the flux.
Definition neumann.f90:199
subroutine neumann_apply_vector(this, x, y, z, n, t, tstep, strong)
Boundary condition apply for a generic Neumann condition to vectors x, y and z.
Definition neumann.f90:148
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
Utilities.
Definition utils.f90:35
Base type for a boundary condition.
Definition bc.f90:57
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:51