64 use json_module,
only : json_file
91 type(
gs_t),
pointer :: gs_xh
99 class(
ksp_t),
allocatable :: ksp
101 integer :: ksp_maxiter
103 integer :: projection_dim
105 integer :: projection_activ_step
107 class(pc_t),
allocatable :: pc
124 integer :: n_dir_bcs = 0
126 integer :: n_neumann_bcs = 0
132 type(json_file),
pointer :: params
138 real(kind=
rp),
pointer :: lambda
140 real(kind=
rp),
pointer :: rho
142 real(kind=
rp),
pointer :: cp
144 character(len=NEKO_MSH_MAX_ZLBL_LEN),
allocatable :: bc_labels(:)
176 type(
mesh_t),
target,
intent(inout) :: msh
177 type(
coef_t),
target,
intent(inout) :: coef
178 type(
gs_t),
target,
intent(inout) :: gs
179 type(json_file),
target,
intent(inout) :: params
180 type(
user_t),
target,
intent(in) :: user
192 real(kind=
rp) :: dtlag(10), tlag(10)
213 real(kind=
rp),
intent(inout) :: t
214 integer,
intent(inout) :: tstep
215 real(kind=
rp),
intent(in) :: dt
229 type(
facet_zone_t),
intent(inout) :: zones(NEKO_MSH_MAX_ZLBLS)
230 character(len=NEKO_MSH_MAX_ZLBL_LEN),
intent(in) :: bc_labels(:)
231 character(len=NEKO_MSH_MAX_ZLBL_LEN) :: bc_label
232 integer :: i, j, bc_idx
233 real(kind=
rp) :: dir_value, flux_value
236 do i = 1,
size(bc_labels)
237 bc_label = trim(bc_labels(i))
238 if (bc_label(1:2) .eq.
'd=')
then
254 this%n_dir_bcs = this%n_dir_bcs + 1
255 call this%dir_bcs(this%n_dir_bcs)%init(this%c_Xh)
256 call this%dir_bcs(this%n_dir_bcs)%mark_zone(zones(i))
257 read(bc_label(3:), *) dir_value
258 call this%dir_bcs(this%n_dir_bcs)%set_g(dir_value)
262 if (bc_label(1:2) .eq.
'n=')
then
263 this%n_neumann_bcs = this%n_neumann_bcs + 1
264 call this%neumann_bcs(this%n_neumann_bcs)%init(this%c_Xh)
265 call this%neumann_bcs(this%n_neumann_bcs)%mark_zone(zones(i))
266 read(bc_label(3:), *) flux_value
267 call this%neumann_bcs(this%n_neumann_bcs)%init_neumann(flux_value)
271 if (bc_label(1:4) .eq.
'user')
then
272 call this%user_bc%mark_zone(zones(i))
277 do i = 1, this%n_dir_bcs
278 call this%dir_bcs(i)%finalize()
279 call bc_list_add(this%bclst_dirichlet, this%dir_bcs(i))
283 call bc_list_init(this%bclst_neumann, this%n_neumann_bcs)
284 do i=1, this%n_neumann_bcs
285 call this%neumann_bcs(i)%finalize()
286 call bc_list_add(this%bclst_neumann, this%neumann_bcs(i))
301 type(
mesh_t),
target,
intent(inout) :: msh
302 type(
coef_t),
target,
intent(inout) :: c_Xh
303 type(
gs_t),
target,
intent(inout) :: gs_Xh
304 type(json_file),
target,
intent(inout) :: params
305 character(len=*),
intent(in) :: scheme
306 type(
user_t),
target,
intent(in) :: user
309 character(len=LOG_SIZE) :: log_buf
311 logical :: logical_val
312 real(kind=
rp) :: real_val, solver_abstol
313 integer :: integer_val, ierr
314 character(len=:),
allocatable :: solver_type, solver_precon
321 call json_get(params,
'case.fluid.velocity_solver.type', solver_type)
322 call json_get(params,
'case.fluid.velocity_solver.preconditioner',&
324 call json_get(params,
'case.fluid.velocity_solver.absolute_tolerance',&
334 write(log_buf,
'(A,ES13.6)')
'rho :', this%rho
336 write(log_buf,
'(A,ES13.6)')
'lambda :', this%lambda
338 write(log_buf,
'(A,ES13.6)')
'cp :', this%cp
342 'case.fluid.velocity_solver.projection_space_size',&
343 this%projection_dim, 20)
345 'case.fluid.velocity_solver.projection_hold_steps',&
346 this%projection_activ_step, 5)
349 write(log_buf,
'(A, A)')
'Type : ', trim(scheme)
351 call neko_log%message(
'Ksp scalar : ('// trim(solver_type) // &
352 ', ' // trim(solver_precon) //
')')
353 write(log_buf,
'(A,ES13.6)')
' `-abs tol :', solver_abstol
357 this%dm_Xh => this%u%dof
358 this%params => params
365 call this%slag%init(this%s, 2)
375 call this%user_bc%init(this%c_Xh)
381 this%bc_labels =
"not"
383 if (params%valid_path(
'case.scalar.boundary_types'))
then
385 'case.scalar.boundary_types', &
394 call this%f_Xh%init(this%dm_Xh, fld_name=
"scalar_rhs")
397 call this%source_term%init(params, this%f_Xh, this%c_Xh, user)
402 call this%user_bc%mark_zone(msh%wall)
403 call this%user_bc%mark_zone(msh%inlet)
404 call this%user_bc%mark_zone(msh%outlet)
405 call this%user_bc%mark_zone(msh%outlet_normal)
406 call this%user_bc%mark_zone(msh%sympln)
407 call this%user_bc%finalize()
408 call this%user_bc%set_coef(this%c_Xh)
409 if (this%user_bc%msk(0) .gt. 0)
call bc_list_add(this%bclst_dirichlet,&
413 call this%field_dir_bc%init(this%c_Xh)
414 call this%field_dir_bc%mark_zones_from_list(msh%labeled_zones, &
415 'd_s', this%bc_labels)
416 call this%field_dir_bc%finalize()
417 call mpi_allreduce(this%field_dir_bc%msk(0), integer_val, 1, &
419 if (integer_val .gt. 0)
call this%field_dir_bc%init_field(
'd_s')
421 call bc_list_add(this%bclst_dirichlet, this%field_dir_bc)
426 this%dirichlet_update_ => user%user_dirichlet_update
431 allocate(this%field_dirichlet_fields%fields(1))
432 this%field_dirichlet_fields%fields(1)%f => &
433 this%field_dir_bc%field_bc
436 call bc_list_add(this%field_dirichlet_bcs, this%field_dir_bc)
443 solver_type, integer_val, solver_abstol)
445 this%c_Xh, this%dm_Xh, this%gs_Xh, this%bclst_dirichlet, solver_precon)
462 if (
allocated(this%ksp))
then
467 if (
allocated(this%pc))
then
472 if (
allocated(this%bc_labels))
then
473 deallocate(this%bc_labels)
476 call this%source_term%free()
481 call this%slag%free()
484 call this%field_dirichlet_fields%free()
486 call this%field_dir_bc%field_bc%free()
487 call this%field_dir_bc%free()
488 if (
associated(this%dirichlet_update_))
then
489 this%dirichlet_update_ => null()
499 if ( (.not.
allocated(this%u%x)) .or. &
500 (.not.
allocated(this%v%x)) .or. &
501 (.not.
allocated(this%w%x)) .or. &
502 (.not.
allocated(this%s%x)))
then
506 if (.not.
allocated(this%ksp))
then
507 call neko_error(
'No Krylov solver for velocity defined')
510 if (.not.
associated(this%Xh))
then
514 if (.not.
associated(this%dm_Xh))
then
518 if (.not.
associated(this%c_Xh))
then
522 if (.not.
associated(this%f_Xh))
then
526 if (.not.
associated(this%params))
then
541 class(
ksp_t),
allocatable,
target,
intent(inout) :: ksp
542 integer,
intent(in),
value :: n
543 integer,
intent(in) :: max_iter
544 character(len=*),
intent(in) :: solver
545 real(kind=
rp) :: abstol
553 class(pc_t),
allocatable,
target,
intent(inout) :: pc
554 class(
ksp_t),
target,
intent(inout) :: ksp
555 type(
coef_t),
target,
intent(inout) :: coef
556 type(
dofmap_t),
target,
intent(inout) :: dof
557 type(
gs_t),
target,
intent(inout) :: gs
558 type(
bc_list_t),
target,
intent(inout) :: bclst
559 character(len=*) :: pctype
563 select type(pcp => pc)
565 call pcp%init(coef, dof, gs)
567 call pcp%init(coef, dof, gs)
569 call pcp%init(coef, dof, gs)
571 if (len_trim(pctype) .gt. 4)
then
572 if (index(pctype,
'+') .eq. 5)
then
573 call pcp%init(dof%msh, dof%Xh, coef, dof, gs, &
574 bclst, trim(pctype(6:)))
579 call pcp%init(dof%msh, dof%Xh, coef, dof, gs, bclst)
593 call this%user_bc%set_eval(usr_eval)
Abstract interface defining a dirichlet condition on a list of fields.
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 to dealocate a scalar formulation.
Abstract interface to initialize a scalar formulation.
Abstract interface to restart a scalar formulation.
Abstract interface to compute a time-step.
Abstract interface defining a user defined scalar boundary condition (pointwise) Just imitating inflo...
Defines a boundary condition.
subroutine, public bc_list_add(bclst, bc)
Add a condition to a list of boundary conditions.
subroutine, public bc_list_init(bclst, size)
Constructor for a list of boundary conditions.
subroutine, public bc_list_free(bclst)
Destructor for a list of boundary conditions.
subroutine, public bc_list_apply_scalar(bclst, x, n, t, tstep)
Apply a list of boundary conditions to a scalar field.
type(mpi_comm) neko_comm
MPI communicator.
Jacobi preconditioner accelerator backend.
Defines a dirichlet boundary condition.
Defines a mapping of the degrees of freedom.
Defines a zone as a subset of facets in a mesh.
Defines inflow dirichlet conditions.
Defines a registry for storing solution fields.
type(field_registry_t), target, public neko_field_registry
Global field registry.
Utilities for retrieving parameters from the case files.
subroutine, public krylov_solver_destroy(ksp)
Destroy an interative Krylov solver.
subroutine, public krylov_solver_factory(ksp, n, solver, max_iter, abstol, M)
Initialize an iterative Krylov solver.
Implements the base abstract type for Krylov solvers plus helper types.
type(log_t), public neko_log
Global log stream.
integer, parameter, public log_size
Implements material_properties_t type.
integer, parameter, public neko_msh_max_zlbls
Max num. zone labels.
integer, parameter, public neko_msh_max_zlbl_len
Max length of a zone label.
Defines a Neumann boundary condition.
integer, parameter, public rp
Global precision used in computations.
subroutine precon_destroy(pc)
Destroy a preconditioner.
subroutine precon_factory(pc, pctype)
Create a preconditioner.
Contains the scalar_scheme_t type.
subroutine scalar_scheme_init(this, msh, c_Xh, gs_Xh, params, scheme, user, material_properties)
Initialize all related components of the current scheme.
subroutine scalar_scheme_add_bcs(this, zones, bc_labels)
Initialize boundary conditions.
subroutine scalar_scheme_free(this)
Deallocate a scalar formulation.
subroutine scalar_scheme_set_user_bc(this, usr_eval)
Initialize a user defined scalar bc.
subroutine scalar_scheme_validate(this)
Validate that all fields, solvers etc necessary for performing time-stepping are defined.
subroutine scalar_scheme_precon_factory(pc, ksp, coef, dof, gs, bclst, pctype)
Initialize a Krylov preconditioner.
subroutine scalar_scheme_solver_factory(ksp, n, solver, max_iter, abstol)
Initialize a linear solver.
Implements the scalar_source_term_t type.
Implements the source_term_t type and a wrapper source_term_wrapper_t.
Defines a function space.
Jacobi preconditioner SX-Aurora backend.
Compound scheme for the advection and diffusion operators in a transport equation.
Implements type time_step_controller.
Interfaces for user interaction with NEKO.
Defines dirichlet conditions for scalars.
A list of boundary conditions.
Base type for a boundary condition.
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 .
User defined dirichlet condition, for which the user can work with an entire field....
field_list_t, To be able to group fields together
Defines a jacobi preconditioner.
Base abstract type for a canonical Krylov method, solving .
Contains all the material properties necessary in the simulation.
Generic Neumann boundary condition. This sets the flux of the field to the chosen value.
Base type for a scalar advection-diffusion solver.
Wrapper contaning and executing the scalar source terms.
The function space for the SEM solution fields.
Defines a jacobi preconditioner for SX-Aurora.
Implements the logic to compute the time coefficients for the advection and diffusion operators in a ...
Provides a tool to set time step dt.
User defined dirichlet condition for scalars.