Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
blasius.f90
Go to the documentation of this file.
1! Copyright (c) 2021-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 blasius
35 use num_types, only : rp
36 use coefs, only : coef_t
37 use utils, only : nonlinear_index
40 use flow_profile
41 use, intrinsic :: iso_fortran_env
42 use, intrinsic :: iso_c_binding
43 use bc, only : bc_t
44 use json_module, only : json_file
45 use json_utils, only : json_get
46 use time_state, only : time_state_t
47 implicit none
48 private
49
53 type, public, extends(bc_t) :: blasius_t
54 real(kind=rp), dimension(3) :: uinf = [0d0, 0d0, 0d0]
55 real(kind=rp) :: delta
56 procedure(blasius_profile), nopass, pointer :: bla => null()
57 type(c_ptr), private :: blax_d = c_null_ptr
58 type(c_ptr), private :: blay_d = c_null_ptr
59 type(c_ptr), private :: blaz_d = c_null_ptr
60 contains
61 procedure, pass(this) :: apply_scalar => blasius_apply_scalar
62 procedure, pass(this) :: apply_vector => blasius_apply_vector
63 procedure, pass(this) :: apply_scalar_dev => blasius_apply_scalar_dev
64 procedure, pass(this) :: apply_vector_dev => blasius_apply_vector_dev
65 procedure, pass(this) :: set_params => blasius_set_params
67 procedure, pass(this) :: init => blasius_init
69 procedure, pass(this) :: init_from_components => &
72 procedure, pass(this) :: free => blasius_free
74 procedure, pass(this) :: finalize => blasius_finalize
75 end type blasius_t
76
77contains
78
82 subroutine blasius_init(this, coef, json)
83 class(blasius_t), intent(inout), target :: this
84 type(coef_t), target, intent(in) :: coef
85 type(json_file), intent(inout) :: json
86 real(kind=rp) :: delta
87 real(kind=rp), allocatable :: uinf(:)
88 character(len=:), allocatable :: approximation
89
90 call this%init_base(coef)
91
92 call json_get(json, 'delta', delta)
93 call json_get(json, 'approximation', approximation)
94 call json_get(json, 'freestream_velocity', uinf)
95
96 if (size(uinf) .ne. 3) then
97 call neko_error("The uinf keyword for the blasius profile should be an &
98 & array of 3 reals")
99 end if
100
101 call this%init_from_components(coef, delta, uinf, approximation)
102
103 end subroutine blasius_init
104
110 subroutine blasius_init_from_components(this, coef, delta, uinf, &
111 approximation)
112 class(blasius_t), intent(inout), target :: this
113 type(coef_t), target, intent(in) :: coef
114 real(kind=rp) :: delta
115 real(kind=rp) :: uinf(3)
116 character(len=*) :: approximation
117
118 call this%init_base(coef)
119
120 this%delta = delta
121 this%uinf = uinf
122
123 select case (trim(approximation))
124 case ('linear')
125 this%bla => blasius_linear
126 case ('quadratic')
127 this%bla => blasius_quadratic
128 case ('cubic')
129 this%bla => blasius_cubic
130 case ('quartic')
131 this%bla => blasius_quartic
132 case ('sin')
133 this%bla => blasius_sin
134 case ('tanh')
135 this%bla => blasius_tanh
136 case default
137 call neko_error('Invalid Blasius approximation')
138 end select
139 end subroutine blasius_init_from_components
140
141 subroutine blasius_free(this)
142 class(blasius_t), target, intent(inout) :: this
143
144 call this%free_base()
145 nullify(this%bla)
146
147 if (c_associated(this%blax_d)) then
148 call device_free(this%blax_d)
149 end if
150
151 if (c_associated(this%blay_d)) then
152 call device_free(this%blay_d)
153 end if
154
155 if (c_associated(this%blaz_d)) then
156 call device_free(this%blaz_d)
157 end if
158
159 end subroutine blasius_free
160
162 subroutine blasius_apply_scalar(this, x, n, time, strong)
163 class(blasius_t), intent(inout) :: this
164 integer, intent(in) :: n
165 real(kind=rp), intent(inout), dimension(n) :: x
166 type(time_state_t), intent(in), optional :: time
167 logical, intent(in), optional :: strong
168 end subroutine blasius_apply_scalar
169
171 subroutine blasius_apply_scalar_dev(this, x_d, time, strong, strm)
172 class(blasius_t), intent(inout), target :: this
173 type(c_ptr), intent(inout) :: x_d
174 type(time_state_t), intent(in), optional :: time
175 logical, intent(in), optional :: strong
176 type(c_ptr), intent(inout) :: strm
177 end subroutine blasius_apply_scalar_dev
178
180 subroutine blasius_apply_vector(this, x, y, z, n, time, strong)
181 class(blasius_t), intent(inout) :: this
182 integer, intent(in) :: n
183 real(kind=rp), intent(inout), dimension(n) :: x
184 real(kind=rp), intent(inout), dimension(n) :: y
185 real(kind=rp), intent(inout), dimension(n) :: z
186 type(time_state_t), intent(in), optional :: time
187 logical, intent(in), optional :: strong
188 integer :: i, m, k, idx(4), facet
189 logical :: strong_
190
191 if (present(strong)) then
192 strong_ = strong
193 else
194 strong_ = .true.
195 end if
196
197 associate(xc => this%coef%dof%x, yc => this%coef%dof%y, &
198 zc => this%coef%dof%z, nx => this%coef%nx, ny => this%coef%ny, &
199 nz => this%coef%nz, lx => this%coef%Xh%lx)
200 m = this%msk(0)
201 if (strong_) then
202 do i = 1, m
203 k = this%msk(i)
204 facet = this%facet(i)
205 idx = nonlinear_index(k, lx, lx, lx)
206 select case (facet)
207 case (1, 2)
208 x(k) = this%bla(zc(idx(1), idx(2), idx(3), idx(4)), &
209 this%delta, this%uinf(1))
210 y(k) = 0.0_rp
211 z(k) = 0.0_rp
212 case (3, 4)
213 x(k) = 0.0_rp
214 y(k) = this%bla(xc(idx(1), idx(2), idx(3), idx(4)), &
215 this%delta, this%uinf(2))
216 z(k) = 0.0_rp
217 case (5, 6)
218 x(k) = 0.0_rp
219 y(k) = 0.0_rp
220 z(k) = this%bla(yc(idx(1), idx(2), idx(3), idx(4)), &
221 this%delta, this%uinf(3))
222 end select
223 end do
224 end if
225 end associate
226 end subroutine blasius_apply_vector
227
229 subroutine blasius_apply_vector_dev(this, x_d, y_d, z_d, time, strong, strm)
230 class(blasius_t), intent(inout), target :: this
231 type(c_ptr), intent(inout) :: x_d
232 type(c_ptr), intent(inout) :: y_d
233 type(c_ptr), intent(inout) :: z_d
234 type(time_state_t), intent(in), optional :: time
235 logical, intent(in), optional :: strong
236 integer :: i, m, k, idx(4), facet
237 integer(c_size_t) :: s
238 real(kind=rp), allocatable :: bla_x(:), bla_y(:), bla_z(:)
239 logical :: strong_
240 type(c_ptr), intent(inout) :: strm
241
242 if (present(strong)) then
243 strong_ = strong
244 else
245 strong_ = .true.
246 end if
247
248 associate(xc => this%coef%dof%x, yc => this%coef%dof%y, &
249 zc => this%coef%dof%z, nx => this%coef%nx, ny => this%coef%ny, &
250 nz => this%coef%nz, lx => this%coef%Xh%lx , &
251 blax_d => this%blax_d, blay_d => this%blay_d, &
252 blaz_d => this%blaz_d)
253
254 m = this%msk(0)
255
256
257 ! Pretabulate values during first call to apply
258 if (.not. c_associated(blax_d) .and. strong_ .and. this%msk(0) .gt. 0) then
259 allocate(bla_x(m), bla_y(m), bla_z(m)) ! Temp arrays
260
261 if (rp .eq. real32) then
262 s = m * 4
263 else if (rp .eq. real64) then
264 s = m * 8
265 end if
266
267 call device_alloc(blax_d, s)
268 call device_alloc(blay_d, s)
269 call device_alloc(blaz_d, s)
270
271 do i = 1, m
272 k = this%msk(i)
273 facet = this%facet(i)
274 idx = nonlinear_index(k, lx, lx, lx)
275 select case (facet)
276 case (1,2)
277 bla_x(i) = this%bla(zc(idx(1), idx(2), idx(3), idx(4)), &
278 this%delta, this%uinf(1))
279 bla_y(i) = 0.0_rp
280 bla_z(i) = 0.0_rp
281 case (3,4)
282 bla_x(i) = 0.0_rp
283 bla_y(i) = this%bla(xc(idx(1), idx(2), idx(3), idx(4)), &
284 this%delta, this%uinf(2))
285 bla_z(i) = 0.0_rp
286 case (5,6)
287 bla_x(i) = 0.0_rp
288 bla_y(i) = 0.0_rp
289 bla_z(i) = this%bla(yc(idx(1), idx(2), idx(3), idx(4)), &
290 this%delta, this%uinf(3))
291 end select
292 end do
293
294 call device_memcpy(bla_x, blax_d, m, host_to_device, sync = .false.)
295 call device_memcpy(bla_y, blay_d, m, host_to_device, sync = .false.)
296 call device_memcpy(bla_z, blaz_d, m, host_to_device, sync = .true.)
297
298 deallocate(bla_x, bla_y, bla_z)
299 end if
300
301 if (strong_ .and. this%msk(0) .gt. 0) then
302 call device_inhom_dirichlet_apply_vector(this%msk_d, x_d, y_d, z_d, &
303 blax_d, blay_d, blaz_d, m, strm)
304 end if
305
306 end associate
307
308 end subroutine blasius_apply_vector_dev
309
311 subroutine blasius_set_params(this, uinf, delta, type)
312 class(blasius_t), intent(inout) :: this
313 real(kind=rp), intent(in) :: uinf(3)
314 real(kind=rp), intent(in) :: delta
315 character(len=*) :: type
316 this%delta = delta
317 this%uinf = uinf
318
319 select case (trim(type))
320 case ('linear')
321 this%bla => blasius_linear
322 case ('quadratic')
323 this%bla => blasius_quadratic
324 case ('cubic')
325 this%bla => blasius_cubic
326 case ('quartic')
327 this%bla => blasius_quartic
328 case ('sin')
329 this%bla => blasius_sin
330 case ('tanh')
331 this%bla => blasius_tanh
332 case default
333 call neko_error('Invalid Blasius approximation')
334 end select
335 end subroutine blasius_set_params
336
338 subroutine blasius_finalize(this, only_facets)
339 class(blasius_t), target, intent(inout) :: this
340 logical, optional, intent(in) :: only_facets
341 logical :: only_facets_
342
343 if (present(only_facets)) then
344 only_facets_ = only_facets
345 else
346 only_facets_ = .false.
347 end if
348
349 call this%finalize_base(only_facets_)
350 end subroutine blasius_finalize
351end module blasius
__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
Abstract interface for computing a Blasius flow profile.
Retrieves a parameter by name or throws an error.
Defines a boundary condition.
Definition bc.f90:34
Defines a Blasius profile dirichlet condition.
Definition blasius.f90:34
subroutine blasius_finalize(this, only_facets)
Finalize.
Definition blasius.f90:339
subroutine blasius_apply_scalar_dev(this, x_d, time, strong, strm)
No-op scalar apply (device version)
Definition blasius.f90:172
subroutine blasius_init(this, coef, json)
Constructor.
Definition blasius.f90:83
subroutine blasius_init_from_components(this, coef, delta, uinf, approximation)
Constructor from components.
Definition blasius.f90:112
subroutine blasius_apply_vector(this, x, y, z, n, time, strong)
Apply blasius conditions (vector valued)
Definition blasius.f90:181
subroutine blasius_apply_scalar(this, x, n, time, strong)
No-op scalar apply.
Definition blasius.f90:163
subroutine blasius_apply_vector_dev(this, x_d, y_d, z_d, time, strong, strm)
Apply blasius conditions (vector valued) (device version)
Definition blasius.f90:230
subroutine blasius_free(this)
Definition blasius.f90:142
subroutine blasius_set_params(this, uinf, delta, type)
Set Blasius parameters.
Definition blasius.f90:312
Coefficients.
Definition coef.f90:34
Device abstraction, common interface for various accelerators.
Definition device.F90:34
integer, parameter, public host_to_device
Definition device.F90:47
subroutine, public device_free(x_d)
Deallocate memory on the device.
Definition device.F90:214
subroutine, public device_alloc(x_d, s)
Allocate memory on the device.
Definition device.F90:187
Defines a flow profile.
real(kind=rp) function, public blasius_quadratic(y, delta, u)
Quadratic approximate Blasius Profile .
real(kind=rp) function, public blasius_quartic(y, delta, u)
Quartic approximate Blasius Profile .
real(kind=rp) function, public blasius_sin(y, delta, u)
Sinusoidal approximate Blasius Profile .
real(kind=rp) function, public blasius_cubic(y, delta, u)
Cubic approximate Blasius Profile .
real(kind=rp) function, public blasius_tanh(y, delta, u)
Hyperbolic tangent approximate Blasius Profile from O. Savas (2012) where is the 99 percent thickne...
real(kind=rp) function, public blasius_linear(y, delta, u)
Linear approximate Blasius profile .
Utilities for retrieving parameters from the case files.
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
Blasius profile for inlet (vector valued).
Definition blasius.f90:53
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.