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
104 class(
ksp_t),
allocatable :: ksp_vel
105 class(
ksp_t),
allocatable :: ksp_prs
106 class(
pc_t),
allocatable :: pc_vel
107 class(
pc_t),
allocatable :: pc_prs
108 integer :: vel_projection_dim
109 integer :: pr_projection_dim
110 integer :: vel_projection_activ_step
111 integer :: pr_projection_activ_step
113 class(
bc_t),
allocatable :: bc_inflow
124 type(json_file),
pointer :: params
130 logical :: forced_flow_rate = .false.
131 logical :: freeze = .false.
133 real(kind=
rp),
pointer :: mu => null()
137 character(len=:),
allocatable :: nut_field_name
139 logical :: variable_material_properties = .false.
141 real(kind=
rp),
pointer :: rho => null()
146 character(len=NEKO_MSH_MAX_ZLBL_LEN),
allocatable :: bc_labels(:)
172 procedure,
private, pass(this) :: set_bc_type_output => &
175 procedure, pass(this) :: update_material_properties => &
182 material_properties, time_scheme)
190 type(
mesh_t),
target,
intent(inout) :: msh
191 integer,
intent(inout) :: lx
192 type(json_file),
target,
intent(inout) :: params
193 type(
user_t),
intent(in) :: user
195 target,
intent(inout) :: material_properties
217 real(kind=
rp),
intent(inout) :: t
218 integer,
intent(inout) :: tstep
219 real(kind=
rp),
intent(in) :: dt
231 real(kind=
rp) :: dtlag(10), tlag(10)
238 module subroutine fluid_scheme_factory(object, type_name)
239 class(fluid_scheme_t),
intent(inout),
allocatable :: object
240 character(len=*) :: type_name
241 end subroutine fluid_scheme_factory
244 public :: fluid_scheme_t, fluid_scheme_factory
249 subroutine fluid_scheme_init_common(this, msh, lx, params, scheme, user, &
250 material_properties, kspv_init)
252 class(fluid_scheme_t),
target,
intent(inout) :: this
253 type(
mesh_t),
target,
intent(inout) :: msh
254 integer,
intent(inout) :: lx
255 character(len=*),
intent(in) :: scheme
256 type(json_file),
target,
intent(inout) :: params
257 type(
user_t),
target,
intent(in) :: user
259 logical,
intent(in) :: kspv_init
261 character(len=LOG_SIZE) :: log_buf
262 real(kind=
rp),
allocatable :: real_vec(:)
263 real(kind=
rp) :: real_val
264 logical :: logical_val
265 integer :: integer_val, ierr
266 character(len=:),
allocatable :: string_val1, string_val2
275 if (msh%gdim .eq. 2)
then
276 call this%Xh%init(
gll, lx, lx)
278 call this%Xh%init(
gll, lx, lx, lx)
283 call this%gs_Xh%init(this%dm_Xh)
285 call this%c_Xh%init(this%gs_Xh)
291 this%params => params
303 if (params%valid_path(
'case.fluid.nut_field'))
then
304 call json_get(params,
'case.fluid.nut_field', this%nut_field_name)
305 this%variable_material_properties = .true.
307 this%nut_field_name =
""
312 call this%mu_field%init(this%dm_Xh,
"mu")
313 call this%rho_field%init(this%dm_Xh,
"mu")
314 call field_cfill(this%mu_field, this%mu, this%mu_field%size())
315 call field_cfill(this%rho_field, this%rho, this%mu_field%size())
321 call cfill(this%mu_field%x, this%mu, this%mu_field%size())
322 call cfill(this%rho_field%x, this%rho, this%rho_field%size())
328 'case.fluid.velocity_solver.projection_space_size', &
329 this%vel_projection_dim, 20)
331 'case.fluid.pressure_solver.projection_space_size', &
332 this%pr_projection_dim, 20)
334 'case.fluid.velocity_solver.projection_hold_steps', &
335 this%vel_projection_activ_step, 5)
337 'case.fluid.pressure_solver.projection_hold_steps', &
338 this%pr_projection_activ_step, 5)
343 if (params%valid_path(
"case.fluid.flow_rate_force"))
then
344 this%forced_flow_rate = .true.
353 write(log_buf,
'(A, A)')
'Type : ', trim(scheme)
356 write(log_buf,
'(A, I1)')
'Poly order : ', lx-1
357 else if (lx .ge. 10)
then
358 write(log_buf,
'(A, I2)')
'Poly order : ', lx-1
360 write(log_buf,
'(A, I3)')
'Poly order : ', lx-1
363 write(log_buf,
'(A, I0)')
'DoFs : ', this%dm_Xh%size()
366 write(log_buf,
'(A,ES13.6)')
'rho :', this%rho
368 write(log_buf,
'(A,ES13.6)')
'mu :', this%mu
371 call json_get(params,
'case.numerics.dealias', logical_val)
372 write(log_buf,
'(A, L1)')
'Dealias : ', logical_val
377 write(log_buf,
'(A, L1)')
'Save bdry : ', logical_val
379 if (logical_val)
then
380 write(log_buf,
'(A)')
'bdry keys: '
382 write(log_buf,
'(A)')
'No-slip wall ''w'' = 1'
384 write(log_buf,
'(A)')
'Inlet/velocity dirichlet ''v'' = 2'
386 write(log_buf,
'(A)')
'Outlet ''o'' = 3'
388 write(log_buf,
'(A)')
'Symmetry ''sym'' = 4'
390 write(log_buf,
'(A)')
'Outlet-normal ''on'' = 5'
392 write(log_buf,
'(A)')
'Periodic = 6'
394 write(log_buf,
'(A)')
'Dirichlet on velocity components: '
396 write(log_buf,
'(A)')
' ''d_vel_u, d_vel_v, d_vel_w'' = 7'
398 write(log_buf,
'(A)')
'Pressure dirichlet ''d_pres'' = 8'
400 write(log_buf,
'(A)')
'''d_vel_(u,v,w)'' and ''d_pres'' = 8'
402 write(log_buf,
'(A)')
'No boundary condition set = 0'
411 this%bc_labels =
"not"
413 if (params%valid_path(
'case.fluid.boundary_types'))
then
415 'case.fluid.boundary_types', &
421 call this%bc_sym%init_base(this%c_Xh)
422 call this%bc_sym%mark_zone(msh%sympln)
423 call this%bc_sym%mark_zones_from_list(msh%labeled_zones, &
424 'sym', this%bc_labels)
425 call this%bc_sym%finalize()
426 call this%bc_sym%init(this%c_Xh)
432 if (params%valid_path(
'case.fluid.inflow_condition'))
then
433 call json_get(params,
'case.fluid.inflow_condition.type', string_val1)
434 if (trim(string_val1) .eq.
"uniform")
then
436 else if (trim(string_val1) .eq.
"blasius")
then
438 else if (trim(string_val1) .eq.
"user")
then
441 call neko_error(
'Invalid inflow condition '//string_val1)
444 call this%bc_inflow%init_base(this%c_Xh)
445 call this%bc_inflow%mark_zone(msh%inlet)
446 call this%bc_inflow%mark_zones_from_list(msh%labeled_zones, &
448 call this%bc_inflow%finalize()
451 if (trim(string_val1) .eq.
"uniform")
then
452 call json_get(params,
'case.fluid.inflow_condition.value', real_vec)
453 select type (bc_if => this%bc_inflow)
455 call bc_if%set_inflow(real_vec)
457 else if (trim(string_val1) .eq.
"blasius")
then
458 select type (bc_if => this%bc_inflow)
460 call json_get(params,
'case.fluid.blasius.delta', real_val)
461 call json_get(params,
'case.fluid.blasius.approximation', &
463 call json_get(params,
'case.fluid.blasius.freestream_velocity', &
466 call bc_if%set_params(real_vec, real_val, string_val2)
469 else if (trim(string_val1) .eq.
"user")
then
473 call this%bc_wall%init_base(this%c_Xh)
474 call this%bc_wall%mark_zone(msh%wall)
475 call this%bc_wall%mark_zones_from_list(msh%labeled_zones, &
477 call this%bc_wall%finalize()
481 call this%user_field_bc_vel%bc_u%init_base(this%c_Xh)
482 call this%user_field_bc_vel%bc_u%mark_zones_from_list(msh%labeled_zones, &
483 'd_vel_u', this%bc_labels)
484 call this%user_field_bc_vel%bc_u%finalize()
486 call mpi_allreduce(this%user_field_bc_vel%bc_u%msk(0), integer_val, 1, &
488 if (integer_val .gt. 0)
then
489 call this%user_field_bc_vel%bc_u%init_field(
'd_vel_u')
493 call this%user_field_bc_vel%bc_v%init_base(this%c_Xh)
494 call this%user_field_bc_vel%bc_v%mark_zones_from_list(msh%labeled_zones, &
495 'd_vel_v', this%bc_labels)
496 call this%user_field_bc_vel%bc_v%finalize()
498 call mpi_allreduce(this%user_field_bc_vel%bc_v%msk(0), integer_val, 1, &
500 if (integer_val .gt. 0)
then
501 call this%user_field_bc_vel%bc_v%init_field(
'd_vel_v')
505 call this%user_field_bc_vel%bc_w%init_base(this%c_Xh)
506 call this%user_field_bc_vel%bc_w%mark_zones_from_list(msh%labeled_zones, &
507 'd_vel_w', this%bc_labels)
508 call this%user_field_bc_vel%bc_w%finalize()
510 call mpi_allreduce(this%user_field_bc_vel%bc_w%msk(0), integer_val, 1, &
512 if (integer_val .gt. 0)
then
513 call this%user_field_bc_vel%bc_w%init_field(
'd_vel_w')
517 call this%user_field_bc_vel%init_base(this%c_Xh)
518 call this%user_field_bc_vel%mark_zones_from_list(msh%labeled_zones, &
519 'd_vel_u', this%bc_labels)
520 call this%user_field_bc_vel%mark_zones_from_list(msh%labeled_zones, &
521 'd_vel_v', this%bc_labels)
522 call this%user_field_bc_vel%mark_zones_from_list(msh%labeled_zones, &
523 'd_vel_w', this%bc_labels)
524 call this%user_field_bc_vel%finalize()
527 call bc_list_add(this%bclst_vel, this%user_field_bc_vel)
532 this%user_field_bc_vel%update => user%user_dirichlet_update
539 call this%user_field_bc_vel%field_list%init(4)
540 call this%user_field_bc_vel%field_list%assign_to_field(1, &
541 this%user_field_bc_vel%bc_u%field_bc)
542 call this%user_field_bc_vel%field_list%assign_to_field(2, &
543 this%user_field_bc_vel%bc_v%field_bc)
544 call this%user_field_bc_vel%field_list%assign_to_field(3, &
545 this%user_field_bc_vel%bc_w%field_bc)
546 call this%user_field_bc_vel%field_list%assign_to_field(4, &
547 this%user_field_bc_prs%field_bc)
549 call bc_list_init(this%user_field_bc_vel%bc_list,
size = 4)
552 this%user_field_bc_vel%bc_u)
554 this%user_field_bc_vel%bc_v)
556 this%user_field_bc_vel%bc_w)
560 call fluid_scheme_set_bc_type_output(this, params)
568 call this%f_x%init(this%dm_Xh, fld_name =
"fluid_rhs_x")
569 call this%f_y%init(this%dm_Xh, fld_name =
"fluid_rhs_y")
570 call this%f_z%init(this%dm_Xh, fld_name =
"fluid_rhs_z")
573 call this%source_term%init(this%f_x, this%f_y, this%f_z, this%c_Xh, user)
574 call this%source_term%add(params,
'case.fluid.source_term')
578 call this%set_bc_type_output(params)
582 call neko_log%section(
"Velocity solver")
584 'case.fluid.velocity_solver.max_iterations', &
586 call json_get(params,
'case.fluid.velocity_solver.type', string_val1)
587 call json_get(params,
'case.fluid.velocity_solver.preconditioner', &
589 call json_get(params,
'case.fluid.velocity_solver.absolute_tolerance', &
592 'case.fluid.velocity_solver.monitor', &
593 logical_val, .false.)
595 call neko_log%message(
'Type : ('// trim(string_val1) // &
596 ', ' // trim(string_val2) //
')')
598 write(log_buf,
'(A,ES13.6)')
'Abs tol :', real_val
600 call fluid_scheme_solver_factory(this%ksp_vel, this%dm_Xh%size(), &
601 string_val1, integer_val, real_val, logical_val)
602 call fluid_scheme_precon_factory(this%pc_vel, this%ksp_vel, &
603 this%c_Xh, this%dm_Xh, this%gs_Xh, this%bclst_vel, string_val2)
616 call this%ulag%init(this%u, 2)
617 call this%vlag%init(this%v, 2)
618 call this%wlag%init(this%w, 2)
621 end subroutine fluid_scheme_init_common
624 subroutine fluid_scheme_init_all(this, msh, lx, params, kspv_init, &
625 kspp_init, scheme, user, &
628 class(fluid_scheme_t),
target,
intent(inout) :: this
629 type(
mesh_t),
target,
intent(inout) :: msh
630 integer,
intent(inout) :: lx
631 type(json_file),
target,
intent(inout) :: params
632 type(
user_t),
target,
intent(in) :: user
636 character(len=*),
intent(in) :: scheme
637 real(kind=
rp) :: abs_tol
638 integer :: integer_val, ierr
639 logical :: logical_val
640 character(len=:),
allocatable :: solver_type, precon_type
641 character(len=LOG_SIZE) :: log_buf
643 call fluid_scheme_init_common(this, msh, lx, params, scheme, user, &
653 call this%bc_prs%init_base(this%c_Xh)
654 call this%bc_prs%mark_zones_from_list(msh%labeled_zones, &
656 call this%bc_prs%mark_zones_from_list(msh%labeled_zones, &
657 'on', this%bc_labels)
660 call this%user_field_bc_prs%init_base(this%c_Xh)
661 call this%user_field_bc_prs%mark_zones_from_list(msh%labeled_zones, &
662 'd_pres', this%bc_labels)
663 call this%user_field_bc_prs%finalize()
664 call mpi_allreduce(this%user_field_bc_prs%msk(0), integer_val, 1, &
667 if (integer_val .gt. 0)
call this%user_field_bc_prs%init_field(
'd_pres')
668 call bc_list_add(this%bclst_prs, this%user_field_bc_prs)
669 call bc_list_add(this%user_field_bc_vel%bc_list, this%user_field_bc_prs)
671 if (msh%outlet%size .gt. 0)
then
672 call this%bc_prs%mark_zone(msh%outlet)
674 if (msh%outlet_normal%size .gt. 0)
then
675 call this%bc_prs%mark_zone(msh%outlet_normal)
678 call this%bc_prs%finalize()
679 call this%bc_prs%set_g(0.0_rp)
681 call this%bc_dong%init_base(this%c_Xh)
682 call this%bc_dong%mark_zones_from_list(msh%labeled_zones, &
683 'o+dong', this%bc_labels)
684 call this%bc_dong%mark_zones_from_list(msh%labeled_zones, &
685 'on+dong', this%bc_labels)
686 call this%bc_dong%finalize()
689 call this%bc_dong%init(this%c_Xh, params)
695 call neko_log%section(
"Pressure solver")
698 'case.fluid.pressure_solver.max_iterations', &
700 call json_get(params,
'case.fluid.pressure_solver.type', solver_type)
701 call json_get(params,
'case.fluid.pressure_solver.preconditioner', &
703 call json_get(params,
'case.fluid.pressure_solver.absolute_tolerance', &
706 'case.fluid.pressure_solver.monitor', &
707 logical_val, .false.)
708 call neko_log%message(
'Type : ('// trim(solver_type) // &
709 ', ' // trim(precon_type) //
')')
710 write(log_buf,
'(A,ES13.6)')
'Abs tol :', abs_tol
713 call fluid_scheme_solver_factory(this%ksp_prs, this%dm_Xh%size(), &
714 solver_type, integer_val, abs_tol, logical_val)
715 call fluid_scheme_precon_factory(this%pc_prs, this%ksp_prs, &
716 this%c_Xh, this%dm_Xh, this%gs_Xh, this%bclst_prs, precon_type)
725 end subroutine fluid_scheme_init_all
728 subroutine fluid_scheme_free(this)
729 class(fluid_scheme_t),
intent(inout) :: this
731 call this%bdry%free()
733 if (
allocated(this%bc_inflow))
then
734 call this%bc_inflow%free()
737 call this%bc_wall%free()
738 call this%bc_sym%free()
743 call this%user_field_bc_prs%field_bc%free()
744 call this%user_field_bc_prs%free()
745 call this%user_field_bc_vel%bc_u%field_bc%free()
746 call this%user_field_bc_vel%bc_v%field_bc%free()
747 call this%user_field_bc_vel%bc_w%field_bc%free()
748 call this%user_field_bc_vel%free()
752 if (
allocated(this%ksp_vel))
then
753 call krylov_solver_destroy(this%ksp_vel)
754 deallocate(this%ksp_vel)
757 if (
allocated(this%ksp_prs))
then
758 call krylov_solver_destroy(this%ksp_prs)
759 deallocate(this%ksp_prs)
762 if (
allocated(this%pc_vel))
then
763 call precon_destroy(this%pc_vel)
764 deallocate(this%pc_vel)
767 if (
allocated(this%pc_prs))
then
768 call precon_destroy(this%pc_prs)
769 deallocate(this%pc_prs)
772 if (
allocated(this%bc_labels))
then
773 deallocate(this%bc_labels)
776 call this%source_term%free()
778 call this%gs_Xh%free()
780 call this%c_Xh%free()
784 call this%scratch%free()
793 call this%ulag%free()
794 call this%vlag%free()
795 call this%wlag%free()
798 if (
associated(this%f_x))
then
802 if (
associated(this%f_y))
then
806 if (
associated(this%f_z))
then
814 call this%mu_field%free()
817 end subroutine fluid_scheme_free
821 subroutine fluid_scheme_validate(this)
822 class(fluid_scheme_t),
target,
intent(inout) :: this
824 logical :: logical_val
826 if ( (.not.
associated(this%u)) .or. &
827 (.not.
associated(this%v)) .or. &
828 (.not.
associated(this%w)) .or. &
829 (.not.
associated(this%p)))
then
833 if ( (.not.
allocated(this%u%x)) .or. &
834 (.not.
allocated(this%v%x)) .or. &
835 (.not.
allocated(this%w%x)) .or. &
836 (.not.
allocated(this%p%x)))
then
840 if (.not.
allocated(this%ksp_vel))
then
841 call neko_error(
'No Krylov solver for velocity defined')
844 if (.not.
allocated(this%ksp_prs))
then
845 call neko_error(
'No Krylov solver for pressure defined')
848 if (.not.
associated(this%params))
then
852 if (
allocated(this%bc_inflow))
then
853 select type (ip => this%bc_inflow)
862 call this%chkp%init(this%u, this%v, this%w, this%p)
867 if (this%params%valid_path(
'case.statistics'))
then
870 if (logical_val)
then
871 call this%mean%init(this%u, this%v, this%w, this%p)
872 call this%stats%init(this%c_Xh, this%mean%u, &
873 this%mean%v, this%mean%w, this%mean%p)
877 end subroutine fluid_scheme_validate
883 subroutine fluid_scheme_bc_apply_vel(this, t, tstep)
884 class(fluid_scheme_t),
intent(inout) :: this
885 real(kind=
rp),
intent(in) :: t
886 integer,
intent(in) :: tstep
889 this%u%x, this%v%x, this%w%x, this%dm_Xh%size(), t, tstep)
891 end subroutine fluid_scheme_bc_apply_vel
895 subroutine fluid_scheme_bc_apply_prs(this, t, tstep)
896 class(fluid_scheme_t),
intent(inout) :: this
897 real(kind=
rp),
intent(in) :: t
898 integer,
intent(in) :: tstep
901 this%p%dof%size(), t, tstep)
903 end subroutine fluid_scheme_bc_apply_prs
907 subroutine fluid_scheme_solver_factory(ksp, n, solver, &
908 max_iter, abstol, monitor)
909 class(
ksp_t),
allocatable,
target,
intent(inout) :: ksp
910 integer,
intent(in),
value :: n
911 character(len=*),
intent(in) :: solver
912 integer,
intent(in) :: max_iter
913 real(kind=
rp),
intent(in) :: abstol
914 logical,
intent(in) :: monitor
916 call krylov_solver_factory(ksp, n, solver, max_iter, abstol, &
919 end subroutine fluid_scheme_solver_factory
922 subroutine fluid_scheme_precon_factory(pc, ksp, coef, dof, gs, bclst, &
924 class(
pc_t),
allocatable,
target,
intent(inout) :: pc
925 class(
ksp_t),
target,
intent(inout) :: ksp
926 type(
coef_t),
target,
intent(inout) :: coef
927 type(
dofmap_t),
target,
intent(inout) :: dof
928 type(
gs_t),
target,
intent(inout) :: gs
929 type(
bc_list_t),
target,
intent(inout) :: bclst
930 character(len=*) :: pctype
932 call precon_factory(pc, pctype)
934 select type (pcp => pc)
936 call pcp%init(coef, dof, gs)
938 call pcp%init(coef, dof, gs)
940 call pcp%init(coef, dof, gs)
942 if (len_trim(pctype) .gt. 4)
then
943 if (index(pctype,
'+') .eq. 5)
then
944 call pcp%init(dof%msh, dof%Xh, coef, dof, gs, bclst, &
950 call pcp%init(dof%msh, dof%Xh, coef, dof, gs, bclst)
956 end subroutine fluid_scheme_precon_factory
959 subroutine fluid_scheme_set_usr_inflow(this, usr_eval)
960 class(fluid_scheme_t),
intent(inout) :: this
963 select type (bc_if => this%bc_inflow)
965 call bc_if%set_eval(usr_eval)
967 call neko_error(
"Not a user defined inflow condition")
969 end subroutine fluid_scheme_set_usr_inflow
972 function fluid_compute_cfl(this, dt)
result(c)
973 class(fluid_scheme_t),
intent(in) :: this
974 real(kind=
rp),
intent(in) :: dt
977 c =
cfl(dt, this%u%x, this%v%x, this%w%x, &
978 this%Xh, this%c_Xh, this%msh%nelv, this%msh%gdim)
980 end function fluid_compute_cfl
984 subroutine fluid_scheme_set_bc_type_output(this, params)
985 class(fluid_scheme_t),
target,
intent(inout) :: this
986 type(json_file),
intent(inout) :: params
996 call this%bdry%init(this%dm_Xh,
'bdry')
999 call bdry_mask%init_base(this%c_Xh)
1000 call bdry_mask%mark_zone(this%msh%wall)
1001 call bdry_mask%mark_zones_from_list(this%msh%labeled_zones, &
1002 'w', this%bc_labels)
1003 call bdry_mask%finalize()
1004 call bdry_mask%set_g(1.0_rp)
1005 call bdry_mask%apply_scalar(this%bdry%x, this%dm_Xh%size())
1006 call bdry_mask%free()
1008 call bdry_mask%init_base(this%c_Xh)
1009 call bdry_mask%mark_zone(this%msh%inlet)
1010 call bdry_mask%mark_zones_from_list(this%msh%labeled_zones, &
1011 'v', this%bc_labels)
1012 call bdry_mask%finalize()
1013 call bdry_mask%set_g(2.0_rp)
1014 call bdry_mask%apply_scalar(this%bdry%x, this%dm_Xh%size())
1015 call bdry_mask%free()
1017 call bdry_mask%init_base(this%c_Xh)
1018 call bdry_mask%mark_zone(this%msh%outlet)
1019 call bdry_mask%mark_zones_from_list(this%msh%labeled_zones, &
1020 'o', this%bc_labels)
1021 call bdry_mask%mark_zones_from_list(this%msh%labeled_zones, &
1022 'o+dong', this%bc_labels)
1023 call bdry_mask%finalize()
1024 call bdry_mask%set_g(3.0_rp)
1025 call bdry_mask%apply_scalar(this%bdry%x, this%dm_Xh%size())
1026 call bdry_mask%free()
1028 call bdry_mask%init_base(this%c_Xh)
1029 call bdry_mask%mark_zone(this%msh%sympln)
1030 call bdry_mask%mark_zones_from_list(this%msh%labeled_zones, &
1031 'sym', this%bc_labels)
1032 call bdry_mask%finalize()
1033 call bdry_mask%set_g(4.0_rp)
1034 call bdry_mask%apply_scalar(this%bdry%x, this%dm_Xh%size())
1035 call bdry_mask%free()
1037 call bdry_mask%init_base(this%c_Xh)
1038 call bdry_mask%mark_zone(this%msh%outlet_normal)
1039 call bdry_mask%mark_zones_from_list(this%msh%labeled_zones, &
1040 'on', this%bc_labels)
1041 call bdry_mask%mark_zones_from_list(this%msh%labeled_zones, &
1042 'on+dong', this%bc_labels)
1043 call bdry_mask%finalize()
1044 call bdry_mask%set_g(5.0_rp)
1045 call bdry_mask%apply_scalar(this%bdry%x, this%dm_Xh%size())
1046 call bdry_mask%free()
1048 call bdry_mask%init_base(this%c_Xh)
1049 call bdry_mask%mark_zone(this%msh%periodic)
1050 call bdry_mask%finalize()
1051 call bdry_mask%set_g(6.0_rp)
1052 call bdry_mask%apply_scalar(this%bdry%x, this%dm_Xh%size())
1053 call bdry_mask%free()
1055 call bdry_mask%init_base(this%c_Xh)
1056 call bdry_mask%mark_zones_from_list(this%msh%labeled_zones, &
1057 'd_vel_u', this%bc_labels)
1058 call bdry_mask%mark_zones_from_list(this%msh%labeled_zones, &
1059 'd_vel_v', this%bc_labels)
1060 call bdry_mask%mark_zones_from_list(this%msh%labeled_zones, &
1061 'd_vel_w', this%bc_labels)
1062 call bdry_mask%finalize()
1063 call bdry_mask%set_g(7.0_rp)
1064 call bdry_mask%apply_scalar(this%bdry%x, this%dm_Xh%size())
1065 call bdry_mask%free()
1067 call bdry_mask%init_base(this%c_Xh)
1068 call bdry_mask%mark_zones_from_list(this%msh%labeled_zones, &
1069 'd_pres', this%bc_labels)
1070 call bdry_mask%finalize()
1071 call bdry_mask%set_g(8.0_rp)
1072 call bdry_mask%apply_scalar(this%bdry%x, this%dm_Xh%size())
1073 call bdry_mask%free()
1077 end subroutine fluid_scheme_set_bc_type_output
1080 subroutine fluid_scheme_update_material_properties(this)
1081 class(fluid_scheme_t),
intent(inout) :: this
1085 if (this%variable_material_properties)
then
1093 call cfill(this%mu_field%x, this%mu, n)
1094 call add2s2(this%mu_field%x, nut%x, this%rho, n)
1098 end subroutine fluid_scheme_update_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 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_init_all(this, msh, lx, params, kspv_init, kspp_init, scheme, user, material_properties)
Initialize all components of the current scheme.
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_set_usr_inflow(this, usr_eval)
Initialize a user defined inflow condition.
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_init_common(this, msh, lx, params, scheme, user, material_properties, kspv_init)
Initialise a fluid scheme.
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 = ...
Defines inflow dirichlet conditions.
Utilities for retrieving parameters from the case files.
Implements the base abstract type for Krylov solvers plus helper types.
integer, parameter, public ksp_max_iter
Maximum number of iters.
type(log_t), public neko_log
Global log stream.
integer, parameter, public log_size
Implements material_properties_t type.
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...
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.
Defines inflow dirichlet conditions.
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.
Dirichlet condition for inlet (vector valued)
Defines a jacobi preconditioner.
Base abstract type for a canonical Krylov method, solving .
Contains all the material properties necessary in the simulation.
Defines a canonical Krylov preconditioner.
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.