38 pnpn_prs_res_factory, pnpn_vel_res_factory, &
39 pnpn_prs_res_stress_factory, pnpn_vel_res_stress_factory
42 rhs_maker_ext_fctry, rhs_maker_oifs_fctry
53 use json_module,
only : json_file
76 type(
field_t) :: p_res, u_res, v_res, w_res
81 class(
ax_t),
allocatable :: ax_vel
83 class(
ax_t),
allocatable :: ax_prs
146 material_properties, time_scheme)
148 type(
mesh_t),
target,
intent(inout) :: msh
149 integer,
intent(inout) :: lx
150 type(json_file),
target,
intent(inout) :: params
151 type(
user_t),
intent(in) :: user
154 character(len=15),
parameter :: scheme =
'Modular (Pn/Pn)'
159 call this%scheme_init(msh, lx, params, .true., .true., scheme, user, &
162 if (this%variable_material_properties .eqv. .true.)
then
164 call ax_helm_factory(this%Ax_vel, full_formulation = .true.)
167 call pnpn_prs_res_stress_factory(this%prs_res)
170 call pnpn_vel_res_stress_factory(this%vel_res)
173 call ax_helm_factory(this%Ax_vel, full_formulation = .false.)
176 call pnpn_prs_res_factory(this%prs_res)
179 call pnpn_vel_res_factory(this%vel_res)
183 call ax_helm_factory(this%Ax_prs, full_formulation = .false.)
187 call rhs_maker_sumab_fctry(this%sumab)
190 call rhs_maker_ext_fctry(this%makeabf)
193 call rhs_maker_bdf_fctry(this%makebdf)
196 call rhs_maker_oifs_fctry(this%makeoifs)
199 associate(xh_lx => this%Xh%lx, xh_ly => this%Xh%ly, xh_lz => this%Xh%lz, &
200 dm_xh => this%dm_Xh, nelv => this%msh%nelv)
202 call this%p_res%init(dm_xh,
"p_res")
203 call this%u_res%init(dm_xh,
"u_res")
204 call this%v_res%init(dm_xh,
"v_res")
205 call this%w_res%init(dm_xh,
"w_res")
206 call this%abx1%init(dm_xh,
"abx1")
207 call this%aby1%init(dm_xh,
"aby1")
208 call this%abz1%init(dm_xh,
"abz1")
209 call this%abx2%init(dm_xh,
"abx2")
210 call this%aby2%init(dm_xh,
"aby2")
211 call this%abz2%init(dm_xh,
"abz2")
212 call this%advx%init(dm_xh,
"advx")
213 call this%advy%init(dm_xh,
"advy")
214 call this%advz%init(dm_xh,
"advz")
225 call this%du%init(dm_xh,
'du')
226 call this%dv%init(dm_xh,
'dv')
227 call this%dw%init(dm_xh,
'dw')
228 call this%dp%init(dm_xh,
'dp')
233 call this%bc_prs_surface%init_base(this%c_Xh)
234 call this%bc_prs_surface%mark_zone(msh%inlet)
235 call this%bc_prs_surface%mark_zones_from_list(msh%labeled_zones,&
239 call this%bc_prs_surface%mark_zones_from_list(msh%labeled_zones,&
240 'd_vel_u', this%bc_labels)
241 call this%bc_prs_surface%mark_zones_from_list(msh%labeled_zones,&
242 'd_vel_v', this%bc_labels)
243 call this%bc_prs_surface%mark_zones_from_list(msh%labeled_zones,&
244 'd_vel_w', this%bc_labels)
245 call this%bc_prs_surface%finalize()
247 call this%bc_sym_surface%init_base(this%c_Xh)
248 call this%bc_sym_surface%mark_zone(msh%sympln)
249 call this%bc_sym_surface%mark_zones_from_list(msh%labeled_zones,&
250 'sym', this%bc_labels)
252 call this%bc_sym_surface%finalize()
254 call this%bc_vel_res_non_normal%init_base(this%c_Xh)
255 call this%bc_vel_res_non_normal%mark_zone(msh%outlet_normal)
256 call this%bc_vel_res_non_normal%mark_zones_from_list(msh%labeled_zones,&
257 'on', this%bc_labels)
258 call this%bc_vel_res_non_normal%mark_zones_from_list(msh%labeled_zones,&
261 call this%bc_vel_res_non_normal%finalize()
262 call this%bc_vel_res_non_normal%init(this%c_Xh)
264 call this%bc_field_dirichlet_p%init_base(this%c_Xh)
265 call this%bc_field_dirichlet_p%mark_zones_from_list(msh%labeled_zones, &
266 'on+dong', this%bc_labels)
267 call this%bc_field_dirichlet_p%mark_zones_from_list(msh%labeled_zones, &
268 'o+dong', this%bc_labels)
269 call this%bc_field_dirichlet_p%mark_zones_from_list(msh%labeled_zones, &
270 'd_pres', this%bc_labels)
271 call this%bc_field_dirichlet_p%finalize()
272 call this%bc_field_dirichlet_p%set_g(0.0_rp)
274 call bc_list_add(this%bclst_dp, this%bc_field_dirichlet_p)
278 call this%bc_field_dirichlet_u%init_base(this%c_Xh)
279 call this%bc_field_dirichlet_u%mark_zones_from_list( &
280 msh%labeled_zones,
'd_vel_u', this%bc_labels)
281 call this%bc_field_dirichlet_u%finalize()
282 call this%bc_field_dirichlet_u%set_g(0.0_rp)
284 call this%bc_field_dirichlet_v%init_base(this%c_Xh)
285 call this%bc_field_dirichlet_v%mark_zones_from_list(msh%labeled_zones, &
288 call this%bc_field_dirichlet_v%finalize()
289 call this%bc_field_dirichlet_v%set_g(0.0_rp)
291 call this%bc_field_dirichlet_w%init_base(this%c_Xh)
292 call this%bc_field_dirichlet_w%mark_zones_from_list(msh%labeled_zones, &
295 call this%bc_field_dirichlet_w%finalize()
296 call this%bc_field_dirichlet_w%set_g(0.0_rp)
298 call this%bc_vel_res%init_base(this%c_Xh)
299 call this%bc_vel_res%mark_zone(msh%inlet)
300 call this%bc_vel_res%mark_zone(msh%wall)
301 call this%bc_vel_res%mark_zones_from_list(msh%labeled_zones, &
303 call this%bc_vel_res%mark_zones_from_list(msh%labeled_zones, &
305 call this%bc_vel_res%finalize()
306 call this%bc_vel_res%set_g(0.0_rp)
308 call bc_list_add(this%bclst_vel_res, this%bc_vel_res)
309 call bc_list_add(this%bclst_vel_res, this%bc_vel_res_non_normal)
315 call bc_list_add(this%bclst_du, this%bc_vel_res_non_normal%bc_x)
317 call bc_list_add(this%bclst_du, this%bc_field_dirichlet_u)
321 call bc_list_add(this%bclst_dv, this%bc_vel_res_non_normal%bc_y)
323 call bc_list_add(this%bclst_dv, this%bc_field_dirichlet_v)
327 call bc_list_add(this%bclst_dw, this%bc_vel_res_non_normal%bc_z)
329 call bc_list_add(this%bclst_dw, this%bc_field_dirichlet_w)
333 if (this%variable_material_properties .and. &
334 this%vel_projection_dim .gt. 0)
then
335 call neko_error(
"Velocity projection not available for full stress &
340 call this%proj_prs%init(this%dm_Xh%size(), this%pr_projection_dim, &
341 this%pr_projection_activ_step)
343 call this%proj_u%init(this%dm_Xh%size(), this%vel_projection_dim, &
344 this%vel_projection_activ_step)
345 call this%proj_v%init(this%dm_Xh%size(), this%vel_projection_dim, &
346 this%vel_projection_activ_step)
347 call this%proj_w%init(this%dm_Xh%size(), this%vel_projection_dim, &
348 this%vel_projection_activ_step)
352 call this%chkp%add_lag(this%ulag, this%vlag, this%wlag)
358 call advection_factory(this%adv, params, this%c_Xh, &
359 this%ulag, this%vlag, this%wlag, &
362 if (params%valid_path(
'case.fluid.flow_rate_force'))
then
363 call this%vol_flow%init(this%dm_Xh, params)
370 real(kind=rp) :: dtlag(10), tlag(10)
371 type(field_t) :: u_temp, v_temp, w_temp
374 n = this%u%dof%size()
378 call col2(this%u%x, this%c_Xh%mult, this%u%dof%size())
379 call col2(this%v%x, this%c_Xh%mult, this%u%dof%size())
380 call col2(this%w%x, this%c_Xh%mult, this%u%dof%size())
381 call col2(this%p%x, this%c_Xh%mult, this%u%dof%size())
382 do i = 1, this%ulag%size()
383 call col2(this%ulag%lf(i)%x, this%c_Xh%mult, this%u%dof%size())
384 call col2(this%vlag%lf(i)%x, this%c_Xh%mult, this%u%dof%size())
385 call col2(this%wlag%lf(i)%x, this%c_Xh%mult, this%u%dof%size())
388 if (neko_bcknd_device .eq. 1)
then
389 associate(u => this%u, v => this%v, w => this%w, &
390 ulag => this%ulag, vlag => this%vlag, wlag => this%wlag,&
392 call device_memcpy(u%x, u%x_d, u%dof%size(), &
393 host_to_device, sync = .false.)
394 call device_memcpy(v%x, v%x_d, v%dof%size(), &
395 host_to_device, sync = .false.)
396 call device_memcpy(w%x, w%x_d, w%dof%size(), &
397 host_to_device, sync = .false.)
398 call device_memcpy(p%x, p%x_d, p%dof%size(), &
399 host_to_device, sync = .false.)
400 call device_memcpy(ulag%lf(1)%x, ulag%lf(1)%x_d, &
401 u%dof%size(), host_to_device, sync = .false.)
402 call device_memcpy(ulag%lf(2)%x, ulag%lf(2)%x_d, &
403 u%dof%size(), host_to_device, sync = .false.)
405 call device_memcpy(vlag%lf(1)%x, vlag%lf(1)%x_d, &
406 v%dof%size(), host_to_device, sync = .false.)
407 call device_memcpy(vlag%lf(2)%x, vlag%lf(2)%x_d, &
408 v%dof%size(), host_to_device, sync = .false.)
410 call device_memcpy(wlag%lf(1)%x, wlag%lf(1)%x_d, &
411 w%dof%size(), host_to_device, sync = .false.)
412 call device_memcpy(wlag%lf(2)%x, wlag%lf(2)%x_d, &
413 w%dof%size(), host_to_device, sync = .false.)
414 call device_memcpy(this%abx1%x, this%abx1%x_d, &
415 w%dof%size(), host_to_device, sync = .false.)
416 call device_memcpy(this%abx2%x, this%abx2%x_d, &
417 w%dof%size(), host_to_device, sync = .false.)
418 call device_memcpy(this%aby1%x, this%aby1%x_d, &
419 w%dof%size(), host_to_device, sync = .false.)
420 call device_memcpy(this%aby2%x, this%aby2%x_d, &
421 w%dof%size(), host_to_device, sync = .false.)
422 call device_memcpy(this%abz1%x, this%abz1%x_d, &
423 w%dof%size(), host_to_device, sync = .false.)
424 call device_memcpy(this%abz2%x, this%abz2%x_d, &
425 w%dof%size(), host_to_device, sync = .false.)
426 call device_memcpy(this%advx%x, this%advx%x_d, &
427 w%dof%size(), host_to_device, sync = .false.)
428 call device_memcpy(this%advy%x, this%advy%x_d, &
429 w%dof%size(), host_to_device, sync = .false.)
430 call device_memcpy(this%advz%x, this%advz%x_d, &
431 w%dof%size(), host_to_device, sync = .false.)
436 call this%gs_Xh%op(this%u, gs_op_add)
437 call this%gs_Xh%op(this%v, gs_op_add)
438 call this%gs_Xh%op(this%w, gs_op_add)
439 call this%gs_Xh%op(this%p, gs_op_add)
441 do i = 1, this%ulag%size()
442 call this%gs_Xh%op(this%ulag%lf(i), gs_op_add)
443 call this%gs_Xh%op(this%vlag%lf(i), gs_op_add)
444 call this%gs_Xh%op(this%wlag%lf(i), gs_op_add)
509 call this%scheme_free()
511 call this%bc_prs_surface%free()
512 call this%bc_sym_surface%free()
513 call bc_list_free(this%bclst_vel_res)
514 call bc_list_free(this%bclst_dp)
515 call this%proj_prs%free()
516 call this%proj_u%free()
517 call this%proj_v%free()
518 call this%proj_w%free()
520 call this%p_res%free()
521 call this%u_res%free()
522 call this%v_res%free()
523 call this%w_res%free()
530 call this%abx1%free()
531 call this%aby1%free()
532 call this%abz1%free()
534 call this%abx2%free()
535 call this%aby2%free()
536 call this%abz2%free()
538 call this%advx%free()
539 call this%advy%free()
540 call this%advz%free()
542 if (
allocated(this%Ax_vel))
then
543 deallocate(this%Ax_vel)
546 if (
allocated(this%Ax_prs))
then
547 deallocate(this%Ax_prs)
550 if (
allocated(this%prs_res))
then
551 deallocate(this%prs_res)
554 if (
allocated(this%vel_res))
then
555 deallocate(this%vel_res)
558 if (
allocated(this%sumab))
then
559 deallocate(this%sumab)
562 if (
allocated(this%makeabf))
then
563 deallocate(this%makeabf)
566 if (
allocated(this%makebdf))
then
567 deallocate(this%makebdf)
570 if (
allocated(this%makeoifs))
then
571 deallocate(this%makeoifs)
574 call this%vol_flow%free()
586 real(kind=rp),
intent(inout) :: t
587 integer,
intent(inout) :: tstep
588 real(kind=rp),
intent(in) :: dt
589 type(time_scheme_controller_t),
intent(inout) :: ext_bdf
590 type(time_step_controller_t),
intent(in) :: dt_controller
594 type(ksp_monitor_t) :: ksp_results(4)
596 type(field_t),
pointer :: u_e, v_e, w_e
598 integer :: temp_indices(3)
600 if (this%freeze)
return
602 n = this%dm_Xh%size()
604 call profiler_start_region(
'Fluid', 1)
605 associate(u => this%u, v => this%v, w => this%w, p => this%p, &
606 du => this%du, dv => this%dv, dw => this%dw, dp => this%dp, &
607 u_res => this%u_res, v_res => this%v_res, w_res => this%w_res, &
608 p_res => this%p_res, ax_vel => this%Ax_vel, ax_prs => this%Ax_prs, &
610 c_xh => this%c_Xh, dm_xh => this%dm_Xh, gs_xh => this%gs_Xh, &
611 ulag => this%ulag, vlag => this%vlag, wlag => this%wlag, &
612 msh => this%msh, prs_res => this%prs_res, &
613 source_term => this%source_term, vel_res => this%vel_res, &
614 sumab => this%sumab, makeoifs => this%makeoifs, &
615 makeabf => this%makeabf, makebdf => this%makebdf, &
616 vel_projection_dim => this%vel_projection_dim, &
617 pr_projection_dim => this%pr_projection_dim, &
618 rho => this%rho, mu => this%mu, oifs => this%oifs, &
619 rho_field => this%rho_field, mu_field => this%mu_field, &
620 f_x => this%f_x, f_y => this%f_y, f_z => this%f_z, &
621 if_variable_dt => dt_controller%if_variable_dt, &
622 dt_last_change => dt_controller%dt_last_change)
625 call this%scratch%request_field(u_e, temp_indices(1))
626 call this%scratch%request_field(v_e, temp_indices(2))
627 call this%scratch%request_field(w_e, temp_indices(3))
628 call sumab%compute_fluid(u_e, v_e, w_e, u, v, w, &
629 ulag, vlag, wlag, ext_bdf%advection_coeffs, ext_bdf%nadv)
632 call this%source_term%compute(t, tstep)
635 if (neko_bcknd_device .eq. 1)
then
636 call device_opcolv(f_x%x_d, f_y%x_d, f_z%x_d, c_xh%B_d, msh%gdim, n)
638 call opcolv(f_x%x, f_y%x, f_z%x, c_xh%B, msh%gdim, n)
643 call this%adv%compute(u, v, w, &
644 this%advx, this%advy, this%advz, &
645 xh, this%c_Xh, dm_xh%size(), dt)
651 call makeabf%compute_fluid(this%abx1, this%aby1, this%abz1,&
652 this%abx2, this%aby2, this%abz2, &
653 f_x%x, f_y%x, f_z%x, &
654 rho, ext_bdf%advection_coeffs, n)
657 call makeoifs%compute_fluid(this%advx%x, this%advy%x, this%advz%x, &
658 f_x%x, f_y%x, f_z%x, &
662 call this%adv%compute(u, v, w, &
664 xh, this%c_Xh, dm_xh%size())
670 call makeabf%compute_fluid(this%abx1, this%aby1, this%abz1,&
671 this%abx2, this%aby2, this%abz2, &
672 f_x%x, f_y%x, f_z%x, &
673 rho, ext_bdf%advection_coeffs, n)
676 call makebdf%compute_fluid(ulag, vlag, wlag, f_x%x, f_y%x, f_z%x, &
677 u, v, w, c_xh%B, rho, dt, &
678 ext_bdf%diffusion_coeffs, ext_bdf%ndiff, n)
689 call this%user_field_bc_vel%update(this%user_field_bc_vel%field_list, &
690 this%user_field_bc_vel%bc_list, this%c_Xh, t, tstep,
"fluid")
692 call this%bc_apply_vel(t, tstep)
693 call this%bc_apply_prs(t, tstep)
696 call this%update_material_properties()
699 call profiler_start_region(
'Pressure_residual', 18)
700 call prs_res%compute(p, p_res,&
705 this%bc_prs_surface, this%bc_sym_surface,&
706 ax_prs, ext_bdf%diffusion_coeffs(1), dt, &
709 call gs_xh%op(p_res, gs_op_add)
710 call bc_list_apply_scalar(this%bclst_dp, p_res%x, p%dof%size(), t, tstep)
711 call profiler_end_region(
'Pressure_residual', 18)
713 call this%proj_prs%pre_solving(p_res%x, tstep, c_xh, n, dt_controller, &
716 call this%pc_prs%update()
717 call profiler_start_region(
'Pressure_solve', 3)
719 this%ksp_prs%solve(ax_prs, dp, p_res%x, n, c_xh, this%bclst_dp, gs_xh)
721 call profiler_end_region(
'Pressure_solve', 3)
723 call this%proj_prs%post_solving(dp%x, ax_prs, c_xh, &
724 this%bclst_dp, gs_xh, n, tstep, dt_controller)
726 call field_add2(p, dp, n)
729 call profiler_start_region(
'Velocity_residual', 19)
730 call vel_res%compute(ax_vel, u, v, w, &
731 u_res, v_res, w_res, &
735 mu_field, rho_field, ext_bdf%diffusion_coeffs(1), &
738 call gs_xh%op(u_res, gs_op_add)
739 call gs_xh%op(v_res, gs_op_add)
740 call gs_xh%op(w_res, gs_op_add)
742 call bc_list_apply_vector(this%bclst_vel_res,&
743 u_res%x, v_res%x, w_res%x, dm_xh%size(),&
748 if (neko_bcknd_device .eq. 1)
then
749 call this%bc_field_dirichlet_u%apply_scalar_dev(u_res%x_d, t, tstep)
750 call this%bc_field_dirichlet_v%apply_scalar_dev(v_res%x_d, t, tstep)
751 call this%bc_field_dirichlet_w%apply_scalar_dev(w_res%x_d, t, tstep)
753 call this%bc_field_dirichlet_u%apply_scalar(u_res%x, n, t, tstep)
754 call this%bc_field_dirichlet_v%apply_scalar(v_res%x, n, t, tstep)
755 call this%bc_field_dirichlet_w%apply_scalar(w_res%x, n, t, tstep)
758 call profiler_end_region(
'Velocity_residual', 19)
760 call this%proj_u%pre_solving(u_res%x, tstep, c_xh, n, dt_controller)
761 call this%proj_v%pre_solving(v_res%x, tstep, c_xh, n, dt_controller)
762 call this%proj_w%pre_solving(w_res%x, tstep, c_xh, n, dt_controller)
764 call this%pc_vel%update()
766 call profiler_start_region(
"Velocity_solve", 4)
767 ksp_results(2:4) = this%ksp_vel%solve_coupled(ax_vel, du, dv, dw, &
768 u_res%x, v_res%x, w_res%x, n, c_xh, &
769 this%bclst_du, this%bclst_dv, this%bclst_dw, gs_xh)
770 call profiler_end_region(
"Velocity_solve", 4)
772 call this%proj_u%post_solving(du%x, ax_vel, c_xh, &
773 this%bclst_du, gs_xh, n, tstep, dt_controller)
774 call this%proj_v%post_solving(dv%x, ax_vel, c_xh, &
775 this%bclst_dv, gs_xh, n, tstep, dt_controller)
776 call this%proj_w%post_solving(dw%x, ax_vel, c_xh, &
777 this%bclst_dw, gs_xh, n, tstep, dt_controller)
779 if (neko_bcknd_device .eq. 1)
then
780 call device_opadd2cm(u%x_d, v%x_d, w%x_d, &
781 du%x_d, dv%x_d, dw%x_d, 1.0_rp, n, msh%gdim)
783 call opadd2cm(u%x, v%x, w%x, du%x, dv%x, dw%x, 1.0_rp, n, msh%gdim)
786 if (this%forced_flow_rate)
then
787 call this%vol_flow%adjust( u, v, w, p, u_res, v_res, w_res, p_res, &
788 c_xh, gs_xh, ext_bdf, rho, mu,&
789 dt, this%bclst_dp, this%bclst_du, this%bclst_dv, &
790 this%bclst_dw, this%bclst_vel_res, ax_vel, this%ksp_prs, &
791 this%ksp_vel, this%pc_prs, this%pc_vel, this%ksp_prs%max_iter, &
792 this%ksp_vel%max_iter)
795 call fluid_step_info(tstep, t, dt, ksp_results)
797 call this%scratch%relinquish_field(temp_indices)
800 call profiler_end_region(
'Fluid', 1)
Copy data between host and device (or device and device)
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.
Subroutines to add advection terms to the RHS of a transport equation.
Defines a Matrix-vector product.
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.
subroutine, public device_opadd2cm(a1_d, a2_d, a3_d, b1_d, b2_d, b3_d, c, n, gdim)
subroutine, public device_opcolv(a1_d, a2_d, a3_d, c_d, gdim, n)
Device abstraction, common interface for various accelerators.
integer, parameter, public host_to_device
Defines a dirichlet boundary condition.
Dirichlet condition applied in the facet normal direction.
subroutine, public field_add2(a, b, n)
Vector addition .
Auxiliary routines for fluid solvers.
subroutine, public fluid_step_info(step, t, dt, ksp_results)
Prints for prs, velx, vely, velz the following: Number of iterations, start residual,...
Modular version of the Classic Nek5000 Pn/Pn formulation for fluids.
subroutine fluid_pnpn_step(this, t, tstep, dt, ext_bdf, dt_controller)
Advance fluid simulation in time.
subroutine fluid_pnpn_free(this)
subroutine fluid_pnpn_init(this, msh, lx, params, user, material_properties, time_scheme)
subroutine fluid_pnpn_restart(this, dtlag, tlag)
Defines Gather-scatter operations.
integer, parameter, public gs_op_add
Utilities for retrieving parameters from the case files.
Implements the base abstract type for Krylov solvers plus helper types.
Implements material_properties_t type.
subroutine, public col2(a, b, n)
Vector multiplication .
Collection of vector field operations operating on and . Note that in general the indices and ....
subroutine, public opcolv(a1, a2, a3, c, gdim, n)
subroutine, public opadd2cm(a1, a2, a3, b1, b2, b3, c, n, gdim)
integer, parameter neko_bcknd_device
Dirichlet condition on axis aligned plane in the non normal direction.
integer, parameter, public rp
Global precision used in computations.
Defines Pressure and velocity residuals in the Pn-Pn formulation.
subroutine, public profiler_start_region(name, region_id)
Started a named (name) profiler region.
subroutine, public profiler_end_region(name, region_id)
End the most recently started profiler region.
Project x onto X, the space of old solutions and back again.
Routines to generate the right-hand sides for the convection-diffusion equation. Employs the EXT/BDF ...
Implements the source_term_t type and a wrapper source_term_wrapper_t.
Compound scheme for the advection and diffusion operators in a transport equation.
Base class for time integration schemes.
Implements type time_step_controller.
Interfaces for user interaction with NEKO.
Base abstract type for computing the advection operator.
Base type for a matrix-vector product providing .
A list of boundary conditions.
Generic Dirichlet boundary condition on .
Dirichlet condition in facet normal direction.
Base type of all fluid formulations.
Type for storing initial and final residuals in a Krylov solver.
Contains all the material properties necessary in the simulation.
Dirichlet condition in non normal direction of a plane.
Abstract type to compute pressure residual.
Abstract type to compute velocity residual.
Abstract type to add contributions to F from lagged BD terms.
Abstract type to sum up contributions to kth order extrapolation scheme.
Abstract type to add contributions of kth order OIFS scheme.
Abstract type to compute extrapolated velocity field for the pressure equation.
Implements the logic to compute the time coefficients for the advection and diffusion operators in a ...
Provides a tool to set time step dt.