22#ifndef INCLUDE_EXADG_INCOMPRESSIBLE_NAVIER_STOKES_TIME_INTEGRATION_CREATE_TIME_INTEGRATOR_H_
23#define INCLUDE_EXADG_INCOMPRESSIBLE_NAVIER_STOKES_TIME_INTEGRATION_CREATE_TIME_INTEGRATOR_H_
25#include <exadg/incompressible_navier_stokes/time_integration/time_int_bdf.h>
26#include <exadg/incompressible_navier_stokes/time_integration/time_int_bdf_coupled_solver.h>
27#include <exadg/incompressible_navier_stokes/time_integration/time_int_bdf_dual_splitting.h>
28#include <exadg/incompressible_navier_stokes/time_integration/time_int_bdf_pressure_correction.h>
37template<
int dim,
typename Number>
38std::shared_ptr<TimeIntBDF<dim, Number>>
39create_time_integrator(std::shared_ptr<SpatialOperatorBase<dim, Number>> pde_operator,
40 std::shared_ptr<HelpersALE<dim, Number>
const> helpers_ale,
41 std::shared_ptr<PostProcessorInterface<Number>> postprocessor,
42 Parameters
const & parameters,
43 MPI_Comm
const & mpi_comm,
46 std::shared_ptr<TimeIntBDF<dim, Number>> time_integrator;
48 if(parameters.temporal_discretization == TemporalDiscretization::BDFCoupledSolution)
50 std::shared_ptr<OperatorCoupled<dim, Number>> operator_coupled =
51 std::dynamic_pointer_cast<OperatorCoupled<dim, Number>>(pde_operator);
53 time_integrator = std::make_shared<IncNS::TimeIntBDFCoupled<dim, Number>>(
54 operator_coupled, helpers_ale, postprocessor, parameters, mpi_comm, is_test);
56 else if(parameters.temporal_discretization == TemporalDiscretization::BDFDualSplittingScheme)
58 std::shared_ptr<OperatorDualSplitting<dim, Number>> operator_dual_splitting =
59 std::dynamic_pointer_cast<OperatorDualSplitting<dim, Number>>(pde_operator);
61 time_integrator = std::make_shared<IncNS::TimeIntBDFDualSplitting<dim, Number>>(
62 operator_dual_splitting, helpers_ale, postprocessor, parameters, mpi_comm, is_test);
64 else if(parameters.temporal_discretization == TemporalDiscretization::BDFPressureCorrection)
66 std::shared_ptr<OperatorPressureCorrection<dim, Number>> operator_pressure_correction =
67 std::dynamic_pointer_cast<OperatorPressureCorrection<dim, Number>>(pde_operator);
69 time_integrator = std::make_shared<IncNS::TimeIntBDFPressureCorrection<dim, Number>>(
70 operator_pressure_correction, helpers_ale, postprocessor, parameters, mpi_comm, is_test);
74 AssertThrow(
false, dealii::ExcMessage(
"Not implemented."));
77 return time_integrator;