52 use precon,
only :
pc_t, precon_factory, precon_destroy
62 use json_module,
only : json_file
78 class(
ksp_t),
allocatable :: ksp_vel
79 class(
ksp_t),
allocatable :: ksp_prs
80 class(
pc_t),
allocatable :: pc_vel
81 class(
pc_t),
allocatable :: pc_prs
82 integer :: vel_projection_dim
83 integer :: pr_projection_dim
84 integer :: vel_projection_activ_step
85 integer :: pr_projection_activ_step
86 logical :: strict_convergence
93 logical :: forced_flow_rate = .false.
96 character(len=:),
allocatable :: nut_field_name
99 type(
field_t),
pointer :: mu_tot => null()
102 integer(kind=i8) :: glb_n_points
104 integer(kind=i8) :: glb_unique_points
119 procedure, pass(this) :: set_material_properties => &
133 module subroutine fluid_scheme_factory(object, type_name)
135 character(len=*) :: type_name
136 end subroutine fluid_scheme_factory
139 public :: fluid_scheme_incompressible_t, fluid_scheme_factory
144 subroutine fluid_scheme_init_base(this, msh, lx, params, scheme, user, &
147 class(fluid_scheme_incompressible_t),
target,
intent(inout) :: this
148 type(
mesh_t),
target,
intent(inout) :: msh
149 integer,
intent(in) :: lx
150 character(len=*),
intent(in) :: scheme
151 type(json_file),
target,
intent(inout) :: params
152 type(
user_t),
target,
intent(in) :: user
153 logical,
intent(in) :: kspv_init
155 character(len=LOG_SIZE) :: log_buf
156 real(kind=
rp),
allocatable :: real_vec(:)
157 real(kind=
rp) :: real_val, kappa, b, z0
158 logical :: logical_val
159 integer :: integer_val, ierr
160 type(json_file) :: wm_json
161 character(len=:),
allocatable :: string_val1, string_val2
162 real(kind=
rp) :: gjp_param_a, gjp_param_b
163 type(json_file) :: json_subdict
171 if (msh%gdim .eq. 2)
then
172 call this%Xh%init(
gll, lx, lx)
174 call this%Xh%init(
gll, lx, lx, lx)
177 call this%dm_Xh%init(msh, this%Xh)
179 call this%gs_Xh%init(this%dm_Xh)
181 call this%c_Xh%init(this%gs_Xh)
184 call this%scratch%init(this%dm_Xh, 10, 2)
194 write(log_buf,
'(A, A)')
'Type : ', trim(scheme)
196 write(log_buf,
'(A, A)')
'Name : ', trim(this%name)
210 call this%set_material_properties(params,
user)
214 'case.fluid.velocity_solver.projection_space_size', &
215 this%vel_projection_dim, 0)
217 'case.fluid.pressure_solver.projection_space_size', &
218 this%pr_projection_dim, 0)
220 'case.fluid.velocity_solver.projection_hold_steps', &
221 this%vel_projection_activ_step, 5)
223 'case.fluid.pressure_solver.projection_hold_steps', &
224 this%pr_projection_activ_step, 5)
229 if (params%valid_path(
"case.fluid.flow_rate_force"))
then
230 this%forced_flow_rate = .true.
235 write(log_buf,
'(A, I1)')
'Poly order : ', lx-1
236 else if (lx .ge. 10)
then
237 write(log_buf,
'(A, I2)')
'Poly order : ', lx-1
239 write(log_buf,
'(A, I3)')
'Poly order : ', lx-1
242 this%glb_n_points = int(this%msh%glb_nelv,
i8)*int(this%Xh%lxyz,
i8)
243 this%glb_unique_points = int(
glsum(this%c_Xh%mult, this%dm_Xh%size()),
i8)
245 write(log_buf,
'(A, I0)')
'GLL points : ', this%glb_n_points
247 write(log_buf,
'(A, I0)')
'Unique pts.: ', this%glb_unique_points
251 call json_get(params,
'case.numerics.dealias', logical_val)
252 write(log_buf,
'(A, L1)')
'Dealias : ', logical_val
258 write(log_buf,
'(A, L1)')
'Save bdry : ', logical_val
262 logical_val, .false.)
263 write(log_buf,
'(A, L1)')
'Full stress: ', logical_val
273 call this%f_x%init(this%dm_Xh, fld_name =
"fluid_rhs_x")
274 call this%f_y%init(this%dm_Xh, fld_name =
"fluid_rhs_y")
275 call this%f_z%init(this%dm_Xh, fld_name =
"fluid_rhs_z")
279 call neko_log%section(
"Velocity solver")
281 'case.fluid.velocity_solver.max_iterations', &
283 call json_get(params,
'case.fluid.velocity_solver.type', string_val1)
284 call json_get(params,
'case.fluid.velocity_solver.preconditioner.type', &
287 'case.fluid.velocity_solver.preconditioner', json_subdict)
288 call json_get(params,
'case.fluid.velocity_solver.absolute_tolerance', &
291 'case.fluid.velocity_solver.monitor', &
292 logical_val, .false.)
294 call neko_log%message(
'Type : ('// trim(string_val1) // &
295 ', ' // trim(string_val2) //
')')
297 write(log_buf,
'(A,ES13.6)')
'Abs tol :', real_val
299 call this%solver_factory(this%ksp_vel, this%dm_Xh%size(), &
300 string_val1, integer_val, real_val, logical_val)
301 call this%precon_factory_(this%pc_vel, this%ksp_vel, &
302 this%c_Xh, this%dm_Xh, this%gs_Xh, this%bcs_vel, &
303 string_val2, json_subdict)
309 this%strict_convergence, .false.)
313 call this%ulag%init(this%u, 2)
314 call this%vlag%init(this%v, 2)
315 call this%wlag%init(this%w, 2)
325 call neko_log%section(
'Fluid Source term')
326 call this%source_term%init(this%f_x, this%f_y, this%f_z, this%c_Xh,
user, &
328 call this%source_term%add(params,
'case.fluid.source_terms')
331 end subroutine fluid_scheme_init_base
333 subroutine fluid_scheme_free(this)
334 class(fluid_scheme_incompressible_t),
intent(inout) :: this
338 if (
allocated(this%ksp_vel))
then
339 call this%ksp_vel%free()
340 deallocate(this%ksp_vel)
343 if (
allocated(this%ksp_prs))
then
344 call this%ksp_prs%free()
345 deallocate(this%ksp_prs)
348 if (
allocated(this%pc_vel))
then
349 call precon_destroy(this%pc_vel)
350 deallocate(this%pc_vel)
353 if (
allocated(this%pc_prs))
then
354 call precon_destroy(this%pc_prs)
355 deallocate(this%pc_prs)
358 call this%source_term%free()
360 call this%gs_Xh%free()
362 call this%c_Xh%free()
364 call this%scratch%free()
375 call this%ulag%free()
376 call this%vlag%free()
377 call this%wlag%free()
380 if (
associated(this%f_x))
then
384 if (
associated(this%f_y))
then
388 if (
associated(this%f_z))
then
399 end subroutine fluid_scheme_free
403 subroutine fluid_scheme_validate(this)
404 class(fluid_scheme_incompressible_t),
target,
intent(inout) :: this
406 logical :: logical_val
408 if ( (.not.
associated(this%u)) .or. &
409 (.not.
associated(this%v)) .or. &
410 (.not.
associated(this%w)) .or. &
411 (.not.
associated(this%p)))
then
415 if ( (.not.
allocated(this%u%x)) .or. &
416 (.not.
allocated(this%v%x)) .or. &
417 (.not.
allocated(this%w%x)) .or. &
418 (.not.
allocated(this%p%x)))
then
422 if (.not.
allocated(this%ksp_vel))
then
423 call neko_error(
'No Krylov solver for velocity defined')
426 if (.not.
allocated(this%ksp_prs))
then
427 call neko_error(
'No Krylov solver for pressure defined')
430 end subroutine fluid_scheme_validate
436 subroutine fluid_scheme_bc_apply_vel(this, time, strong)
437 class(fluid_scheme_incompressible_t),
intent(inout) :: this
439 logical,
intent(in) :: strong
441 class(
bc_t),
pointer :: b
444 call this%bcs_vel%apply_vector(&
445 this%u%x, this%v%x, this%w%x, this%dm_Xh%size(), time, strong)
454 call this%bcs_vel%apply_vector(&
455 this%u%x, this%v%x, this%w%x, this%dm_Xh%size(), time, strong)
463 do i = 1, this%bcs_vel%size()
464 b => this%bcs_vel%get(i)
469 end subroutine fluid_scheme_bc_apply_vel
473 subroutine fluid_scheme_bc_apply_prs(this, time)
474 class(fluid_scheme_incompressible_t),
intent(inout) :: this
478 class(
bc_t),
pointer :: b
481 call this%bcs_prs%apply(this%p, time)
485 call this%bcs_prs%apply(this%p, time)
489 do i = 1, this%bcs_prs%size()
490 b => this%bcs_prs%get(i)
495 end subroutine fluid_scheme_bc_apply_prs
499 subroutine fluid_scheme_solver_factory(ksp, n, solver, &
500 max_iter, abstol, monitor)
501 class(
ksp_t),
allocatable,
target,
intent(inout) :: ksp
502 integer,
intent(in),
value :: n
503 character(len=*),
intent(in) :: solver
504 integer,
intent(in) :: max_iter
505 real(kind=
rp),
intent(in) :: abstol
506 logical,
intent(in) :: monitor
508 call krylov_solver_factory(ksp, n, solver, max_iter, abstol, &
511 end subroutine fluid_scheme_solver_factory
514 subroutine fluid_scheme_precon_factory(this, pc, ksp, coef, dof, gs, bclst, &
516 class(fluid_scheme_incompressible_t),
intent(inout) :: this
517 class(
pc_t),
allocatable,
target,
intent(inout) :: pc
518 class(
ksp_t),
target,
intent(inout) :: ksp
519 type(
coef_t),
target,
intent(in) :: coef
520 type(
dofmap_t),
target,
intent(in) :: dof
521 type(
gs_t),
target,
intent(inout) :: gs
522 type(
bc_list_t),
target,
intent(inout) :: bclst
523 character(len=*) :: pctype
524 type(json_file),
intent(inout) :: pcparams
526 call precon_factory(pc, pctype)
528 select type (pcp => pc)
530 call pcp%init(coef, dof, gs)
532 call pcp%init(coef, dof, gs)
534 call pcp%init(coef, dof, gs)
536 call pcp%init(coef, bclst, pcparams)
538 call pcp%init(coef, bclst, pcparams)
543 end subroutine fluid_scheme_precon_factory
546 function fluid_compute_cfl(this, dt)
result(c)
547 class(fluid_scheme_incompressible_t),
intent(in) :: this
548 real(kind=
rp),
intent(in) :: dt
551 c =
cfl(dt, this%u%x, this%v%x, this%w%x, &
552 this%Xh, this%c_Xh, this%msh%nelv, this%msh%gdim)
554 end function fluid_compute_cfl
561 subroutine fluid_scheme_update_material_properties(this, time)
562 class(fluid_scheme_incompressible_t),
intent(inout) :: this
566 call this%user_material_properties(this%name, this%material_properties, &
569 if (len(trim(this%nut_field_name)) > 0)
then
582 call device_memcpy(this%rho%x, this%rho%x_d, this%rho%size(), &
585 end subroutine fluid_scheme_update_material_properties
590 subroutine fluid_scheme_set_material_properties(this, params, user)
591 class(fluid_scheme_incompressible_t),
target,
intent(inout) :: this
592 type(json_file),
intent(inout) :: params
593 type(
user_t),
target,
intent(in) :: user
594 character(len=LOG_SIZE) :: log_buf
597 logical :: nondimensional
598 real(kind=
rp) :: dummy_lambda, dummy_cp
599 real(kind=
rp) :: const_mu, const_rho
612 call this%material_properties%init(2)
613 call this%material_properties%assign(1, this%rho)
614 call this%material_properties%assign(2, this%mu)
616 if (.not.
associated(
user%material_properties, dummy_mp_ptr))
then
618 write(log_buf,
'(A)')
'Material properties must be set in the user' // &
621 this%user_material_properties =>
user%material_properties
623 call user%material_properties(this%name, this%material_properties, &
629 if (params%valid_path(
'case.fluid.Re') .and. &
630 (params%valid_path(
'case.fluid.mu') .or. &
631 params%valid_path(
'case.fluid.rho')))
then
632 call neko_error(
"To set the material properties for the fluid, " // &
633 "either provide Re OR mu and rho in the case file.")
635 else if (params%valid_path(
'case.fluid.Re'))
then
637 write(log_buf,
'(A)')
'Non-dimensional fluid material properties &
640 write(log_buf,
'(A)')
'Density will be set to 1, dynamic viscosity to&
645 call json_get(params,
'case.fluid.Re', const_mu)
646 write(log_buf,
'(A)')
'Read non-dimensional material properties'
648 write(log_buf,
'(A,ES13.6)')
'Re :', const_mu
654 const_mu = 1.0_rp/const_mu
657 call json_get(params,
'case.fluid.mu', const_mu)
658 call json_get(params,
'case.fluid.rho', const_rho)
664 if (
associated(
user%material_properties, dummy_mp_ptr))
then
671 write(log_buf,
'(A,ES13.6)')
'rho :', const_rho
673 write(log_buf,
'(A,ES13.6)')
'mu :', const_mu
685 call device_memcpy(this%rho%x, this%rho%x_d, this%rho%size(), &
689 call device_memcpy(this%mu_tot%x, this%mu_tot%x_d, this%mu%size(), &
692 end subroutine fluid_scheme_set_material_properties
Copy data between host and device (or device and device)
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...
Retrieves a parameter by name or throws an error.
Abstract interface for setting material properties.
Defines a boundary condition.
Jacobi preconditioner accelerator backend.
Device abstraction, common interface for various accelerators.
subroutine, public device_event_sync(event)
Synchronize an event.
integer, parameter, public device_to_host
type(c_ptr), bind(C), public glb_cmd_event
Event for the global command queue.
Defines a dirichlet boundary condition.
Defines a mapping of the degrees of freedom.
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_copy(a, b, n)
Copy a vector .
Defines a registry for storing solution fields.
type(field_registry_t), target, public neko_field_registry
Global field registry.
subroutine fluid_scheme_set_material_properties(this, params, user)
Sets rho and mu.
subroutine fluid_scheme_update_material_properties(this, time)
Call user material properties routine and update the values of mu if necessary.
subroutine fluid_scheme_precon_factory(this, pc, ksp, coef, dof, gs, bclst, pctype, pcparams)
Initialize a Krylov preconditioner.
subroutine fluid_scheme_free(this)
subroutine fluid_scheme_validate(this)
Validate that all fields, solvers etc necessary for performing time-stepping are defined.
real(kind=rp) function fluid_compute_cfl(this, dt)
Compute CFL.
subroutine fluid_scheme_bc_apply_vel(this, time, strong)
Apply all boundary conditions defined for velocity Here we perform additional gs operations to take c...
subroutine fluid_scheme_solver_factory(ksp, n, solver, max_iter, abstol, monitor)
Initialize a linear solver.
subroutine fluid_scheme_bc_apply_prs(this, time)
Apply all boundary conditions defined for pressure.
subroutine fluid_scheme_init_base(this, msh, lx, params, scheme, user, kspv_init)
Initialise a fluid scheme.
Implements the fluid_source_term_t type.
Computes various statistics for the fluid fields. We use the Reynolds decomposition for a field u = ...
Utilities for retrieving parameters from the case files.
subroutine, public json_extract_object(json, name, object)
Extract object as a separate JSON dictionary.
Implements the base abstract type for Krylov solvers plus helper types.
integer, parameter, public ksp_max_iter
Maximum number of iters.
integer, parameter, public neko_log_verbose
Verbose log level.
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(dt, u, v, w, xh, coef, nelv, gdim)
Hybrid ph-multigrid preconditioner.
Defines a registry for storing and requesting temporary fields This can be used when you have a funct...
Implements the source_term_t type and a wrapper source_term_wrapper_t.
Defines a function space.
integer, parameter, public gll
Defines a container for all statistics.
Jacobi preconditioner SX-Aurora backend.
Module with things related to the simulation time.
Interfaces for user interaction with NEKO.
subroutine, public dummy_user_material_properties(scheme_name, properties, time)
Base type for a boundary condition.
A list of allocatable `bc_t`. Follows the standard interface of lists.
Coefficients defined on a given (mesh, ) tuple. Arrays use indices (i,j,k,e): element e,...
Defines a jacobi preconditioner.
Generic Dirichlet boundary condition on .
Base type of all fluid formulations.
Base type of all fluid formulations.
Wrapper contaning and executing the fluid source terms.
Defines a jacobi preconditioner.
Base abstract type for a canonical Krylov method, solving .
Defines a canonical Krylov preconditioner.
Defines a jacobi preconditioner for SX-Aurora.
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...