46 use krylov,
only :
ksp_t, krylov_solver_factory, krylov_solver_destroy, &
62 use precon,
only :
pc_t, precon_factory, precon_destroy
74 use json_module,
only : json_file
107 class(
ksp_t),
allocatable :: ksp_vel
108 class(
ksp_t),
allocatable :: ksp_prs
109 class(
pc_t),
allocatable :: pc_vel
110 class(
pc_t),
allocatable :: pc_prs
111 integer :: vel_projection_dim
112 integer :: pr_projection_dim
113 integer :: vel_projection_activ_step
114 integer :: pr_projection_activ_step
116 class(
bc_t),
allocatable :: bc_inflow
119 logical :: if_gradient_jump_penalty
135 type(json_file),
pointer :: params
141 logical :: forced_flow_rate = .false.
142 logical :: freeze = .false.
148 character(len=:),
allocatable :: nut_field_name
150 logical :: variable_material_properties = .false.
157 character(len=NEKO_MSH_MAX_ZLBL_LEN),
allocatable :: bc_labels(:)
176 procedure, pass(this) :: set_material_properties => &
186 procedure,
private, pass(this) :: set_bc_type_output => &
189 procedure, pass(this) :: update_material_properties => &
203 type(
mesh_t),
target,
intent(inout) :: msh
204 integer,
intent(inout) :: lx
205 type(json_file),
target,
intent(inout) :: params
206 type(
user_t),
target,
intent(in) :: user
228 real(kind=
rp),
intent(inout) :: t
229 integer,
intent(inout) :: tstep
230 real(kind=
rp),
intent(in) :: dt
242 real(kind=
rp) :: dtlag(10), tlag(10)
249 module subroutine fluid_scheme_factory(object, type_name)
250 class(fluid_scheme_t),
intent(inout),
allocatable :: object
251 character(len=*) :: type_name
252 end subroutine fluid_scheme_factory
255 public :: fluid_scheme_t, fluid_scheme_factory
260 subroutine fluid_scheme_init_common(this, msh, lx, params, scheme, user, &
263 class(fluid_scheme_t),
target,
intent(inout) :: this
264 type(
mesh_t),
target,
intent(inout) :: msh
265 integer,
intent(inout) :: lx
266 character(len=*),
intent(in) :: scheme
267 type(json_file),
target,
intent(inout) :: params
268 type(
user_t),
target,
intent(in) :: user
269 logical,
intent(in) :: kspv_init
271 character(len=LOG_SIZE) :: log_buf
272 real(kind=
rp),
allocatable :: real_vec(:)
273 real(kind=
rp) :: real_val, kappa, b, z0
274 logical :: logical_val
275 integer :: integer_val, ierr
276 type(json_file) :: wm_json
277 character(len=:),
allocatable :: string_val1, string_val2
285 if (msh%gdim .eq. 2)
then
286 call this%Xh%init(
gll, lx, lx)
288 call this%Xh%init(
gll, lx, lx, lx)
291 call this%dm_Xh%init(msh, this%Xh)
293 call this%gs_Xh%init(this%dm_Xh)
295 call this%c_Xh%init(this%gs_Xh)
301 this%params => params
309 write(log_buf,
'(A, A)')
'Type : ', trim(scheme)
315 call this%set_material_properties(params, user)
320 if (params%valid_path(
'case.fluid.nut_field'))
then
321 call json_get(params,
'case.fluid.nut_field', this%nut_field_name)
322 this%variable_material_properties = .true.
324 this%nut_field_name =
""
329 call this%mu_field%init(this%dm_Xh,
"mu")
330 call this%rho_field%init(this%dm_Xh,
"mu")
331 call field_cfill(this%mu_field, this%mu, this%mu_field%size())
332 call field_cfill(this%rho_field, this%rho, this%mu_field%size())
338 call cfill(this%mu_field%x, this%mu, this%mu_field%size())
339 call cfill(this%rho_field%x, this%rho, this%rho_field%size())
345 'case.fluid.velocity_solver.projection_space_size', &
346 this%vel_projection_dim, 20)
348 'case.fluid.pressure_solver.projection_space_size', &
349 this%pr_projection_dim, 20)
351 'case.fluid.velocity_solver.projection_hold_steps', &
352 this%vel_projection_activ_step, 5)
354 'case.fluid.pressure_solver.projection_hold_steps', &
355 this%pr_projection_activ_step, 5)
360 if (params%valid_path(
"case.fluid.flow_rate_force"))
then
361 this%forced_flow_rate = .true.
366 write(log_buf,
'(A, I1)')
'Poly order : ', lx-1
367 else if (lx .ge. 10)
then
368 write(log_buf,
'(A, I2)')
'Poly order : ', lx-1
370 write(log_buf,
'(A, I3)')
'Poly order : ', lx-1
373 write(log_buf,
'(A, I0)')
'DoFs : ', this%dm_Xh%size()
376 write(log_buf,
'(A,ES13.6)')
'rho :', this%rho
378 write(log_buf,
'(A,ES13.6)')
'mu :', this%mu
381 call json_get(params,
'case.numerics.dealias', logical_val)
382 write(log_buf,
'(A, L1)')
'Dealias : ', logical_val
387 write(log_buf,
'(A, L1)')
'Save bdry : ', logical_val
389 if (logical_val)
then
390 write(log_buf,
'(A)')
'bdry keys: '
392 write(log_buf,
'(A)')
'No-slip wall ''w'' = 1'
394 write(log_buf,
'(A)')
'Inlet/velocity dirichlet ''v'' = 2'
396 write(log_buf,
'(A)')
'Outlet ''o'' = 3'
398 write(log_buf,
'(A)')
'Symmetry ''sym'' = 4'
400 write(log_buf,
'(A)')
'Outlet-normal ''on'' = 5'
402 write(log_buf,
'(A)')
'Periodic = 6'
404 write(log_buf,
'(A)')
'Dirichlet on velocity components: '
406 write(log_buf,
'(A)')
' ''d_vel_u, d_vel_v, d_vel_w'' = 7'
408 write(log_buf,
'(A)')
'Pressure dirichlet ''d_pres'' = 8'
410 write(log_buf,
'(A)')
'''d_vel_(u,v,w)'' and ''d_pres'' = 8'
412 write(log_buf,
'(A)')
'Shear stress ''sh'' = 9'
414 write(log_buf,
'(A)')
'Wall modelling ''wm'' = 10'
416 write(log_buf,
'(A)')
'No boundary condition set = 0'
425 this%bc_labels =
"not"
427 if (params%valid_path(
'case.fluid.boundary_types'))
then
429 'case.fluid.boundary_types', &
435 call this%bc_sym%init_base(this%c_Xh)
436 call this%bc_sym%mark_zone(msh%sympln)
437 call this%bc_sym%mark_zones_from_list(msh%labeled_zones, &
438 'sym', this%bc_labels)
439 call this%bc_sym%finalize()
440 call this%bc_sym%init(this%c_Xh)
444 call this%bc_sh%init_base(this%c_Xh)
445 call this%bc_sh%mark_zones_from_list(msh%labeled_zones, &
446 'sh', this%bc_labels)
447 call this%bc_sh%finalize()
450 call this%bc_sh%init_shear_stress(this%c_Xh)
452 call bc_list_add(this%bclst_vel, this%bc_sh%symmetry)
455 if (this%bc_sh%msk(0) .gt. 0)
then
456 call params%get(
'case.fluid.shear_stress.value', real_vec, logical_val)
457 if (.not. logical_val .and. this%bc_sh%msk(0) .gt. 0)
then
458 call neko_warning(
"No stress values provided for sh boundaries, &
459 & defaulting to 0. Use fluid.shear_stress.value to set.")
460 allocate(real_vec(3))
462 else if (
size(real_vec) .ne. 3)
then
463 call neko_error (
"The shear stress vector provided in &
464 &fluid.shear_stress.value should have 3 components.")
466 call this%bc_sh%set_stress(real_vec(1), real_vec(2), real_vec(3))
470 call bc_list_add(this%bclst_vel_neumann, this%bc_sh)
475 if (params%valid_path(
'case.fluid.inflow_condition'))
then
476 call json_get(params,
'case.fluid.inflow_condition.type', string_val1)
477 if (trim(string_val1) .eq.
"uniform")
then
479 else if (trim(string_val1) .eq.
"blasius")
then
481 else if (trim(string_val1) .eq.
"user")
then
484 call neko_error(
'Invalid inflow condition '//string_val1)
487 call this%bc_inflow%init_base(this%c_Xh)
488 call this%bc_inflow%mark_zone(msh%inlet)
489 call this%bc_inflow%mark_zones_from_list(msh%labeled_zones, &
491 call this%bc_inflow%finalize()
494 if (trim(string_val1) .eq.
"uniform")
then
495 call json_get(params,
'case.fluid.inflow_condition.value', real_vec)
496 select type (bc_if => this%bc_inflow)
498 call bc_if%set_inflow(real_vec)
500 else if (trim(string_val1) .eq.
"blasius")
then
501 select type (bc_if => this%bc_inflow)
503 call json_get(params,
'case.fluid.blasius.delta', real_val)
504 call json_get(params,
'case.fluid.blasius.approximation', &
506 call json_get(params,
'case.fluid.blasius.freestream_velocity', &
509 call bc_if%set_params(real_vec, real_val, string_val2)
512 else if (trim(string_val1) .eq.
"user")
then
520 call this%bc_wallmodel%init_base(this%c_Xh)
521 call this%bc_wallmodel%mark_zones_from_list(msh%labeled_zones,&
522 'wm', this%bc_labels)
523 call this%bc_wallmodel%finalize()
525 if (this%bc_wallmodel%msk(0) .gt. 0)
then
527 call this%bc_wallmodel%init_wall_model_bc(wm_json, this%mu / this%rho)
529 call this%bc_wallmodel%shear_stress_t%init_shear_stress(this%c_Xh)
532 call bc_list_add(this%bclst_vel, this%bc_wallmodel%symmetry)
533 call bc_list_add(this%bclst_vel_neumann, this%bc_wallmodel)
535 call this%bc_wall%init_base(this%c_Xh)
536 call this%bc_wall%mark_zone(msh%wall)
537 call this%bc_wall%mark_zones_from_list(msh%labeled_zones, &
539 call this%bc_wall%finalize()
543 call this%user_field_bc_vel%bc_u%init_base(this%c_Xh)
544 call this%user_field_bc_vel%bc_u%mark_zones_from_list(msh%labeled_zones, &
545 'd_vel_u', this%bc_labels)
546 call this%user_field_bc_vel%bc_u%finalize()
548 call mpi_allreduce(this%user_field_bc_vel%bc_u%msk(0), integer_val, 1, &
550 if (integer_val .gt. 0)
then
551 call this%user_field_bc_vel%bc_u%init_field(
'd_vel_u')
555 call this%user_field_bc_vel%bc_v%init_base(this%c_Xh)
556 call this%user_field_bc_vel%bc_v%mark_zones_from_list(msh%labeled_zones, &
557 'd_vel_v', this%bc_labels)
558 call this%user_field_bc_vel%bc_v%finalize()
560 call mpi_allreduce(this%user_field_bc_vel%bc_v%msk(0), integer_val, 1, &
562 if (integer_val .gt. 0)
then
563 call this%user_field_bc_vel%bc_v%init_field(
'd_vel_v')
567 call this%user_field_bc_vel%bc_w%init_base(this%c_Xh)
568 call this%user_field_bc_vel%bc_w%mark_zones_from_list(msh%labeled_zones, &
569 'd_vel_w', this%bc_labels)
570 call this%user_field_bc_vel%bc_w%finalize()
572 call mpi_allreduce(this%user_field_bc_vel%bc_w%msk(0), integer_val, 1, &
574 if (integer_val .gt. 0)
then
575 call this%user_field_bc_vel%bc_w%init_field(
'd_vel_w')
579 call this%user_field_bc_vel%init_base(this%c_Xh)
580 call this%user_field_bc_vel%mark_zones_from_list(msh%labeled_zones, &
581 'd_vel_u', this%bc_labels)
582 call this%user_field_bc_vel%mark_zones_from_list(msh%labeled_zones, &
583 'd_vel_v', this%bc_labels)
584 call this%user_field_bc_vel%mark_zones_from_list(msh%labeled_zones, &
585 'd_vel_w', this%bc_labels)
586 call this%user_field_bc_vel%finalize()
589 call bc_list_add(this%bclst_vel, this%user_field_bc_vel)
594 this%user_field_bc_vel%update => user%user_dirichlet_update
601 call this%user_field_bc_vel%field_list%init(4)
602 call this%user_field_bc_vel%field_list%assign_to_field(1, &
603 this%user_field_bc_vel%bc_u%field_bc)
604 call this%user_field_bc_vel%field_list%assign_to_field(2, &
605 this%user_field_bc_vel%bc_v%field_bc)
606 call this%user_field_bc_vel%field_list%assign_to_field(3, &
607 this%user_field_bc_vel%bc_w%field_bc)
608 call this%user_field_bc_vel%field_list%assign_to_field(4, &
609 this%user_field_bc_prs%field_bc)
611 call bc_list_init(this%user_field_bc_vel%bc_list,
size = 4)
614 this%user_field_bc_vel%bc_u)
616 this%user_field_bc_vel%bc_v)
618 this%user_field_bc_vel%bc_w)
622 call fluid_scheme_set_bc_type_output(this, params)
630 call this%f_x%init(this%dm_Xh, fld_name =
"fluid_rhs_x")
631 call this%f_y%init(this%dm_Xh, fld_name =
"fluid_rhs_y")
632 call this%f_z%init(this%dm_Xh, fld_name =
"fluid_rhs_z")
635 call this%source_term%init(this%f_x, this%f_y, this%f_z, this%c_Xh, user)
636 call this%source_term%add(params,
'case.fluid.source_terms')
640 call this%set_bc_type_output(params)
644 call neko_log%section(
"Velocity solver")
646 'case.fluid.velocity_solver.max_iterations', &
648 call json_get(params,
'case.fluid.velocity_solver.type', string_val1)
649 call json_get(params,
'case.fluid.velocity_solver.preconditioner', &
651 call json_get(params,
'case.fluid.velocity_solver.absolute_tolerance', &
654 'case.fluid.velocity_solver.monitor', &
655 logical_val, .false.)
657 call neko_log%message(
'Type : ('// trim(string_val1) // &
658 ', ' // trim(string_val2) //
')')
660 write(log_buf,
'(A,ES13.6)')
'Abs tol :', real_val
662 call fluid_scheme_solver_factory(this%ksp_vel, this%dm_Xh%size(), &
663 string_val1, integer_val, real_val, logical_val)
664 call fluid_scheme_precon_factory(this%pc_vel, this%ksp_vel, &
665 this%c_Xh, this%dm_Xh, this%gs_Xh, this%bclst_vel, string_val2)
678 call this%ulag%init(this%u, 2)
679 call this%vlag%init(this%v, 2)
680 call this%wlag%init(this%w, 2)
683 end subroutine fluid_scheme_init_common
686 subroutine fluid_scheme_init_all(this, msh, lx, params, kspv_init, &
687 kspp_init, scheme, user)
689 class(fluid_scheme_t),
target,
intent(inout) :: this
690 type(
mesh_t),
target,
intent(inout) :: msh
691 integer,
intent(inout) :: lx
692 type(json_file),
target,
intent(inout) :: params
693 type(
user_t),
target,
intent(in) :: user
696 character(len=*),
intent(in) :: scheme
697 real(kind=
rp) :: abs_tol
698 integer :: integer_val, ierr
699 logical :: logical_val
700 character(len=:),
allocatable :: solver_type, precon_type
701 character(len=LOG_SIZE) :: log_buf
702 real(kind=
rp) :: gjp_param_a, gjp_param_b
704 call fluid_scheme_init_common(this, msh, lx, params, scheme, user, &
714 call this%bc_prs%init_base(this%c_Xh)
715 call this%bc_prs%mark_zones_from_list(msh%labeled_zones, &
717 call this%bc_prs%mark_zones_from_list(msh%labeled_zones, &
718 'on', this%bc_labels)
721 call this%user_field_bc_prs%init_base(this%c_Xh)
722 call this%user_field_bc_prs%mark_zones_from_list(msh%labeled_zones, &
723 'd_pres', this%bc_labels)
724 call this%user_field_bc_prs%finalize()
725 call mpi_allreduce(this%user_field_bc_prs%msk(0), integer_val, 1, &
728 if (integer_val .gt. 0)
call this%user_field_bc_prs%init_field(
'd_pres')
729 call bc_list_add(this%bclst_prs, this%user_field_bc_prs)
730 call bc_list_add(this%user_field_bc_vel%bc_list, this%user_field_bc_prs)
732 if (msh%outlet%size .gt. 0)
then
733 call this%bc_prs%mark_zone(msh%outlet)
735 if (msh%outlet_normal%size .gt. 0)
then
736 call this%bc_prs%mark_zone(msh%outlet_normal)
739 call this%bc_prs%finalize()
740 call this%bc_prs%set_g(0.0_rp)
742 call this%bc_dong%init_base(this%c_Xh)
743 call this%bc_dong%mark_zones_from_list(msh%labeled_zones, &
744 'o+dong', this%bc_labels)
745 call this%bc_dong%mark_zones_from_list(msh%labeled_zones, &
746 'on+dong', this%bc_labels)
747 call this%bc_dong%finalize()
750 call this%bc_dong%init(this%c_Xh, params)
756 call neko_log%section(
"Pressure solver")
759 'case.fluid.pressure_solver.max_iterations', &
761 call json_get(params,
'case.fluid.pressure_solver.type', solver_type)
762 call json_get(params,
'case.fluid.pressure_solver.preconditioner', &
764 call json_get(params,
'case.fluid.pressure_solver.absolute_tolerance', &
767 'case.fluid.pressure_solver.monitor', &
768 logical_val, .false.)
769 call neko_log%message(
'Type : ('// trim(solver_type) // &
770 ', ' // trim(precon_type) //
')')
771 write(log_buf,
'(A,ES13.6)')
'Abs tol :', abs_tol
774 call fluid_scheme_solver_factory(this%ksp_prs, this%dm_Xh%size(), &
775 solver_type, integer_val, abs_tol, logical_val)
776 call fluid_scheme_precon_factory(this%pc_prs, this%ksp_prs, &
777 this%c_Xh, this%dm_Xh, this%gs_Xh, this%bclst_prs, precon_type)
785 'case.fluid.gradient_jump_penalty.enabled',&
786 this%if_gradient_jump_penalty, .false.)
788 if (this%if_gradient_jump_penalty .eqv. .true.)
then
789 if ((this%dm_Xh%xh%lx - 1) .eq. 1)
then
791 'case.fluid.gradient_jump_penalty.tau',&
792 gjp_param_a, 0.02_rp)
796 'case.fluid.gradient_jump_penalty.scaling_factor',&
799 'case.fluid.gradient_jump_penalty.scaling_exponent',&
802 call this%gradient_jump_penalty_u%init(params, this%dm_Xh, this%c_Xh, &
803 gjp_param_a, gjp_param_b)
804 call this%gradient_jump_penalty_v%init(params, this%dm_Xh, this%c_Xh, &
805 gjp_param_a, gjp_param_b)
806 call this%gradient_jump_penalty_w%init(params, this%dm_Xh, this%c_Xh, &
807 gjp_param_a, gjp_param_b)
812 end subroutine fluid_scheme_init_all
815 subroutine fluid_scheme_free(this)
816 class(fluid_scheme_t),
intent(inout) :: this
818 call this%bdry%free()
820 if (
allocated(this%bc_inflow))
then
821 call this%bc_inflow%free()
824 call this%bc_wall%free()
825 call this%bc_sym%free()
826 call this%bc_sh%free()
831 call this%user_field_bc_prs%field_bc%free()
832 call this%user_field_bc_prs%free()
833 call this%user_field_bc_vel%bc_u%field_bc%free()
834 call this%user_field_bc_vel%bc_v%field_bc%free()
835 call this%user_field_bc_vel%bc_w%field_bc%free()
836 call this%user_field_bc_vel%free()
840 if (
allocated(this%ksp_vel))
then
841 call krylov_solver_destroy(this%ksp_vel)
842 deallocate(this%ksp_vel)
845 if (
allocated(this%ksp_prs))
then
846 call krylov_solver_destroy(this%ksp_prs)
847 deallocate(this%ksp_prs)
850 if (
allocated(this%pc_vel))
then
851 call precon_destroy(this%pc_vel)
852 deallocate(this%pc_vel)
855 if (
allocated(this%pc_prs))
then
856 call precon_destroy(this%pc_prs)
857 deallocate(this%pc_prs)
860 if (
allocated(this%bc_labels))
then
861 deallocate(this%bc_labels)
864 call this%source_term%free()
866 call this%gs_Xh%free()
868 call this%c_Xh%free()
872 call this%scratch%free()
881 call this%ulag%free()
882 call this%vlag%free()
883 call this%wlag%free()
886 if (
associated(this%f_x))
then
890 if (
associated(this%f_y))
then
894 if (
associated(this%f_z))
then
902 call this%mu_field%free()
905 if (this%if_gradient_jump_penalty .eqv. .true.)
then
906 call this%gradient_jump_penalty_u%free()
907 call this%gradient_jump_penalty_v%free()
908 call this%gradient_jump_penalty_w%free()
912 end subroutine fluid_scheme_free
916 subroutine fluid_scheme_validate(this)
917 class(fluid_scheme_t),
target,
intent(inout) :: this
919 logical :: logical_val
921 if ( (.not.
associated(this%u)) .or. &
922 (.not.
associated(this%v)) .or. &
923 (.not.
associated(this%w)) .or. &
924 (.not.
associated(this%p)))
then
928 if ( (.not.
allocated(this%u%x)) .or. &
929 (.not.
allocated(this%v%x)) .or. &
930 (.not.
allocated(this%w%x)) .or. &
931 (.not.
allocated(this%p%x)))
then
935 if (.not.
allocated(this%ksp_vel))
then
936 call neko_error(
'No Krylov solver for velocity defined')
939 if (.not.
allocated(this%ksp_prs))
then
940 call neko_error(
'No Krylov solver for pressure defined')
943 if (.not.
associated(this%params))
then
947 if (
allocated(this%bc_inflow))
then
948 select type (ip => this%bc_inflow)
957 call this%chkp%init(this%u, this%v, this%w, this%p)
959 end subroutine fluid_scheme_validate
965 subroutine fluid_scheme_bc_apply_vel(this, t, tstep)
966 class(fluid_scheme_t),
intent(inout) :: this
967 real(kind=
rp),
intent(in) :: t
968 integer,
intent(in) :: tstep
971 this%u%x, this%v%x, this%w%x, this%dm_Xh%size(), t, tstep)
973 end subroutine fluid_scheme_bc_apply_vel
977 subroutine fluid_scheme_bc_apply_prs(this, t, tstep)
978 class(fluid_scheme_t),
intent(inout) :: this
979 real(kind=
rp),
intent(in) :: t
980 integer,
intent(in) :: tstep
983 this%p%dof%size(), t, tstep)
985 end subroutine fluid_scheme_bc_apply_prs
989 subroutine fluid_scheme_solver_factory(ksp, n, solver, &
990 max_iter, abstol, monitor)
991 class(
ksp_t),
allocatable,
target,
intent(inout) :: ksp
992 integer,
intent(in),
value :: n
993 character(len=*),
intent(in) :: solver
994 integer,
intent(in) :: max_iter
995 real(kind=
rp),
intent(in) :: abstol
996 logical,
intent(in) :: monitor
998 call krylov_solver_factory(ksp, n, solver, max_iter, abstol, &
1001 end subroutine fluid_scheme_solver_factory
1004 subroutine fluid_scheme_precon_factory(pc, ksp, coef, dof, gs, bclst, &
1006 class(
pc_t),
allocatable,
target,
intent(inout) :: pc
1007 class(
ksp_t),
target,
intent(inout) :: ksp
1008 type(
coef_t),
target,
intent(inout) :: coef
1009 type(
dofmap_t),
target,
intent(inout) :: dof
1010 type(
gs_t),
target,
intent(inout) :: gs
1011 type(
bc_list_t),
target,
intent(inout) :: bclst
1012 character(len=*) :: pctype
1014 call precon_factory(pc, pctype)
1016 select type (pcp => pc)
1018 call pcp%init(coef, dof, gs)
1020 call pcp%init(coef, dof, gs)
1022 call pcp%init(coef, dof, gs)
1024 if (len_trim(pctype) .gt. 4)
then
1025 if (index(pctype,
'+') .eq. 5)
then
1026 call pcp%init(dof%msh, dof%Xh, coef, dof, gs, bclst, &
1029 call neko_error(
'Unknown coarse grid solver')
1032 call pcp%init(dof%msh, dof%Xh, coef, dof, gs, bclst)
1038 end subroutine fluid_scheme_precon_factory
1041 subroutine fluid_scheme_set_usr_inflow(this, usr_eval)
1042 class(fluid_scheme_t),
intent(inout) :: this
1045 select type (bc_if => this%bc_inflow)
1047 call bc_if%set_eval(usr_eval)
1049 call neko_error(
"Not a user defined inflow condition")
1051 end subroutine fluid_scheme_set_usr_inflow
1054 function fluid_compute_cfl(this, dt)
result(c)
1055 class(fluid_scheme_t),
intent(in) :: this
1056 real(kind=
rp),
intent(in) :: dt
1059 c =
cfl(dt, this%u%x, this%v%x, this%w%x, &
1060 this%Xh, this%c_Xh, this%msh%nelv, this%msh%gdim)
1062 end function fluid_compute_cfl
1066 subroutine fluid_scheme_set_bc_type_output(this, params)
1067 class(fluid_scheme_t),
target,
intent(inout) :: this
1068 type(json_file),
intent(inout) :: params
1078 call this%bdry%init(this%dm_Xh,
'bdry')
1081 call bdry_mask%init_base(this%c_Xh)
1082 call bdry_mask%mark_zone(this%msh%wall)
1083 call bdry_mask%mark_zones_from_list(this%msh%labeled_zones, &
1084 'w', this%bc_labels)
1085 call bdry_mask%finalize()
1086 call bdry_mask%set_g(1.0_rp)
1087 call bdry_mask%apply_scalar(this%bdry%x, this%dm_Xh%size())
1088 call bdry_mask%free()
1090 call bdry_mask%init_base(this%c_Xh)
1091 call bdry_mask%mark_zone(this%msh%inlet)
1092 call bdry_mask%mark_zones_from_list(this%msh%labeled_zones, &
1093 'v', this%bc_labels)
1094 call bdry_mask%finalize()
1095 call bdry_mask%set_g(2.0_rp)
1096 call bdry_mask%apply_scalar(this%bdry%x, this%dm_Xh%size())
1097 call bdry_mask%free()
1099 call bdry_mask%init_base(this%c_Xh)
1100 call bdry_mask%mark_zone(this%msh%outlet)
1101 call bdry_mask%mark_zones_from_list(this%msh%labeled_zones, &
1102 'o', this%bc_labels)
1103 call bdry_mask%mark_zones_from_list(this%msh%labeled_zones, &
1104 'o+dong', this%bc_labels)
1105 call bdry_mask%finalize()
1106 call bdry_mask%set_g(3.0_rp)
1107 call bdry_mask%apply_scalar(this%bdry%x, this%dm_Xh%size())
1108 call bdry_mask%free()
1110 call bdry_mask%init_base(this%c_Xh)
1111 call bdry_mask%mark_zone(this%msh%sympln)
1112 call bdry_mask%mark_zones_from_list(this%msh%labeled_zones, &
1113 'sym', this%bc_labels)
1114 call bdry_mask%finalize()
1115 call bdry_mask%set_g(4.0_rp)
1116 call bdry_mask%apply_scalar(this%bdry%x, this%dm_Xh%size())
1117 call bdry_mask%free()
1119 call bdry_mask%init_base(this%c_Xh)
1120 call bdry_mask%mark_zone(this%msh%outlet_normal)
1121 call bdry_mask%mark_zones_from_list(this%msh%labeled_zones, &
1122 'on', this%bc_labels)
1123 call bdry_mask%mark_zones_from_list(this%msh%labeled_zones, &
1124 'on+dong', this%bc_labels)
1125 call bdry_mask%finalize()
1126 call bdry_mask%set_g(5.0_rp)
1127 call bdry_mask%apply_scalar(this%bdry%x, this%dm_Xh%size())
1128 call bdry_mask%free()
1130 call bdry_mask%init_base(this%c_Xh)
1131 call bdry_mask%mark_zone(this%msh%periodic)
1132 call bdry_mask%finalize()
1133 call bdry_mask%set_g(6.0_rp)
1134 call bdry_mask%apply_scalar(this%bdry%x, this%dm_Xh%size())
1135 call bdry_mask%free()
1137 call bdry_mask%init_base(this%c_Xh)
1138 call bdry_mask%mark_zones_from_list(this%msh%labeled_zones, &
1139 'd_vel_u', this%bc_labels)
1140 call bdry_mask%mark_zones_from_list(this%msh%labeled_zones, &
1141 'd_vel_v', this%bc_labels)
1142 call bdry_mask%mark_zones_from_list(this%msh%labeled_zones, &
1143 'd_vel_w', this%bc_labels)
1144 call bdry_mask%finalize()
1145 call bdry_mask%set_g(7.0_rp)
1146 call bdry_mask%apply_scalar(this%bdry%x, this%dm_Xh%size())
1147 call bdry_mask%free()
1149 call bdry_mask%init_base(this%c_Xh)
1150 call bdry_mask%mark_zones_from_list(this%msh%labeled_zones, &
1151 'd_pres', this%bc_labels)
1152 call bdry_mask%finalize()
1153 call bdry_mask%set_g(8.0_rp)
1154 call bdry_mask%apply_scalar(this%bdry%x, this%dm_Xh%size())
1155 call bdry_mask%free()
1157 call bdry_mask%init_base(this%c_Xh)
1158 call bdry_mask%mark_zones_from_list(this%msh%labeled_zones, &
1159 'sh', this%bc_labels)
1160 call bdry_mask%finalize()
1161 call bdry_mask%set_g(9.0_rp)
1162 call bdry_mask%apply_scalar(this%bdry%x, this%dm_Xh%size())
1163 call bdry_mask%free()
1165 call bdry_mask%init_base(this%c_Xh)
1166 call bdry_mask%mark_zones_from_list(this%msh%labeled_zones, &
1167 'wm', this%bc_labels)
1168 call bdry_mask%finalize()
1169 call bdry_mask%set_g(10.0_rp)
1170 call bdry_mask%apply_scalar(this%bdry%x, this%dm_Xh%size())
1171 call bdry_mask%free()
1175 end subroutine fluid_scheme_set_bc_type_output
1178 subroutine fluid_scheme_update_material_properties(this)
1179 class(fluid_scheme_t),
intent(inout) :: this
1183 if (this%variable_material_properties)
then
1191 call cfill(this%mu_field%x, this%mu, n)
1192 call add2s2(this%mu_field%x, nut%x, this%rho, n)
1196 end subroutine fluid_scheme_update_material_properties
1201 subroutine fluid_scheme_set_material_properties(this, params, user)
1202 class(fluid_scheme_t),
intent(inout) :: this
1203 type(json_file),
intent(inout) :: params
1204 type(
user_t),
target,
intent(in) :: user
1205 character(len=LOG_SIZE) :: log_buf
1208 logical :: nondimensional
1209 real(kind=
rp) :: dummy_lambda, dummy_cp
1213 if (.not.
associated(user%material_properties, dummy_mp_ptr))
then
1215 write(log_buf,
'(A)')
"Material properties must be set in the user&
1218 call user%material_properties(0.0_rp, 0, this%rho, this%mu, &
1219 dummy_cp, dummy_lambda, params)
1222 if (params%valid_path(
'case.fluid.Re') .and. &
1223 (params%valid_path(
'case.fluid.mu') .or. &
1224 params%valid_path(
'case.fluid.rho')))
then
1225 call neko_error(
"To set the material properties for the fluid,&
1226 & either provide Re OR mu and rho in the case file.")
1229 else if (params%valid_path(
'case.fluid.Re'))
then
1231 write(log_buf,
'(A)')
'Non-dimensional fluid material properties &
1234 write(log_buf,
'(A)')
'Density will be set to 1, dynamic viscosity to&
1239 call json_get(params,
'case.fluid.Re', this%mu)
1240 write(log_buf,
'(A)')
'Read non-dimensional material properties'
1242 write(log_buf,
'(A,ES13.6)')
'Re :', this%mu
1248 this%mu = 1.0_rp/this%mu
1251 call json_get(params,
'case.fluid.mu', this%mu)
1252 call json_get(params,
'case.fluid.rho', this%rho)
1256 end subroutine fluid_scheme_set_material_properties
Abstract interface to dealocate a fluid formulation.
Abstract interface to initialize a fluid formulation.
Abstract interface to restart a fluid scheme.
Abstract interface to compute a time-step.
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.
Abstract interface defining a user defined inflow condition (pointwise)
Defines a boundary condition.
subroutine, public bc_list_add(bclst, bc)
Add a condition to a list of boundary conditions.
subroutine, public bc_list_apply_vector(bclst, x, y, z, n, t, tstep)
Apply a list of boundary conditions to a vector field.
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.
Defines a Blasius profile dirichlet condition.
type(mpi_comm) neko_comm
MPI communicator.
Jacobi preconditioner accelerator backend.
subroutine, public device_add2s2(a_d, b_d, c1, n)
Vector addition with scalar multiplication (multiplication on first argument)
subroutine, public device_cfill(a_d, c, n)
Set all elements to a constant c .
Defines a dirichlet boundary condition.
Defines a mapping of the degrees of freedom.
Defines a dong outflow condition.
Defines inflow dirichlet conditions.
Defines inflow dirichlet conditions.
subroutine, public field_cfill(a, c, n)
Set all elements to a constant c .
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_free(this)
Deallocate a fluid formulation.
subroutine fluid_scheme_update_material_properties(this)
Update the values of mu_field if necessary.
subroutine fluid_scheme_set_bc_type_output(this, params)
Set boundary types for the diagnostic output.
real(kind=rp) function fluid_compute_cfl(this, dt)
Compute CFL.
subroutine fluid_scheme_init_all(this, msh, lx, params, kspv_init, kspp_init, scheme, user)
Initialize all components of the current scheme.
subroutine fluid_scheme_set_usr_inflow(this, usr_eval)
Initialize a user defined inflow condition.
subroutine fluid_scheme_init_common(this, msh, lx, params, scheme, user, kspv_init)
Initialise a fluid scheme.
subroutine fluid_scheme_bc_apply_vel(this, t, tstep)
Apply all boundary conditions defined for velocity Here we perform additional gs operations to take c...
subroutine fluid_scheme_validate(this)
Validate that all fields, solvers etc necessary for performing time-stepping are defined.
subroutine fluid_scheme_bc_apply_prs(this, t, tstep)
Apply all boundary conditions defined for pressure.
Implements the fluid_source_term_t type.
Computes various statistics for the fluid fields. We use the Reynolds decomposition for a field u = ...
Implements gradient_jump_penalty_t.
Defines inflow dirichlet conditions.
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
subroutine, public cfill(a, c, n)
Set all elements to a constant c .
subroutine, public add2s2(a, b, c1, n)
Vector addition with scalar multiplication (multiplication on second argument)
Defines a mean flow field.
Defines a mean squared flow field.
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.
integer, parameter neko_bcknd_device
integer, parameter, public rp
Global precision used in computations.
real(kind=rp) function, public cfl(dt, u, v, w, Xh, coef, nelv, gdim)
Defines a registry for storing and requesting temporary fields This can be used when you have a funct...
Defines a shear stress boundary condition for a vector field. Maintainer: Timofey Mukha.
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.
Mixed Dirichlet-Neumann axis aligned symmetry plane.
Compound scheme for the advection and diffusion operators in a transport equation.
Implements type time_step_controller.
Interfaces for user interaction with NEKO.
subroutine, public dummy_user_material_properties(t, tstep, rho, mu, cp, lambda, params)
Defines inflow dirichlet conditions.
subroutine, public neko_warning(warning_msg)
Reports a warning to standard output.
Defines the wall_model_bc_t type. Maintainer: Timofey Mukha.
Defines wall boundary conditions.
A list of boundary conditions.
Base type for a boundary condition.
Blasius profile for inlet (vector valued)
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 .
Dong outflow condition Follows "A Convective-like Energy-Stable Open Boundary Condition for Simulati...
User defined dirichlet condition, for which the user can work with an entire field....
Extension of the user defined dirichlet condition field_dirichlet
Base type of all fluid formulations.
Wrapper contaning and executing the fluid source terms.
Implements the gradient jump penalty.
Dirichlet condition for inlet (vector valued)
Defines a jacobi preconditioner.
Base abstract type for a canonical Krylov method, solving .
Defines a canonical Krylov preconditioner.
A shear stress boundary condition.
The function space for the SEM solution fields.
Defines a jacobi preconditioner for SX-Aurora.
Mixed Dirichlet-Neumann symmetry plane condition.
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 inlet (vector valued)
No-slip Wall boundary condition.
A shear stress boundary condition, computing the stress values using a wall model.