Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
fluid_scheme_compressible.f90
Go to the documentation of this file.
1! Copyright (c) 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!
34 use field, only : field_t
37
40 use json_module, only : json_file
41 use num_types, only : rp
42 use mesh, only : mesh_t
44 use space, only : gll
45 use user_intf, only : user_t
47 use mpi_f08
48 use operators, only : cfl_compressible
54 use time_state, only : time_state_t
55 implicit none
56 private
57
59 type, public, abstract, extends(fluid_scheme_base_t) :: fluid_scheme_compressible_t
61 type(field_t), pointer :: m_x => null()
62 type(field_t), pointer :: m_y => null()
63 type(field_t), pointer :: m_z => null()
64 type(field_t), pointer :: e => null()
65 type(field_t), pointer :: max_wave_speed => null()
66 type(field_t), pointer :: s => null()
67
68 real(kind=rp) :: gamma
69
70 type(scratch_registry_t) :: scratch
71
72 contains
74 procedure, pass(this) :: scheme_init => fluid_scheme_compressible_init
76 procedure, pass(this) :: scheme_free => fluid_scheme_compressible_free
77
79 procedure, pass(this) :: validate => fluid_scheme_compressible_validate
81 procedure, pass(this) :: compute_cfl &
84 procedure, pass(this) :: update_material_properties => &
87 procedure, pass(this) :: compute_entropy => &
90 procedure, pass(this) :: compute_max_wave_speed => &
92
94
95contains
103 subroutine fluid_scheme_compressible_init(this, msh, lx, params, scheme, user)
104 class(fluid_scheme_compressible_t), target, intent(inout) :: this
105 type(mesh_t), target, intent(inout) :: msh
106 integer, intent(in) :: lx
107 character(len=*), intent(in) :: scheme
108 type(json_file), target, intent(inout) :: params
109 type(user_t), target, intent(in) :: user
110
111 !
112 ! SEM simulation fundamentals
113 !
114
115 this%msh => msh
116
117 if (msh%gdim .eq. 2) then
118 call this%Xh%init(gll, lx, lx)
119 else
120 call this%Xh%init(gll, lx, lx, lx)
121 end if
122
123 call this%dm_Xh%init(msh, this%Xh)
124
125 call this%gs_Xh%init(this%dm_Xh)
126
127 call this%c_Xh%init(this%gs_Xh)
128
129 ! Local scratch registry
130 call this%scratch%init(this%dm_Xh, 10, 2)
131
132 ! Case parameters
133 this%params => params
134
135 ! Assign a name
136 call json_get_or_default(params, 'case.fluid.name', this%name, "fluid")
137
138 ! Fill mu and rho field with the physical value
139 call neko_field_registry%add_field(this%dm_Xh, this%name // "_mu")
140 call neko_field_registry%add_field(this%dm_Xh, this%name // "_rho")
141 this%mu => neko_field_registry%get_field(this%name // "_mu")
142 this%rho => neko_field_registry%get_field(this%name // "_rho")
143 call field_cfill(this%mu, 0.0_rp, this%mu%size())
144
145 ! Assign momentum fields
146 call neko_field_registry%add_field(this%dm_Xh, "m_x")
147 call neko_field_registry%add_field(this%dm_Xh, "m_y")
148 call neko_field_registry%add_field(this%dm_Xh, "m_z")
149 this%m_x => neko_field_registry%get_field("m_x")
150 this%m_y => neko_field_registry%get_field("m_y")
151 this%m_z => neko_field_registry%get_field("m_z")
152 call this%m_x%init(this%dm_Xh, "m_x")
153 call this%m_y%init(this%dm_Xh, "m_y")
154 call this%m_z%init(this%dm_Xh, "m_z")
155
156 ! Assign energy field
157 call neko_field_registry%add_field(this%dm_Xh, "E")
158 this%E => neko_field_registry%get_field("E")
159 call this%E%init(this%dm_Xh, "E")
160
161 ! Assign maximum wave speed field
162 call neko_field_registry%add_field(this%dm_Xh, "max_wave_speed")
163 this%max_wave_speed => neko_field_registry%get_field("max_wave_speed")
164 call this%max_wave_speed%init(this%dm_Xh, "max_wave_speed")
165
166 ! Assign entropy field
167 call neko_field_registry%add_field(this%dm_Xh, "S")
168 this%S => neko_field_registry%get_field("S")
169 call this%S%init(this%dm_Xh, "S")
170
171 ! ! Assign velocity fields
172 call neko_field_registry%add_field(this%dm_Xh, "u")
173 call neko_field_registry%add_field(this%dm_Xh, "v")
174 call neko_field_registry%add_field(this%dm_Xh, "w")
175 this%u => neko_field_registry%get_field("u")
176 this%v => neko_field_registry%get_field("v")
177 this%w => neko_field_registry%get_field("w")
178 call this%u%init(this%dm_Xh, "u")
179 call this%v%init(this%dm_Xh, "v")
180 call this%w%init(this%dm_Xh, "w")
181 call neko_field_registry%add_field(this%dm_Xh, 'p')
182 this%p => neko_field_registry%get_field('p')
183 call this%p%init(this%dm_Xh, "p")
184
185 !
186 ! Setup right-hand side fields.
187 !
188 allocate(this%f_x)
189 allocate(this%f_y)
190 allocate(this%f_z)
191 call this%f_x%init(this%dm_Xh, fld_name = "fluid_rhs_x")
192 call this%f_y%init(this%dm_Xh, fld_name = "fluid_rhs_y")
193 call this%f_z%init(this%dm_Xh, fld_name = "fluid_rhs_z")
194
195 ! Compressible parameters
196 call json_get_or_default(params, 'case.fluid.gamma', this%gamma, 1.4_rp)
197 end subroutine fluid_scheme_compressible_init
198
202 class(fluid_scheme_compressible_t), intent(inout) :: this
203 call this%dm_Xh%free()
204 call this%gs_Xh%free()
205 call this%c_Xh%free()
206 call this%Xh%free()
207
208 if (associated(this%m_x)) then
209 call this%m_x%free()
210 end if
211
212 if (associated(this%m_y)) then
213 call this%m_y%free()
214 end if
215
216 if (associated(this%m_z)) then
217 call this%m_z%free()
218 end if
219
220 if (associated(this%E)) then
221 call this%E%free()
222 end if
223
224 if (associated(this%max_wave_speed)) then
225 call this%max_wave_speed%free()
226 end if
227
228 if (associated(this%S)) then
229 call this%S%free()
230 end if
231
232 nullify(this%m_x)
233 nullify(this%m_y)
234 nullify(this%m_z)
235 nullify(this%E)
236 nullify(this%max_wave_speed)
237 nullify(this%S)
238
239 nullify(this%u)
240 nullify(this%v)
241 nullify(this%w)
242 nullify(this%p)
243 nullify(this%rho)
244 nullify(this%mu)
245
246 end subroutine fluid_scheme_compressible_free
247
251 class(fluid_scheme_compressible_t), target, intent(inout) :: this
252 integer :: n
253 type(field_t), pointer :: temp
254 integer :: temp_indices(1)
255
256 n = this%dm_Xh%size()
257 call this%scratch%request_field(temp, temp_indices(1))
258
260 call field_col3(this%m_x, this%u, this%rho)
261 call field_col3(this%m_y, this%v, this%rho)
262 call field_col3(this%m_z, this%w, this%rho)
263
266 call field_cmult2(this%E, this%p, 1.0_rp/(this%gamma - 1.0_rp), n)
267 call field_col3(temp, this%u, this%u, n)
268 call field_addcol3(temp, this%v, this%v, n)
269 call field_addcol3(temp, this%w, this%w, n)
270 call field_col2(temp, this%rho, n)
271 call field_cmult(temp, 0.5_rp, n)
272 call field_add2(this%E, temp, n)
273
274 call this%scratch%relinquish_field(temp_indices)
275
277 call this%compute_max_wave_speed()
278
280
285 function fluid_scheme_compressible_compute_cfl(this, dt) result(c)
286 class(fluid_scheme_compressible_t), intent(in) :: this
287 real(kind=rp), intent(in) :: dt
288 real(kind=rp) :: c
289 integer :: n
290
291 associate(u => this%u, v => this%v, w => this%w, p => this%p, &
292 rho => this%rho, xh => this%Xh, c_xh => this%c_Xh, &
293 msh => this%msh, gamma => this%gamma, &
294 max_wave_speed => this%max_wave_speed)
295
296 n = xh%lx * xh%ly * xh%lz * msh%nelv
297
298 ! Use the compressible CFL function with precomputed maximum wave speed
299 c = cfl_compressible(dt, max_wave_speed%x, xh, c_xh, msh%nelv, msh%gdim)
300 end associate
301
303
307 class(fluid_scheme_compressible_t), intent(inout) :: this
308 type(time_state_t), intent(in) :: time
310
314 class(fluid_scheme_compressible_t), intent(inout) :: this
315 integer :: n
316
317 n = this%S%dof%size()
318
320 if (neko_bcknd_device .eq. 1) then
321 call compressible_ops_device_compute_entropy(this%S, this%p, this%rho, this%gamma, n)
322 else
323 call compressible_ops_cpu_compute_entropy(this%S%x, this%p%x, this%rho%x, this%gamma, n)
324 end if
325
327
331 class(fluid_scheme_compressible_t), intent(inout) :: this
332 integer :: n
333
334 n = this%u%dof%size()
335
337 if (neko_bcknd_device .eq. 1) then
338 call compressible_ops_device_compute_max_wave_speed(this%max_wave_speed, &
339 this%u, this%v, this%w, this%gamma, this%p, this%rho, n)
340 else
341 call compressible_ops_cpu_compute_max_wave_speed(this%max_wave_speed%x, &
342 this%u%x, this%v%x, this%w%x, this%gamma, this%p%x, this%rho%x, n)
343 end if
344
Abstract interface to sets rho and mu.
Retrieves a parameter by name or assigns a provided default value. In the latter case also adds the m...
CPU implementation of compressible flow operations.
subroutine, public compressible_ops_cpu_compute_entropy(s, p, rho, gamma, n)
Compute entropy field S = 1/(gamma-1) * rho * (log(p) - gamma * log(rho)) on CPU.
subroutine, public compressible_ops_cpu_compute_max_wave_speed(max_wave_speed, u, v, w, gamma, p, rho, n)
Compute maximum wave speed for compressible flows on CPU.
Device implementation of compressible flow operations.
subroutine, public compressible_ops_device_compute_entropy(s, p, rho, gamma, n)
Compute entropy field S = 1/(gamma-1) * rho * (log(p) - gamma * log(rho)) on device.
subroutine, public compressible_ops_device_compute_max_wave_speed(max_wave_speed, u, v, w, gamma, p, rho, n)
Compute maximum wave speed for compressible flows on device.
subroutine, public field_col2(a, b, n)
Vector multiplication .
subroutine, public field_cmult2(a, b, c, n)
Multiplication by constant c .
subroutine, public field_cfill(a, c, n)
Set all elements to a constant c .
subroutine, public field_addcol3(a, b, c, n)
Returns .
subroutine, public field_add2(a, b, n)
Vector addition .
subroutine, public field_col3(a, b, c, n)
Vector multiplication with 3 vectors .
subroutine, public field_cmult(a, c, n)
Multiplication by constant c .
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
subroutine fluid_scheme_compressible_free(this)
Free allocated memory and cleanup resources.
subroutine fluid_scheme_compressible_update_material_properties(this, time)
Set rho and mu.
subroutine fluid_scheme_compressible_init(this, msh, lx, params, scheme, user)
Initialize common data for compressible fluid scheme.
real(kind=rp) function fluid_scheme_compressible_compute_cfl(this, dt)
Compute CFL number.
subroutine fluid_scheme_compressible_compute_max_wave_speed(this)
Compute maximum wave speed for compressible flows.
subroutine fluid_scheme_compressible_validate(this)
Validate field initialization and compute derived quantities.
subroutine fluid_scheme_compressible_compute_entropy(this)
Compute entropy field S = 1/(gamma-1) * rho * (log(p) - gamma * log(rho))
Utilities for retrieving parameters from the case files.
Defines a mesh.
Definition mesh.f90:34
Build configurations.
integer, parameter neko_bcknd_device
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
Operators.
Definition operators.f90:34
real(kind=rp) function, public cfl_compressible(dt, max_wave_speed, xh, coef, nelv, gdim)
Defines a registry for storing and requesting temporary fields This can be used when you have a funct...
Defines a function space.
Definition space.f90:34
integer, parameter, public gll
Definition space.f90:48
Module with things related to the simulation time.
Interfaces for user interaction with NEKO.
Definition user_intf.f90:34
Base type of all fluid formulations.
A struct that contains all info about the time, expand as needed.
A type collecting all the overridable user routines and flag to suppress type injection from custom m...