245 subroutine fluid_pnpn_init(this, msh, lx, params, user, chkp)
246 class(fluid_pnpn_t),
target,
intent(inout) :: this
247 type(
mesh_t),
target,
intent(inout) :: msh
248 integer,
intent(in) :: lx
249 type(json_file),
target,
intent(inout) :: params
250 type(
user_t),
target,
intent(in) :: user
251 type(
chkp_t),
target,
intent(inout) :: chkp
252 character(len=15),
parameter :: scheme =
'Modular (Pn/Pn)'
254 class(
bc_t),
pointer :: bc_i, vel_bc
255 real(kind=
rp) :: abs_tol
256 character(len=LOG_SIZE) :: log_buf
257 integer :: ierr, integer_val, solver_maxiter
258 character(len=:),
allocatable :: solver_type, precon_type
259 logical :: monitor, found
261 type(json_file) :: numerics_params, precon_params
266 call this%init_base(msh, lx, params, scheme,
user, .true.)
277 call json_get(params,
'case.numerics.time_order', integer_val)
278 allocate(this%ext_bdf)
279 call this%ext_bdf%init(integer_val)
282 this%full_stress_formulation, .false.)
284 if (this%full_stress_formulation .eqv. .true.)
then
286 call ax_helm_factory(this%Ax_vel, full_formulation = .true.)
289 call pnpn_prs_res_stress_factory(this%prs_res)
292 call pnpn_vel_res_stress_factory(this%vel_res)
295 call ax_helm_factory(this%Ax_vel, full_formulation = .false.)
298 call pnpn_prs_res_factory(this%prs_res)
301 call pnpn_vel_res_factory(this%vel_res)
306 if (params%valid_path(
'case.fluid.nut_field'))
then
307 if (this%full_stress_formulation .eqv. .false.)
then
308 call neko_error(
"You need to set full_stress_formulation to " // &
309 "true for the fluid to have a spatially varying " // &
312 call json_get(params,
'case.fluid.nut_field', this%nut_field_name)
314 this%nut_field_name =
""
318 call ax_helm_factory(this%Ax_prs, full_formulation = .false.)
322 call rhs_maker_sumab_fctry(this%sumab)
325 call rhs_maker_ext_fctry(this%makeabf)
328 call rhs_maker_bdf_fctry(this%makebdf)
331 call rhs_maker_oifs_fctry(this%makeoifs)
334 associate(xh_lx => this%Xh%lx, xh_ly => this%Xh%ly, xh_lz => this%Xh%lz, &
335 dm_xh => this%dm_Xh, nelv => this%msh%nelv)
337 call this%p_res%init(dm_xh,
"p_res")
338 call this%u_res%init(dm_xh,
"u_res")
339 call this%v_res%init(dm_xh,
"v_res")
340 call this%w_res%init(dm_xh,
"w_res")
341 call this%abx1%init(dm_xh,
"abx1")
342 call this%aby1%init(dm_xh,
"aby1")
343 call this%abz1%init(dm_xh,
"abz1")
344 call this%abx2%init(dm_xh,
"abx2")
345 call this%aby2%init(dm_xh,
"aby2")
346 call this%abz2%init(dm_xh,
"abz2")
347 call this%advx%init(dm_xh,
"advx")
348 call this%advy%init(dm_xh,
"advy")
349 call this%advz%init(dm_xh,
"advz")
352 call this%du%init(this%dm_Xh,
'du')
353 call this%dv%init(this%dm_Xh,
'dv')
354 call this%dw%init(this%dm_Xh,
'dw')
355 call this%dp%init(this%dm_Xh,
'dp')
358 call this%setup_bcs(
user, params)
362 if (found)
call this%write_boundary_conditions()
364 call this%proj_prs%init(this%dm_Xh%size(), this%pr_projection_dim, &
365 this%pr_projection_activ_step)
367 call this%proj_vel%init(this%dm_Xh%size(), this%vel_projection_dim, &
368 this%vel_projection_activ_step)
374 if (params%valid_path(
'case.fluid.flow_rate_force'))
then
375 call this%vol_flow%init(this%dm_Xh, params)
379 call neko_log%section(
"Pressure solver")
382 'case.fluid.pressure_solver.max_iterations', &
384 call json_get(params,
'case.fluid.pressure_solver.type', solver_type)
385 call json_get(params,
'case.fluid.pressure_solver.preconditioner.type', &
388 'case.fluid.pressure_solver.preconditioner', precon_params)
389 call json_get(params,
'case.fluid.pressure_solver.absolute_tolerance', &
393 call neko_log%message(
'Type : ('// trim(solver_type) // &
394 ', ' // trim(precon_type) //
')')
395 write(log_buf,
'(A,ES13.6)')
'Abs tol :', abs_tol
398 call this%solver_factory(this%ksp_prs, this%dm_Xh%size(), &
399 solver_type, solver_maxiter, abs_tol, monitor)
400 call this%precon_factory_(this%pc_prs, this%ksp_prs, &
401 this%c_Xh, this%dm_Xh, this%gs_Xh, this%bcs_prs, &
402 precon_type, precon_params)
408 call advection_factory(this%adv, numerics_params, this%c_Xh, &
409 this%ulag, this%vlag, this%wlag, &
410 chkp%dtlag, chkp%tlag, this%ext_bdf, &
416 call this%chkp%init(this%u, this%v, this%w, this%p)
418 this%chkp%abx1 => this%abx1
419 this%chkp%abx2 => this%abx2
420 this%chkp%aby1 => this%aby1
421 this%chkp%aby2 => this%aby2
422 this%chkp%abz1 => this%abz1
423 this%chkp%abz2 => this%abz2
424 call this%chkp%add_lag(this%ulag, this%vlag, this%wlag)
432 subroutine fluid_pnpn_restart(this, chkp)
433 class(fluid_pnpn_t),
target,
intent(inout) :: this
434 type(chkp_t),
intent(inout) :: chkp
435 real(kind=rp) :: dtlag(10), tlag(10)
436 type(field_t) :: u_temp, v_temp, w_temp
442 n = this%u%dof%size()
443 if (
allocated(chkp%previous_mesh%elements) .or. &
444 chkp%previous_Xh%lx .ne. this%Xh%lx)
then
445 associate(u => this%u, v => this%v, w => this%w, p => this%p, &
446 c_xh => this%c_Xh, ulag => this%ulag, vlag => this%vlag, &
448 do concurrent(j = 1:n)
449 u%x(j,1,1,1) = u%x(j,1,1,1) * c_xh%mult(j,1,1,1)
450 v%x(j,1,1,1) = v%x(j,1,1,1) * c_xh%mult(j,1,1,1)
451 w%x(j,1,1,1) = w%x(j,1,1,1) * c_xh%mult(j,1,1,1)
452 p%x(j,1,1,1) = p%x(j,1,1,1) * c_xh%mult(j,1,1,1)
454 do i = 1, this%ulag%size()
455 do concurrent(j = 1:n)
456 ulag%lf(i)%x(j,1,1,1) = ulag%lf(i)%x(j,1,1,1) &
458 vlag%lf(i)%x(j,1,1,1) = vlag%lf(i)%x(j,1,1,1) &
460 wlag%lf(i)%x(j,1,1,1) = wlag%lf(i)%x(j,1,1,1) &
467 if (neko_bcknd_device .eq. 1)
then
468 associate(u => this%u, v => this%v, w => this%w, &
469 ulag => this%ulag, vlag => this%vlag, wlag => this%wlag,&
471 call device_memcpy(u%x, u%x_d, u%dof%size(), &
472 host_to_device, sync = .false.)
473 call device_memcpy(v%x, v%x_d, v%dof%size(), &
474 host_to_device, sync = .false.)
475 call device_memcpy(w%x, w%x_d, w%dof%size(), &
476 host_to_device, sync = .false.)
477 call device_memcpy(p%x, p%x_d, p%dof%size(), &
478 host_to_device, sync = .false.)
479 call device_memcpy(ulag%lf(1)%x, ulag%lf(1)%x_d, &
480 u%dof%size(), host_to_device, sync = .false.)
481 call device_memcpy(ulag%lf(2)%x, ulag%lf(2)%x_d, &
482 u%dof%size(), host_to_device, sync = .false.)
484 call device_memcpy(vlag%lf(1)%x, vlag%lf(1)%x_d, &
485 v%dof%size(), host_to_device, sync = .false.)
486 call device_memcpy(vlag%lf(2)%x, vlag%lf(2)%x_d, &
487 v%dof%size(), host_to_device, sync = .false.)
489 call device_memcpy(wlag%lf(1)%x, wlag%lf(1)%x_d, &
490 w%dof%size(), host_to_device, sync = .false.)
491 call device_memcpy(wlag%lf(2)%x, wlag%lf(2)%x_d, &
492 w%dof%size(), host_to_device, sync = .false.)
493 call device_memcpy(this%abx1%x, this%abx1%x_d, &
494 w%dof%size(), host_to_device, sync = .false.)
495 call device_memcpy(this%abx2%x, this%abx2%x_d, &
496 w%dof%size(), host_to_device, sync = .false.)
497 call device_memcpy(this%aby1%x, this%aby1%x_d, &
498 w%dof%size(), host_to_device, sync = .false.)
499 call device_memcpy(this%aby2%x, this%aby2%x_d, &
500 w%dof%size(), host_to_device, sync = .false.)
501 call device_memcpy(this%abz1%x, this%abz1%x_d, &
502 w%dof%size(), host_to_device, sync = .false.)
503 call device_memcpy(this%abz2%x, this%abz2%x_d, &
504 w%dof%size(), host_to_device, sync = .false.)
505 call device_memcpy(this%advx%x, this%advx%x_d, &
506 w%dof%size(), host_to_device, sync = .false.)
507 call device_memcpy(this%advy%x, this%advy%x_d, &
508 w%dof%size(), host_to_device, sync = .false.)
509 call device_memcpy(this%advz%x, this%advz%x_d, &
510 w%dof%size(), host_to_device, sync = .false.)
517 if (
allocated(chkp%previous_mesh%elements) &
518 .or. chkp%previous_Xh%lx .ne. this%Xh%lx)
then
519 call this%gs_Xh%op(this%u, gs_op_add)
520 call this%gs_Xh%op(this%v, gs_op_add)
521 call this%gs_Xh%op(this%w, gs_op_add)
522 call this%gs_Xh%op(this%p, gs_op_add)
524 do i = 1, this%ulag%size()
525 call this%gs_Xh%op(this%ulag%lf(i), gs_op_add)
526 call this%gs_Xh%op(this%vlag%lf(i), gs_op_add)
527 call this%gs_Xh%op(this%wlag%lf(i), gs_op_add)
614 subroutine fluid_pnpn_step(this, time, dt_controller)
615 class(fluid_pnpn_t),
target,
intent(inout) :: this
616 type(time_state_t),
intent(in) :: time
617 type(time_step_controller_t),
intent(in) :: dt_controller
621 type(ksp_monitor_t) :: ksp_results(4)
623 type(file_t) :: dump_file
624 class(bc_t),
pointer :: bc_i
625 type(non_normal_t),
pointer :: bc_j
627 if (this%freeze)
return
629 n = this%dm_Xh%size()
631 call profiler_start_region(
'Fluid', 1)
632 associate(u => this%u, v => this%v, w => this%w, p => this%p, &
633 u_e => this%u_e, v_e => this%v_e, w_e => this%w_e, &
634 du => this%du, dv => this%dv, dw => this%dw, dp => this%dp, &
635 u_res => this%u_res, v_res => this%v_res, w_res => this%w_res, &
636 p_res => this%p_res, ax_vel => this%Ax_vel, ax_prs => this%Ax_prs, &
638 c_xh => this%c_Xh, dm_xh => this%dm_Xh, gs_xh => this%gs_Xh, &
639 ulag => this%ulag, vlag => this%vlag, wlag => this%wlag, &
640 msh => this%msh, prs_res => this%prs_res, &
641 source_term => this%source_term, vel_res => this%vel_res, &
642 sumab => this%sumab, makeoifs => this%makeoifs, &
643 makeabf => this%makeabf, makebdf => this%makebdf, &
644 vel_projection_dim => this%vel_projection_dim, &
645 pr_projection_dim => this%pr_projection_dim, &
647 rho => this%rho, mu_tot => this%mu_tot, &
648 f_x => this%f_x, f_y => this%f_y, f_z => this%f_z, &
649 t => time%t, tstep => time%tstep, dt => time%dt, &
650 ext_bdf => this%ext_bdf, event => glb_cmd_event)
653 call sumab%compute_fluid(u_e, v_e, w_e, u, v, w, &
654 ulag, vlag, wlag, ext_bdf%advection_coeffs, ext_bdf%nadv)
657 call this%source_term%compute(time)
660 call this%bcs_vel%apply_vector(f_x%x, f_y%x, f_z%x, &
661 this%dm_Xh%size(), time, strong = .false.)
665 call this%adv%compute(u, v, w, &
666 this%advx, this%advy, this%advz, &
667 xh, this%c_Xh, dm_xh%size(), dt)
674 call makeabf%compute_fluid(this%abx1, this%aby1, this%abz1,&
675 this%abx2, this%aby2, this%abz2, &
676 f_x%x, f_y%x, f_z%x, &
677 rho%x(1,1,1,1), ext_bdf%advection_coeffs, n)
680 call makeoifs%compute_fluid(this%advx%x, this%advy%x, this%advz%x, &
681 f_x%x, f_y%x, f_z%x, &
682 rho%x(1,1,1,1), dt, n)
685 call this%adv%compute(u, v, w, &
687 xh, this%c_Xh, dm_xh%size())
693 call makeabf%compute_fluid(this%abx1, this%aby1, this%abz1,&
694 this%abx2, this%aby2, this%abz2, &
695 f_x%x, f_y%x, f_z%x, &
696 rho%x(1,1,1,1), ext_bdf%advection_coeffs, n)
699 call makebdf%compute_fluid(ulag, vlag, wlag, f_x%x, f_y%x, f_z%x, &
700 u, v, w, c_xh%B, rho%x(1,1,1,1), dt, &
701 ext_bdf%diffusion_coeffs, ext_bdf%ndiff, n)
708 call this%bc_apply_vel(time, strong = .true.)
709 call this%bc_apply_prs(time)
712 call this%update_material_properties(time)
715 call profiler_start_region(
'Pressure_residual', 18)
717 call prs_res%compute(p, p_res,&
722 this%bc_prs_surface, this%bc_sym_surface,&
723 ax_prs, ext_bdf%diffusion_coeffs(1), dt, &
727 if (.not. this%prs_dirichlet)
call ortho(p_res%x, this%glb_n_points, n)
729 call gs_xh%op(p_res, gs_op_add, event)
730 call device_event_sync(event)
733 call this%bclst_dp%apply_scalar(p_res%x, p%dof%size(), time)
736 call profiler_end_region(
'Pressure_residual', 18)
739 call this%proj_prs%pre_solving(p_res%x, tstep, c_xh, n, dt_controller, &
742 call this%pc_prs%update()
744 call profiler_start_region(
'Pressure_solve', 3)
748 this%ksp_prs%solve(ax_prs, dp, p_res%x, n, c_xh, &
749 this%bclst_dp, gs_xh)
750 ksp_results(1)%name =
'Pressure'
753 call profiler_end_region(
'Pressure_solve', 3)
755 call this%proj_prs%post_solving(dp%x, ax_prs, c_xh, &
756 this%bclst_dp, gs_xh, n, tstep, dt_controller)
759 call field_add2(p, dp, n)
760 if (.not. this%prs_dirichlet)
call ortho(p%x, this%glb_n_points, n)
763 call profiler_start_region(
'Velocity_residual', 19)
764 call vel_res%compute(ax_vel, u, v, w, &
765 u_res, v_res, w_res, &
769 mu_tot, rho, ext_bdf%diffusion_coeffs(1), &
772 call gs_xh%op(u_res, gs_op_add, event)
773 call device_event_sync(event)
774 call gs_xh%op(v_res, gs_op_add, event)
775 call device_event_sync(event)
776 call gs_xh%op(w_res, gs_op_add, event)
777 call device_event_sync(event)
780 call this%bclst_vel_res%apply(u_res, v_res, w_res, time)
783 call profiler_end_region(
'Velocity_residual', 19)
785 call this%proj_vel%pre_solving(u_res%x, v_res%x, w_res%x, &
786 tstep, c_xh, n, dt_controller,
'Velocity')
788 call this%pc_vel%update()
790 call profiler_start_region(
"Velocity_solve", 4)
791 ksp_results(2:4) = this%ksp_vel%solve_coupled(ax_vel, du, dv, dw, &
792 u_res%x, v_res%x, w_res%x, n, c_xh, &
793 this%bclst_du, this%bclst_dv, this%bclst_dw, gs_xh, &
794 this%ksp_vel%max_iter)
795 call profiler_end_region(
"Velocity_solve", 4)
796 if (this%full_stress_formulation)
then
797 ksp_results(2)%name =
'Momentum'
799 ksp_results(2)%name =
'X-Velocity'
800 ksp_results(3)%name =
'Y-Velocity'
801 ksp_results(4)%name =
'Z-Velocity'
804 call this%proj_vel%post_solving(du%x, dv%x, dw%x, ax_vel, c_xh, &
805 this%bclst_du, this%bclst_dv, this%bclst_dw, gs_xh, n, tstep, &
808 if (neko_bcknd_device .eq. 1)
then
809 call device_opadd2cm(u%x_d, v%x_d, w%x_d, &
810 du%x_d, dv%x_d, dw%x_d, 1.0_rp, n, msh%gdim)
812 call opadd2cm(u%x, v%x, w%x, du%x, dv%x, dw%x, 1.0_rp, n, msh%gdim)
815 if (this%forced_flow_rate)
then
817 call this%vol_flow%adjust( u, v, w, p, u_res, v_res, w_res, p_res, &
818 c_xh, gs_xh, ext_bdf, rho%x(1,1,1,1), mu_tot, &
819 dt, this%bclst_dp, this%bclst_du, this%bclst_dv, &
820 this%bclst_dw, this%bclst_vel_res, ax_vel, ax_prs, this%ksp_prs, &
821 this%ksp_vel, this%pc_prs, this%pc_vel, this%ksp_prs%max_iter, &
822 this%ksp_vel%max_iter)
825 call fluid_step_info(time, ksp_results, &
826 this%full_stress_formulation, this%strict_convergence)
829 call profiler_end_region(
'Fluid', 1)
834 subroutine fluid_pnpn_setup_bcs(this, user, params)
835 class(fluid_pnpn_t),
target,
intent(inout) :: this
836 type(user_t),
target,
intent(in) :: user
837 type(json_file),
intent(inout) :: params
838 integer :: i, n_bcs, zone_index, j, zone_size, global_zone_size, ierr
839 class(bc_t),
pointer :: bc_i
840 type(json_core) :: core
841 type(json_value),
pointer :: bc_object
842 type(json_file) :: bc_subdict
845 logical,
allocatable :: marked_zones(:)
846 integer,
allocatable :: zone_indices(:)
849 call this%bclst_vel_res%init()
850 call this%bclst_du%init()
851 call this%bclst_dv%init()
852 call this%bclst_dw%init()
853 call this%bclst_dp%init()
855 call this%bc_vel_res%init_from_components(this%c_Xh)
856 call this%bc_du%init_from_components(this%c_Xh)
857 call this%bc_dv%init_from_components(this%c_Xh)
858 call this%bc_dw%init_from_components(this%c_Xh)
859 call this%bc_dp%init_from_components(this%c_Xh)
862 call this%bc_prs_surface%init_from_components(this%c_Xh)
863 call this%bc_sym_surface%init_from_components(this%c_Xh)
866 if (params%valid_path(
'case.fluid.boundary_conditions'))
then
867 call params%info(
'case.fluid.boundary_conditions', n_children = n_bcs)
868 call params%get_core(core)
869 call params%get(
'case.fluid.boundary_conditions', bc_object, found)
874 call this%bcs_vel%init(n_bcs)
876 allocate(marked_zones(
size(this%msh%labeled_zones)))
877 marked_zones = .false.
881 call json_extract_item(core, bc_object, i, bc_subdict)
883 call json_get(bc_subdict,
"zone_indices", zone_indices)
888 do j = 1,
size(zone_indices)
889 zone_size = this%msh%labeled_zones(zone_indices(j))%size
890 call mpi_allreduce(zone_size, global_zone_size, 1, &
891 mpi_integer, mpi_max, neko_comm, ierr)
893 if (global_zone_size .eq. 0)
then
894 write(error_unit,
'(A, A, I0, A, A, I0, A)')
"*** ERROR ***: ",&
895 "Zone index ", zone_indices(j), &
896 " is invalid as this zone has 0 size, meaning it ", &
897 "is not in the mesh. Check fluid boundary condition ", &
902 if (marked_zones(zone_indices(j)) .eqv. .true.)
then
903 write(error_unit,
'(A, A, I0, A, A, A, A)')
"*** ERROR ***: ", &
904 "Zone with index ", zone_indices(j), &
905 " has already been assigned a boundary condition. ", &
906 "Please check your boundary_conditions entry for the ", &
907 "fluid and make sure that each zone index appears only ", &
908 "in a single boundary condition."
911 marked_zones(zone_indices(j)) = .true.
916 call velocity_bc_factory(bc_i, this, bc_subdict, this%c_Xh,
user)
920 if (
associated(bc_i))
then
934 call this%bclst_vel_res%append(bc_i)
935 call this%bc_du%mark_facets(bc_i%bc_x%marked_facet)
936 call this%bc_dv%mark_facets(bc_i%bc_y%marked_facet)
937 call this%bc_dw%mark_facets(bc_i%bc_z%marked_facet)
939 call this%bcs_vel%append(bc_i)
941 call this%bc_sym_surface%mark_facets(bc_i%marked_facet)
942 type is (non_normal_t)
945 call this%bclst_vel_res%append(bc_i)
946 call this%bc_du%mark_facets(bc_i%bc_x%marked_facet)
947 call this%bc_dv%mark_facets(bc_i%bc_y%marked_facet)
948 call this%bc_dw%mark_facets(bc_i%bc_z%marked_facet)
949 type is (shear_stress_t)
951 call this%bclst_vel_res%append(bc_i%symmetry)
952 call this%bclst_du%append(bc_i%symmetry%bc_x)
953 call this%bclst_dv%append(bc_i%symmetry%bc_y)
954 call this%bclst_dw%append(bc_i%symmetry%bc_z)
956 call this%bcs_vel%append(bc_i)
957 type is (wall_model_bc_t)
959 call this%bclst_vel_res%append(bc_i%symmetry)
960 call this%bclst_du%append(bc_i%symmetry%bc_x)
961 call this%bclst_dv%append(bc_i%symmetry%bc_y)
962 call this%bclst_dw%append(bc_i%symmetry%bc_z)
964 call this%bcs_vel%append(bc_i)
971 if (bc_i%strong .eqv. .true.)
then
972 call this%bc_vel_res%mark_facets(bc_i%marked_facet)
973 call this%bc_du%mark_facets(bc_i%marked_facet)
974 call this%bc_dv%mark_facets(bc_i%marked_facet)
975 call this%bc_dw%mark_facets(bc_i%marked_facet)
977 call this%bc_prs_surface%mark_facets(bc_i%marked_facet)
980 call this%bcs_vel%append(bc_i)
986 do i = 1,
size(this%msh%labeled_zones)
987 if ((this%msh%labeled_zones(i)%size .gt. 0) .and. &
988 (marked_zones(i) .eqv. .false.))
then
989 write(error_unit,
'(A, A, I0)')
"*** ERROR ***: ", &
990 "No fluid boundary condition assigned to zone ", i
998 call this%bcs_prs%init(n_bcs)
1002 call json_extract_item(core, bc_object, i, bc_subdict)
1004 call pressure_bc_factory(bc_i, this, bc_subdict, this%c_Xh,
user)
1008 if (
associated(bc_i))
then
1009 call this%bcs_prs%append(bc_i)
1012 if (bc_i%strong .eqv. .true.)
then
1013 call this%bc_dp%mark_facets(bc_i%marked_facet)
1021 do i = 1,
size(this%msh%labeled_zones)
1022 if (this%msh%labeled_zones(i)%size .gt. 0)
then
1023 call neko_error(
"No boundary_conditions entry in the case file!")
1029 call this%bcs_vel%init()
1030 call this%bcs_prs%init()
1034 call this%bc_prs_surface%finalize()
1035 call this%bc_sym_surface%finalize()
1037 call this%bc_vel_res%finalize()
1038 call this%bc_du%finalize()
1039 call this%bc_dv%finalize()
1040 call this%bc_dw%finalize()
1041 call this%bc_dp%finalize()
1043 call this%bclst_vel_res%append(this%bc_vel_res)
1044 call this%bclst_du%append(this%bc_du)
1045 call this%bclst_dv%append(this%bc_dv)
1046 call this%bclst_dw%append(this%bc_dw)
1047 call this%bclst_dp%append(this%bc_dp)
1050 this%prs_dirichlet = .not. this%bclst_dp%is_empty()
1051 call mpi_allreduce(mpi_in_place, this%prs_dirichlet, 1, &
1052 mpi_logical, mpi_lor, neko_comm)
1057 subroutine fluid_pnpn_write_boundary_conditions(this)
1063 class(fluid_pnpn_t),
target,
intent(inout) :: this
1064 type(dirichlet_t) :: bdry_mask
1065 type(field_t),
pointer :: bdry_field
1066 type(file_t) :: bdry_file
1067 integer :: temp_index, i
1068 class(bc_t),
pointer :: bci
1069 character(len=LOG_SIZE) :: log_buf
1071 call neko_log%section(
"Fluid boundary conditions")
1072 write(log_buf,
'(A)')
'Marking using integer keys in bdry0.f00000'
1073 call neko_log%message(log_buf)
1074 write(log_buf,
'(A)')
'Condition-value pairs: '
1075 call neko_log%message(log_buf)
1076 write(log_buf,
'(A)')
' no_slip = 1'
1077 call neko_log%message(log_buf)
1078 write(log_buf,
'(A)')
' velocity_value = 2'
1079 call neko_log%message(log_buf)
1080 write(log_buf,
'(A)')
' outflow, normal_outflow (+dong) = 3'
1081 call neko_log%message(log_buf)
1082 write(log_buf,
'(A)')
' symmetry = 4'
1083 call neko_log%message(log_buf)
1084 write(log_buf,
'(A)')
' periodic = 6'
1085 call neko_log%message(log_buf)
1086 write(log_buf,
'(A)')
' user_velocity = 7'
1087 call neko_log%message(log_buf)
1088 write(log_buf,
'(A)')
' user_pressure = 8'
1089 call neko_log%message(log_buf)
1090 write(log_buf,
'(A)')
' shear_stress = 9'
1091 call neko_log%message(log_buf)
1092 write(log_buf,
'(A)')
' wall_modelling = 10'
1093 call neko_log%message(log_buf)
1094 write(log_buf,
'(A)')
' blasius_profile = 11'
1095 call neko_log%message(log_buf)
1096 call neko_log%end_section()
1098 call this%scratch%request_field(bdry_field, temp_index)
1103 call bdry_mask%init_from_components(this%c_Xh, 6.0_rp)
1104 call bdry_mask%mark_zone(this%msh%periodic)
1105 call bdry_mask%finalize()
1106 call bdry_mask%apply_scalar(bdry_field%x, this%dm_Xh%size())
1107 call bdry_mask%free()
1109 do i = 1, this%bcs_prs%size()
1110 bci => this%bcs_prs%get(i)
1111 select type (
bc => bci)
1112 type is (zero_dirichlet_t)
1113 call bdry_mask%init_from_components(this%c_Xh, 3.0_rp)
1114 call bdry_mask%mark_facets(bci%marked_facet)
1115 call bdry_mask%finalize()
1116 call bdry_mask%apply_scalar(bdry_field%x, this%dm_Xh%size())
1117 call bdry_mask%free()
1119 call bdry_mask%init_from_components(this%c_Xh, 3.0_rp)
1120 call bdry_mask%mark_facets(bci%marked_facet)
1121 call bdry_mask%finalize()
1122 call bdry_mask%apply_scalar(bdry_field%x, this%dm_Xh%size())
1123 call bdry_mask%free()
1125 call bdry_mask%init_from_components(this%c_Xh, 8.0_rp)
1126 call bdry_mask%mark_facets(bci%marked_facet)
1127 call bdry_mask%finalize()
1128 call bdry_mask%apply_scalar(bdry_field%x, this%dm_Xh%size())
1129 call bdry_mask%free()
1133 do i = 1, this%bcs_vel%size()
1134 bci => this%bcs_vel%get(i)
1135 select type (
bc => bci)
1136 type is (zero_dirichlet_t)
1137 call bdry_mask%init_from_components(this%c_Xh, 1.0_rp)
1138 call bdry_mask%mark_facets(bci%marked_facet)
1139 call bdry_mask%finalize()
1140 call bdry_mask%apply_scalar(bdry_field%x, this%dm_Xh%size())
1141 call bdry_mask%free()
1143 call bdry_mask%init_from_components(this%c_Xh, 2.0_rp)
1144 call bdry_mask%mark_facets(bci%marked_facet)
1145 call bdry_mask%finalize()
1146 call bdry_mask%apply_scalar(bdry_field%x, this%dm_Xh%size())
1147 call bdry_mask%free()
1148 type is (symmetry_t)
1149 call bdry_mask%init_from_components(this%c_Xh, 4.0_rp)
1150 call bdry_mask%mark_facets(bci%marked_facet)
1151 call bdry_mask%finalize()
1152 call bdry_mask%apply_scalar(bdry_field%x, this%dm_Xh%size())
1153 call bdry_mask%free()
1155 call bdry_mask%init_from_components(this%c_Xh, 7.0_rp)
1156 call bdry_mask%mark_facets(bci%marked_facet)
1157 call bdry_mask%finalize()
1158 call bdry_mask%apply_scalar(bdry_field%x, this%dm_Xh%size())
1159 call bdry_mask%free()
1160 type is (shear_stress_t)
1161 call bdry_mask%init_from_components(this%c_Xh, 9.0_rp)
1162 call bdry_mask%mark_facets(bci%marked_facet)
1163 call bdry_mask%finalize()
1164 call bdry_mask%apply_scalar(bdry_field%x, this%dm_Xh%size())
1165 call bdry_mask%free()
1166 type is (wall_model_bc_t)
1167 call bdry_mask%init_from_components(this%c_Xh, 10.0_rp)
1168 call bdry_mask%mark_facets(bci%marked_facet)
1169 call bdry_mask%finalize()
1170 call bdry_mask%apply_scalar(bdry_field%x, this%dm_Xh%size())
1171 call bdry_mask%free()
1173 call bdry_mask%init_from_components(this%c_Xh, 11.0_rp)
1174 call bdry_mask%mark_facets(bci%marked_facet)
1175 call bdry_mask%finalize()
1176 call bdry_mask%apply_scalar(bdry_field%x, this%dm_Xh%size())
1177 call bdry_mask%free()
1182 call bdry_file%init(
'bdry.fld')
1183 call bdry_file%write(bdry_field)
1185 call this%scratch%relinquish_field(temp_index)