Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
dirichlet.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!
36 use num_types, only : rp
37 use bc, only : bc_t
38 use coefs, only : coef_t
39 use json_module, only : json_file
40 use json_utils, only : json_get
41 use, intrinsic :: iso_c_binding, only : c_ptr
42 use time_state, only : time_state_t
43 implicit none
44 private
45
48 type, public, extends(bc_t) :: dirichlet_t
49 real(kind=rp), private :: g
50 contains
51 procedure, pass(this) :: apply_scalar => dirichlet_apply_scalar
52 procedure, pass(this) :: apply_vector => dirichlet_apply_vector
53 procedure, pass(this) :: apply_scalar_dev => dirichlet_apply_scalar_dev
54 procedure, pass(this) :: apply_vector_dev => dirichlet_apply_vector_dev
55 procedure, pass(this) :: set_g => dirichlet_set_g
57 procedure, pass(this) :: init => dirichlet_init
59 procedure, pass(this) :: init_from_components => &
62 procedure, pass(this) :: free => dirichlet_free
64 procedure, pass(this) :: finalize => dirichlet_finalize
65 end type dirichlet_t
66
67contains
68
72 subroutine dirichlet_init(this, coef, json)
73 class(dirichlet_t), intent(inout), target :: this
74 type(coef_t), target, intent(in) :: coef
75 type(json_file), intent(inout) ::json
76 real(kind=rp) :: g
77
78 call this%init_base(coef)
79 call json_get(json , "value", g)
80
81 this%g = g
82 end subroutine dirichlet_init
83
87 subroutine dirichlet_init_from_components(this, coef, g)
88 class(dirichlet_t), intent(inout), target :: this
89 type(coef_t), target, intent(in) :: coef
90 real(kind=rp), intent(in) :: g
91
92 call this%init_base(coef)
93 this%g = g
95
98 subroutine dirichlet_apply_scalar(this, x, n, time, strong)
99 class(dirichlet_t), intent(inout) :: this
100 integer, intent(in) :: n
101 real(kind=rp), intent(inout), dimension(n) :: x
102 type(time_state_t), intent(in), optional :: time
103 logical, intent(in), optional :: strong
104 integer :: i, m, k
105 logical :: strong_
106
107 if (present(strong)) then
108 strong_ = strong
109 else
110 strong_ = .true.
111 end if
112
113 if (strong_) then
114 m = this%msk(0)
115 do i = 1, m
116 k = this%msk(i)
117 x(k) = this%g
118 end do
119 end if
120 end subroutine dirichlet_apply_scalar
121
124 subroutine dirichlet_apply_vector(this, x, y, z, n, time, strong)
125 class(dirichlet_t), intent(inout) :: this
126 integer, intent(in) :: n
127 real(kind=rp), intent(inout), dimension(n) :: x
128 real(kind=rp), intent(inout), dimension(n) :: y
129 real(kind=rp), intent(inout), dimension(n) :: z
130 type(time_state_t), intent(in), optional :: time
131 logical, intent(in), optional :: strong
132 integer :: i, m, k
133 logical :: strong_
134
135 if (present(strong)) then
136 strong_ = strong
137 else
138 strong_ = .true.
139 end if
140
141 if (strong_) then
142 m = this%msk(0)
143 do i = 1, m
144 k = this%msk(i)
145 x(k) = this%g
146 y(k) = this%g
147 z(k) = this%g
148 end do
149 end if
150
151 end subroutine dirichlet_apply_vector
152
155 subroutine dirichlet_apply_scalar_dev(this, x_d, time, strong, strm)
156 class(dirichlet_t), intent(inout), target :: this
157 type(c_ptr), intent(inout) :: x_d
158 type(time_state_t), intent(in), optional :: time
159 logical, intent(in), optional :: strong
160 type(c_ptr), intent(inout) :: strm
161 logical :: strong_
162
163 if (present(strong)) then
164 strong_ = strong
165 else
166 strong_ = .true.
167 end if
168
169 if (strong_ .and. this%msk(0) .gt. 0) then
170 call device_dirichlet_apply_scalar(this%msk_d, x_d, &
171 this%g, size(this%msk), strm)
172 end if
173
174 end subroutine dirichlet_apply_scalar_dev
175
178 subroutine dirichlet_apply_vector_dev(this, x_d, y_d, z_d, &
179 time, strong, strm)
180 class(dirichlet_t), intent(inout), target :: this
181 type(c_ptr), intent(inout) :: x_d
182 type(c_ptr), intent(inout) :: y_d
183 type(c_ptr), intent(inout) :: z_d
184 type(time_state_t), intent(in), optional :: time
185 logical, intent(in), optional :: strong
186 type(c_ptr), intent(inout) :: strm
187 logical :: strong_
188
189 if (present(strong)) then
190 strong_ = strong
191 else
192 strong_ = .true.
193 end if
194
195 if (strong_ .and. this%msk(0) .gt. 0) then
196 call device_dirichlet_apply_vector(this%msk_d, x_d, y_d, z_d, this%g, &
197 size(this%msk), strm)
198 end if
199
200 end subroutine dirichlet_apply_vector_dev
201
203 subroutine dirichlet_set_g(this, g)
204 class(dirichlet_t), intent(inout) :: this
205 real(kind=rp), intent(in) :: g
206
207 this%g = g
208
209 end subroutine dirichlet_set_g
210
212 subroutine dirichlet_free(this)
213 class(dirichlet_t), target, intent(inout) :: this
214
215 call this%free_base
216
217 end subroutine dirichlet_free
218
220 subroutine dirichlet_finalize(this, only_facets)
221 class(dirichlet_t), target, intent(inout) :: this
222 logical, optional, intent(in) :: only_facets
223 logical :: only_facets_
224
225 if (present(only_facets)) then
226 only_facets_ = only_facets
227 else
228 only_facets_ = .false.
229 end if
230
231 call this%finalize_base(only_facets_)
232 end subroutine dirichlet_finalize
233
234end module dirichlet
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_dirichlet_apply_scalar(msk, x, g, m, strm)
subroutine, public device_dirichlet_apply_vector(msk, x, y, z, g, m, strm)
Defines a dirichlet boundary condition.
Definition dirichlet.f90:34
subroutine dirichlet_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 dirichlet_finalize(this, only_facets)
Finalize.
subroutine dirichlet_apply_scalar_dev(this, x_d, time, strong, strm)
Boundary condition apply for a generic Dirichlet condition to a vector x (device version)
subroutine dirichlet_apply_scalar(this, x, n, time, strong)
Boundary condition apply for a generic Dirichlet condition to a vector x.
Definition dirichlet.f90:99
subroutine dirichlet_apply_vector(this, x, y, z, n, time, strong)
Boundary condition apply for a generic Dirichlet condition to vectors x, y and z.
subroutine dirichlet_free(this)
Destructor.
subroutine dirichlet_init_from_components(this, coef, g)
Constructor from components.
Definition dirichlet.f90:88
subroutine dirichlet_init(this, coef, json)
Constructor from JSON.
Definition dirichlet.f90:73
subroutine dirichlet_set_g(this, g)
Set value of .
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
Generic Dirichlet boundary condition on .
Definition dirichlet.f90:48
A struct that contains all info about the time, expand as needed.