Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
inflow.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!
34module inflow
36 use num_types, only : rp
37 use bc, only : bc_t
38 use, intrinsic :: iso_c_binding, only : c_ptr, c_loc
39 use coefs, only : coef_t
40 use json_module, only : json_file
41 use json_utils, only : json_get
42 use time_state, only : time_state_t
43 implicit none
44 private
45
47 type, public, extends(bc_t) :: inflow_t
48 real(kind=rp), dimension(3) :: x = [0d0, 0d0, 0d0]
49 contains
50 procedure, pass(this) :: apply_scalar => inflow_apply_scalar
51 procedure, pass(this) :: apply_vector => inflow_apply_vector
52 procedure, pass(this) :: apply_scalar_dev => inflow_apply_scalar_dev
53 procedure, pass(this) :: apply_vector_dev => inflow_apply_vector_dev
55 procedure, pass(this) :: init => inflow_init
57 procedure, pass(this) :: free => inflow_free
59 procedure, pass(this) :: finalize => inflow_finalize
60 end type inflow_t
61
62contains
63
67 subroutine inflow_init(this, coef, json)
68 class(inflow_t), intent(inout), target :: this
69 type(coef_t), target, intent(in) :: coef
70 type(json_file), intent(inout) ::json
71 real(kind=rp), allocatable :: x(:)
72
73 call this%init_base(coef)
74 call json_get(json, 'value', x)
75 this%x = x
76 end subroutine inflow_init
77
79 subroutine inflow_apply_scalar(this, x, n, time, strong)
80 class(inflow_t), intent(inout) :: this
81 integer, intent(in) :: n
82 real(kind=rp), intent(inout), dimension(n) :: x
83 type(time_state_t), intent(in), optional :: time
84 logical, intent(in), optional :: strong
85 end subroutine inflow_apply_scalar
86
88 subroutine inflow_apply_scalar_dev(this, x_d, time, strong, strm)
89 class(inflow_t), intent(inout), target :: this
90 type(c_ptr), intent(inout) :: x_d
91 type(time_state_t), intent(in), optional :: time
92 logical, intent(in), optional :: strong
93 type(c_ptr), intent(inout) :: strm
94 end subroutine inflow_apply_scalar_dev
95
97 subroutine inflow_apply_vector(this, x, y, z, n, time, strong)
98 class(inflow_t), intent(inout) :: this
99 integer, intent(in) :: n
100 real(kind=rp), intent(inout), dimension(n) :: x
101 real(kind=rp), intent(inout), dimension(n) :: y
102 real(kind=rp), intent(inout), dimension(n) :: z
103 type(time_state_t), intent(in), optional :: time
104 logical, intent(in), optional :: strong
105 integer :: i, m, k
106 logical :: strong_
107
108 if (present(strong)) then
109 strong_ = strong
110 else
111 strong_ = .true.
112 end if
113
114 m = this%msk(0)
115
116 if (strong_) then
117 do i = 1, m
118 k = this%msk(i)
119 x(k) = this%x(1)
120 y(k) = this%x(2)
121 z(k) = this%x(3)
122 end do
123 end if
124 end subroutine inflow_apply_vector
125
127 subroutine inflow_apply_vector_dev(this, x_d, y_d, z_d, &
128 time, strong, strm)
129 class(inflow_t), intent(inout), target :: this
130 type(c_ptr), intent(inout) :: x_d
131 type(c_ptr), intent(inout) :: y_d
132 type(c_ptr), intent(inout) :: z_d
133 type(time_state_t), intent(in), optional :: time
134 logical, intent(in), optional :: strong
135 type(c_ptr), intent(inout) :: strm
136 logical :: strong_
137
138 if (present(strong)) then
139 strong_ = strong
140 else
141 strong_ = .true.
142 end if
143
144 if (strong_ .and. (this%msk(0) .gt. 0)) then
145 call device_inflow_apply_vector(this%msk_d, x_d, y_d, z_d, &
146 c_loc(this%x), this%msk(0), strm)
147 end if
148
149 end subroutine inflow_apply_vector_dev
150
152 subroutine inflow_free(this)
153 class(inflow_t), target, intent(inout) :: this
154
155 call this%free_base()
156 end subroutine inflow_free
157
159 subroutine inflow_finalize(this, only_facets)
160 class(inflow_t), target, intent(inout) :: this
161 logical, optional, intent(in) :: only_facets
162 logical :: only_facets_
163
164 integer :: i
165
166 if (present(only_facets)) then
167 only_facets_ = only_facets
168 else
169 only_facets_ = .false.
170 end if
171
172 call this%finalize_base(only_facets_)
173 end subroutine inflow_finalize
174
175
176end module inflow
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_inflow_apply_vector(msk, x, y, z, g, m, strm)
Defines inflow dirichlet conditions.
Definition inflow.f90:34
subroutine inflow_apply_scalar_dev(this, x_d, time, strong, strm)
No-op scalar apply (device version)
Definition inflow.f90:89
subroutine inflow_apply_scalar(this, x, n, time, strong)
No-op scalar apply.
Definition inflow.f90:80
subroutine inflow_finalize(this, only_facets)
Finalize.
Definition inflow.f90:160
subroutine inflow_init(this, coef, json)
Constructor.
Definition inflow.f90:68
subroutine inflow_apply_vector(this, x, y, z, n, time, strong)
Apply inflow conditions (vector valued)
Definition inflow.f90:98
subroutine inflow_apply_vector_dev(this, x_d, y_d, z_d, time, strong, strm)
Apply inflow conditions (vector valued) (device version)
Definition inflow.f90:129
subroutine inflow_free(this)
Destructor.
Definition inflow.f90:153
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.
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
Dirichlet condition for inlet (vector valued)
Definition inflow.f90:47
A struct that contains all info about the time, expand as needed.