Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
dong_outflow.f90
Go to the documentation of this file.
1! Copyright (c) 2022-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!
35 use neko_config
36 use dirichlet, only : dirichlet_t
38 use num_types, only : rp, c_rp
39 use bc, only : bc_t
40 use field, only : field_t
41 use dofmap, only : dofmap_t
42 use coefs, only : coef_t
43 use utils, only : nonlinear_index
46 use, intrinsic :: iso_c_binding, only : c_ptr, c_sizeof, c_null_ptr
47 use json_module, only : json_file
49 use utils, only : neko_error
50 use time_state, only : time_state_t
51 implicit none
52 private
53
59 type, public, extends(bc_t) :: dong_outflow_t
60 type(field_t), pointer :: u
61 type(field_t), pointer :: v
62 type(field_t), pointer :: w
63 real(kind=rp) :: delta
64 real(kind=rp) :: uinf
65 type(c_ptr) :: normal_x_d = c_null_ptr
66 type(c_ptr) :: normal_y_d = c_null_ptr
67 type(c_ptr) :: normal_z_d = c_null_ptr
68 contains
69 procedure, pass(this) :: apply_scalar => dong_outflow_apply_scalar
70 procedure, pass(this) :: apply_vector => dong_outflow_apply_vector
71 procedure, pass(this) :: apply_scalar_dev => dong_outflow_apply_scalar_dev
72 procedure, pass(this) :: apply_vector_dev => dong_outflow_apply_vector_dev
74 procedure, pass(this) :: init => dong_outflow_init
76 procedure, pass(this) :: free => dong_outflow_free
78 procedure, pass(this) :: finalize => dong_outflow_finalize
79 end type dong_outflow_t
80
81contains
85 subroutine dong_outflow_init(this, coef, json)
86 class(dong_outflow_t), target, intent(inout) :: this
87 type(coef_t), target, intent(in) :: coef
88 type(json_file), intent(inout) :: json
89 call this%free()
90 call this%init_base(coef)
91
92 call json_get_or_default(json, 'delta', this%delta, 0.01_rp)
93 call json_get_or_default(json, 'velocity_scale', this%uinf, 1.0_rp)
94
95 end subroutine dong_outflow_init
96
99 subroutine dong_outflow_apply_scalar(this, x, n, time, strong)
100 class(dong_outflow_t), intent(inout) :: this
101 integer, intent(in) :: n
102 real(kind=rp), intent(inout), dimension(n) :: x
103 type(time_state_t), intent(in), optional :: time
104 logical, intent(in), optional :: strong
105 integer :: i, m, k, facet, idx(4)
106 real(kind=rp) :: vn, s0, ux, uy, uz, normal_xyz(3)
107 logical :: strong_
108
109 if (present(strong)) then
110 strong_ = strong
111 else
112 strong_ = .true.
113 end if
114
115 !Im actually not sure what to do if one has two dong that share a corner.
116 if (strong_) then
117 m = this%msk(0)
118 do i = 1, m
119 k = this%msk(i)
120 facet = this%facet(i)
121 ux = this%u%x(k,1,1,1)
122 uy = this%v%x(k,1,1,1)
123 uz = this%w%x(k,1,1,1)
124 idx = nonlinear_index(k, this%Xh%lx, this%Xh%lx, this%Xh%lx)
125 normal_xyz = this%coef%get_normal(idx(1), idx(2), idx(3), idx(4), &
126 facet)
127 vn = ux*normal_xyz(1) + uy*normal_xyz(2) + uz*normal_xyz(3)
128 s0 = 0.5_rp*(1.0_rp - tanh(vn / (this%uinf * this%delta)))
129
130 x(k) = -0.5*(ux*ux+uy*uy+uz*uz)*s0
131 end do
132 end if
133 end subroutine dong_outflow_apply_scalar
134
137 subroutine dong_outflow_apply_vector(this, x, y, z, n, time, strong)
138 class(dong_outflow_t), intent(inout) :: this
139 integer, intent(in) :: n
140 real(kind=rp), intent(inout), dimension(n) :: x
141 real(kind=rp), intent(inout), dimension(n) :: y
142 real(kind=rp), intent(inout), dimension(n) :: z
143 type(time_state_t), intent(in), optional :: time
144 logical, intent(in), optional :: strong
145
146 end subroutine dong_outflow_apply_vector
147
150 subroutine dong_outflow_apply_scalar_dev(this, x_d, time, strong, strm)
151 class(dong_outflow_t), intent(inout), target :: this
152 type(c_ptr), intent(inout) :: x_d
153 type(time_state_t), intent(in), optional :: time
154 logical, intent(in), optional :: strong
155 type(c_ptr), intent(inout) :: strm
156 logical :: strong_
157
158 if (present(strong)) then
159 strong_ = strong
160 else
161 strong_ = .true.
162 end if
163
164 if (strong_ .and. this%msk(0) .gt. 0) then
165 call device_dong_outflow_apply_scalar(this%msk_d, x_d, &
166 this%normal_x_d, this%normal_y_d, this%normal_z_d, &
167 this%u%x_d, this%v%x_d, this%w%x_d, &
168 this%uinf, this%delta, &
169 this%msk(0), strm)
170 end if
171
172 end subroutine dong_outflow_apply_scalar_dev
173
176 subroutine dong_outflow_apply_vector_dev(this, x_d, y_d, z_d, time, &
177 strong, strm)
178 class(dong_outflow_t), intent(inout), target :: this
179 type(c_ptr), intent(inout) :: x_d
180 type(c_ptr), intent(inout) :: y_d
181 type(c_ptr), intent(inout) :: z_d
182 type(time_state_t), intent(in), optional :: time
183 logical, intent(in), optional :: strong
184 type(c_ptr), intent(inout) :: strm
185
186 !call device_dong_outflow_apply_vector(this%msk_d, x_d, y_d, z_d, &
187 ! this%g, size(this%msk))
188
189 end subroutine dong_outflow_apply_vector_dev
190
192 subroutine dong_outflow_free(this)
193 class(dong_outflow_t), target, intent(inout) :: this
194
195 call this%free_base
196
197 end subroutine dong_outflow_free
198
200 subroutine dong_outflow_finalize(this, only_facets)
201 class(dong_outflow_t), target, intent(inout) :: this
202 logical, optional, intent(in) :: only_facets
203 real(kind=rp), allocatable :: temp_x(:)
204 real(kind=rp), allocatable :: temp_y(:)
205 real(kind=rp), allocatable :: temp_z(:)
206 real(c_rp) :: dummy
207 integer :: i, m, k, facet, idx(4)
208 real(kind=rp) :: normal_xyz(3)
209
210 if (present(only_facets)) then
211 if (only_facets .eqv. .false.) then
212 call neko_error("For dong_outflow_t, only_facets has to be true.")
213 end if
214 end if
215
216 call this%finalize_base(.true.)
217
218 this%u => neko_field_registry%get_field("u")
219 this%v => neko_field_registry%get_field("v")
220 this%w => neko_field_registry%get_field("w")
221 if ((neko_bcknd_device .eq. 1) .and. (this%msk(0) .gt. 0)) then
222 call device_alloc(this%normal_x_d, c_sizeof(dummy)*this%msk(0))
223 call device_alloc(this%normal_y_d, c_sizeof(dummy)*this%msk(0))
224 call device_alloc(this%normal_z_d, c_sizeof(dummy)*this%msk(0))
225 m = this%msk(0)
226 allocate(temp_x(m))
227 allocate(temp_y(m))
228 allocate(temp_z(m))
229 do i = 1, m
230 k = this%msk(i)
231 facet = this%facet(i)
232 idx = nonlinear_index(k, this%Xh%lx, this%Xh%lx, this%Xh%lx)
233 normal_xyz = &
234 this%coef%get_normal(idx(1), idx(2), idx(3), idx(4), facet)
235 temp_x(i) = normal_xyz(1)
236 temp_y(i) = normal_xyz(2)
237 temp_z(i) = normal_xyz(3)
238 end do
239 call device_memcpy(temp_x, this%normal_x_d, m, host_to_device, &
240 sync = .false.)
241 call device_memcpy(temp_y, this%normal_y_d, m, host_to_device, &
242 sync = .false.)
243 call device_memcpy(temp_z, this%normal_z_d, m, host_to_device, &
244 sync = .true.)
245 deallocate( temp_x, temp_y, temp_z)
246 end if
247 end subroutine dong_outflow_finalize
248
249end module dong_outflow
__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 assigns a provided default value. In the latter case also adds the m...
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_dong_outflow_apply_scalar(msk, x, normal_x, normal_y, normal_z, u, v, w, uinf, delta, m, strm)
Device abstraction, common interface for various accelerators.
Definition device.F90:34
integer, parameter, public host_to_device
Definition device.F90:47
subroutine, public device_alloc(x_d, s)
Allocate memory on the device.
Definition device.F90:187
Defines a dirichlet boundary condition.
Definition dirichlet.f90:34
Defines a mapping of the degrees of freedom.
Definition dofmap.f90:35
Defines a dong outflow condition.
subroutine dong_outflow_apply_vector_dev(this, x_d, y_d, z_d, time, strong, strm)
Boundary condition apply for a generic Dirichlet condition to vectors x, y and z (device version)
subroutine dong_outflow_apply_vector(this, x, y, z, n, time, strong)
Boundary condition apply for a generic Dirichlet condition to vectors x, y and z.
subroutine dong_outflow_apply_scalar_dev(this, x_d, time, strong, strm)
Boundary condition apply for a generic Dirichlet condition to a vector x (device version)
subroutine dong_outflow_free(this)
Destructor.
subroutine dong_outflow_apply_scalar(this, x, n, time, strong)
Boundary condition apply for a generic Dirichlet condition to a vector x.
subroutine dong_outflow_init(this, coef, json)
Constructor.
subroutine dong_outflow_finalize(this, only_facets)
Finalize.
Defines a registry for storing solution fields.
type(field_registry_t), target, public neko_field_registry
Global field registry.
Defines a field.
Definition field.f90:34
Utilities for retrieving parameters from the case files.
Build configurations.
integer, parameter neko_bcknd_device
integer, parameter, public c_rp
Definition num_types.f90:13
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
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
Generic Dirichlet boundary condition on .
Definition dirichlet.f90:48
Dong outflow condition Follows "A Convective-like Energy-Stable Open Boundary Condition for Simulati...
A struct that contains all info about the time, expand as needed.