40 use json_module,
only : json_file
68 type(
field_t),
pointer :: max_wave_speed => null()
71 real(kind=
rp) :: gamma
74 integer(kind=i8) :: glb_n_points
76 integer(kind=i8) :: glb_unique_points
89 procedure, pass(this) :: compute_cfl &
95 procedure, pass(this) :: compute_entropy => &
98 procedure, pass(this) :: compute_max_wave_speed => &
101 procedure, pass(this) :: log_solver_info => &
116 type(
mesh_t),
target,
intent(inout) :: msh
117 integer,
intent(in) :: lx
118 character(len=*),
intent(in) :: scheme
119 type(json_file),
target,
intent(inout) :: params
120 type(
user_t),
target,
intent(in) :: user
128 if (msh%gdim .eq. 2)
then
129 call this%Xh%init(
gll, lx, lx)
131 call this%Xh%init(
gll, lx, lx, lx)
134 call this%dm_Xh%init(msh, this%Xh)
136 call this%gs_Xh%init(this%dm_Xh)
138 call this%c_Xh%init(this%gs_Xh)
141 call this%scratch%init(this%dm_Xh, 10, 2)
144 this%params => params
163 call this%m_x%init(this%dm_Xh,
"m_x")
164 call this%m_y%init(this%dm_Xh,
"m_y")
165 call this%m_z%init(this%dm_Xh,
"m_z")
170 call this%E%init(this%dm_Xh,
"E")
175 call this%max_wave_speed%init(this%dm_Xh,
"max_wave_speed")
180 call this%S%init(this%dm_Xh,
"S")
189 call this%u%init(this%dm_Xh,
"u")
190 call this%v%init(this%dm_Xh,
"v")
191 call this%w%init(this%dm_Xh,
"w")
194 call this%p%init(this%dm_Xh,
"p")
202 call this%f_x%init(this%dm_Xh, fld_name =
"fluid_rhs_x")
203 call this%f_y%init(this%dm_Xh, fld_name =
"fluid_rhs_y")
204 call this%f_z%init(this%dm_Xh, fld_name =
"fluid_rhs_z")
210 this%glb_n_points = int(this%msh%glb_nelv,
i8)*int(this%Xh%lxyz,
i8)
211 this%glb_unique_points = int(
glsum(this%c_Xh%mult, this%dm_Xh%size()),
i8)
216 call this%log_solver_info(params, scheme, lx)
223 call this%dm_Xh%free()
224 call this%gs_Xh%free()
225 call this%c_Xh%free()
228 if (
associated(this%m_x))
then
232 if (
associated(this%m_y))
then
236 if (
associated(this%m_z))
then
240 if (
associated(this%E))
then
244 if (
associated(this%max_wave_speed))
then
245 call this%max_wave_speed%free()
248 if (
associated(this%S))
then
256 nullify(this%max_wave_speed)
274 integer :: temp_indices(1)
276 n = this%dm_Xh%size()
277 call this%scratch%request_field(temp, temp_indices(1))
286 call field_cmult2(this%E, this%p, 1.0_rp/(this%gamma - 1.0_rp), n)
294 call this%scratch%relinquish_field(temp_indices)
297 call this%compute_max_wave_speed()
307 real(kind=
rp),
intent(in) :: dt
311 associate(u => this%u, v => this%v, w => this%w, p => this%p, &
312 rho => this%rho, xh => this%Xh, c_xh => this%c_Xh, &
313 msh => this%msh, gamma => this%gamma, &
314 max_wave_speed => this%max_wave_speed)
316 n = xh%lx * xh%ly * xh%lz * msh%nelv
328 type(time_state_t),
intent(in) :: time
337 n = this%S%dof%size()
340 if (neko_bcknd_device .eq. 1)
then
341 call compressible_ops_device_compute_entropy(this%S, this%p, this%rho, this%gamma, n)
343 call compressible_ops_cpu_compute_entropy(this%S%x, this%p%x, this%rho%x, this%gamma, n)
354 n = this%u%dof%size()
357 if (neko_bcknd_device .eq. 1)
then
358 call compressible_ops_device_compute_max_wave_speed(this%max_wave_speed, &
359 this%u, this%v, this%w, this%gamma, this%p, this%rho, n)
361 call compressible_ops_cpu_compute_max_wave_speed(this%max_wave_speed%x, &
362 this%u%x, this%v%x, this%w%x, this%gamma, this%p%x, this%rho%x, n)
374 type(json_file),
intent(inout) :: params
375 character(len=*),
intent(in) :: scheme
376 integer,
intent(in) :: lx
377 character(len=LOG_SIZE) :: log_buf
378 logical :: logical_val
379 real(kind=rp) :: real_val
380 integer :: integer_val
382 call neko_log%section(
'Fluid')
383 write(log_buf,
'(A, A)')
'Type : ', trim(scheme)
384 call neko_log%message(log_buf)
385 write(log_buf,
'(A, A)')
'Name : ', trim(this%name)
386 call neko_log%message(log_buf)
390 write(log_buf,
'(A, I1)')
'Poly order : ', lx-1
391 else if (lx .ge. 10)
then
392 write(log_buf,
'(A, I2)')
'Poly order : ', lx-1
394 write(log_buf,
'(A, I3)')
'Poly order : ', lx-1
396 call neko_log%message(log_buf)
399 write(log_buf,
'(A, I0)')
'GLL points : ', this%glb_n_points
400 call neko_log%message(log_buf)
401 write(log_buf,
'(A, I0)')
'Unique pts.: ', this%glb_unique_points
402 call neko_log%message(log_buf)
405 write(log_buf,
'(A,ES13.6)')
'gamma :', this%gamma
406 call neko_log%message(log_buf)
409 call json_get_or_default(params,
'case.numerics.c_avisc_low', real_val, 0.5_rp)
410 write(log_buf,
'(A,ES13.6)')
'c_avisc_low:', real_val
411 call neko_log%message(log_buf)
413 call json_get_or_default(params,
'case.numerics.time_order', integer_val, 4)
414 write(log_buf,
'(A, I0)')
'RK order : ', integer_val
415 call neko_log%message(log_buf)
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_log_solver_info(this, params, scheme, lx)
Log comprehensive solver information.
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.
type(log_t), public neko_log
Global log stream.
integer, parameter, public log_size
real(kind=rp) function, public glsum(a, n)
Sum a vector of length n.
integer, parameter neko_bcknd_device
integer, parameter, public i8
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...