40 use json_module,
only : json_file
65 type(
field_t),
pointer :: max_wave_speed => null()
68 real(kind=
rp) :: gamma
81 procedure, pass(this) :: compute_cfl &
87 procedure, pass(this) :: compute_entropy => &
90 procedure, pass(this) :: compute_max_wave_speed => &
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
117 if (msh%gdim .eq. 2)
then
118 call this%Xh%init(
gll, lx, lx)
120 call this%Xh%init(
gll, lx, lx, lx)
123 call this%dm_Xh%init(msh, this%Xh)
125 call this%gs_Xh%init(this%dm_Xh)
127 call this%c_Xh%init(this%gs_Xh)
130 call this%scratch%init(this%dm_Xh, 10, 2)
133 this%params => params
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")
159 call this%E%init(this%dm_Xh,
"E")
164 call this%max_wave_speed%init(this%dm_Xh,
"max_wave_speed")
169 call this%S%init(this%dm_Xh,
"S")
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")
183 call this%p%init(this%dm_Xh,
"p")
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")
203 call this%dm_Xh%free()
204 call this%gs_Xh%free()
205 call this%c_Xh%free()
208 if (
associated(this%m_x))
then
212 if (
associated(this%m_y))
then
216 if (
associated(this%m_z))
then
220 if (
associated(this%E))
then
224 if (
associated(this%max_wave_speed))
then
225 call this%max_wave_speed%free()
228 if (
associated(this%S))
then
236 nullify(this%max_wave_speed)
254 integer :: temp_indices(1)
256 n = this%dm_Xh%size()
257 call this%scratch%request_field(temp, temp_indices(1))
266 call field_cmult2(this%E, this%p, 1.0_rp/(this%gamma - 1.0_rp), n)
274 call this%scratch%relinquish_field(temp_indices)
277 call this%compute_max_wave_speed()
287 real(kind=
rp),
intent(in) :: dt
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)
296 n = xh%lx * xh%ly * xh%lz * msh%nelv
308 type(time_state_t),
intent(in) :: time
317 n = this%S%dof%size()
320 if (neko_bcknd_device .eq. 1)
then
321 call compressible_ops_device_compute_entropy(this%S, this%p, this%rho, this%gamma, n)
323 call compressible_ops_cpu_compute_entropy(this%S%x, this%p%x, this%rho%x, this%gamma, n)
334 n = this%u%dof%size()
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)
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)
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.
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.
integer, parameter neko_bcknd_device
integer, parameter, public rp
Global precision used in computations.
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.
integer, parameter, public gll
Module with things related to the simulation time.
Interfaces for user interaction with NEKO.
Base type of all fluid formulations.
Base type of compressible 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...