17 #ifndef __deal2__time_dependent_h
18 #define __deal2__time_dependent_h
24 #include <deal.II/base/config.h>
26 #include <deal.II/base/subscriptor.h>
27 #include <deal.II/base/smartpointer.h>
36 template <
typename number>
class Vector;
445 const TimeSteppingData &data_dual,
446 const TimeSteppingData &data_postprocess);
612 template <
typename InitFunctionObject,
typename LoopFunctionObject>
613 void do_loop (InitFunctionObject init_function,
614 LoopFunctionObject loop_function,
615 const TimeSteppingData ×tepping_data,
704 std::vector<SmartPointer<TimeStepBase,TimeDependent> >
timesteps;
744 void end_sweep (
const unsigned int begin_timestep,
745 const unsigned int end_timestep);
768 primal_problem = 0x0,
811 virtual void wake_up (
const unsigned int);
830 virtual void sleep (
const unsigned int);
1227 <<
"The parameter " << arg1 <<
" has an invalid value.");
1502 <<
"The following value does not fulfill the requirements: " << arg1);
1548 <<
"The following value does not fulfill the requirements: " << arg1);
1593 grid_refinement = 0x1000
1681 virtual void wake_up (
const unsigned int wakeup_level);
1708 virtual void sleep (
const unsigned int);
1861 template <
typename InitFunctionObject,
typename LoopFunctionObject>
1863 LoopFunctionObject loop_function,
1878 const unsigned int n_timesteps =
timesteps.size();
1882 for (
unsigned int step=0; step<n_timesteps; ++step)
1886 init_function (static_cast<typename InitFunctionObject::argument_type>
1890 init_function (static_cast<typename InitFunctionObject::argument_type>
1891 (&*
timesteps[n_timesteps-@ref step_1
"step-1"]));
1897 for (
int step=-timestepping_data.
look_ahead; step<0; ++step)
1898 for (
int look_ahead=0;
1899 look_ahead<=static_cast<int>(timestepping_data.
look_ahead); ++look_ahead)
1903 if (step+look_ahead >= 0)
1904 timesteps[step+look_ahead]->wake_up(look_ahead);
1907 if (n_timesteps-(step+look_ahead) < n_timesteps)
1908 timesteps[n_timesteps-(step+look_ahead)]->wake_up(look_ahead);
1913 for (
unsigned int step=0; step<n_timesteps; ++step)
1917 for (
unsigned int look_ahead=0;
1918 look_ahead<=timestepping_data.
look_ahead; ++look_ahead)
1922 if (step+look_ahead < n_timesteps)
1923 timesteps[step+look_ahead]->wake_up(look_ahead);
1926 if (n_timesteps > (step+look_ahead))
1927 timesteps[n_timesteps-(step+look_ahead)-1]->wake_up(look_ahead);
1936 loop_function (static_cast<typename LoopFunctionObject::argument_type>
1940 loop_function (static_cast<typename LoopFunctionObject::argument_type>
1941 (&*
timesteps[n_timesteps-@ref step_1
"step-1"]));
1946 for (
unsigned int look_back=0;
1947 look_back<=timestepping_data.
look_back; ++look_back)
1951 if (step>=look_back)
1952 timesteps[step-look_back]->sleep(look_back);
1955 if (n_timesteps-(step-look_back) <= n_timesteps)
1956 timesteps[n_timesteps-(step-look_back)-1]->sleep(look_back);
1962 for (
int step=n_timesteps;
1963 step<static_cast<int>(n_timesteps+timestepping_data.
look_back); ++step)
1964 for (
int look_back=0;
1965 look_back<=static_cast<int>(timestepping_data.
look_back); ++look_back)
1969 if ((step-look_back >= 0)
1971 (step-look_back < static_cast<int>(n_timesteps)))
1972 timesteps[step-look_back]->sleep(look_back);
1975 if ((step-look_back >= 0)
1977 (step-look_back < static_cast<int>(n_timesteps)))
1978 timesteps[n_timesteps-(step-look_back)-1]->sleep(look_back);
1983 DEAL_II_NAMESPACE_CLOSE
void set_sweep_no(const unsigned int sweep_no)
virtual void solve_primal_problem()=0
const unsigned int look_back
std::size_t memory_consumption() const
virtual void start_sweep()
void refine_grid(const RefinementData data)
virtual void solve_dual_problem()
TimeDependent(const TimeSteppingData &data_primal, const TimeSteppingData &data_dual, const TimeSteppingData &data_postprocess)
void set_timestep_no(const unsigned int step_no)
void solve_primal_problem()
void add_timestep(TimeStepBase *new_timestep)
void solve_dual_problem()
DeclException1(ExcInvalidValue, double,<< "The following value does not fulfill the requirements: "<< arg1)
const TimeStepBase * next_timestep
const TimeSteppingData timestepping_data_postprocess
const unsigned int wakeup_level_to_build_grid
DeclException0(ExcGridNotDeleted)
virtual void init_for_primal_problem()
void set_next_timestep(const TimeStepBase *next)
virtual ~TimeStepBase_Tria()
const TimeStepBase * previous_timestep
const TimeSteppingData timestepping_data_dual
virtual void postprocess_timestep()
static CorrectionRelaxations default_correction_relaxations
double get_forward_timestep() const
const unsigned int cell_number_correction_steps
DeclException0(ExcInvalidPosition)
const unsigned int max_refinement_level
const double cell_number_corridor_top
void set_previous_timestep(const TimeStepBase *previous)
void do_loop(InitFunctionObject init_function, LoopFunctionObject loop_function, const TimeSteppingData ×tepping_data, const Direction direction)
virtual void sleep(const unsigned int)
virtual void init_for_refinement()
const RefinementFlags refinement_flags
const std::vector< std::vector< std::pair< unsigned int, double > > > correction_relaxations
virtual void init_for_dual_problem()
virtual std::size_t memory_consumption() const
SmartPointer< const Triangulation< dim, dim >, TimeStepBase_Tria< dim > > coarse_grid
RefinementData(const double refinement_threshold, const double coarsening_threshold=0)
double get_backward_timestep() const
TimeStepBase & operator=(const TimeStepBase &)
virtual std::size_t memory_consumption() const
virtual void get_tria_refinement_criteria(Vector< float > &criteria) const =0
const TimeSteppingData timestepping_data_primal
DeclException1(ExcInvalidParameter, int,<< "The parameter "<< arg1<< " has an invalid value.")
BlockCompressedSparsityPattern CompressedBlockSparsityPattern DEAL_II_DEPRECATED
virtual void init_for_postprocessing()
const unsigned int first_sweep_with_correction
const unsigned int sleep_level_to_delete_grid
DeclException1(ExcInvalidValue, double,<< "The following value does not fulfill the requirements: "<< arg1)
const double coarsening_threshold
virtual void sleep(const unsigned int)
void delete_timestep(const unsigned int position)
TimeStepBase_Tria_Flags::Flags< dim > Flags
unsigned int get_timestep_no() const
RefinementFlags(const unsigned int max_refinement_level=0, const unsigned int first_sweep_with_correction=0, const unsigned int min_cells_for_correction=0, const double cell_number_corridor_top=(1<< dim), const double cell_number_corridor_bottom=1, const CorrectionRelaxations &correction_relaxations=CorrectionRelaxations(), const unsigned int cell_number_correction_steps=0, const bool mirror_flags_to_previous_grid=false, const bool adapt_grids=false)
std::vector< SmartPointer< TimeStepBase, TimeDependent > > timesteps
const double cell_number_corridor_bottom
const unsigned int look_ahead
const bool mirror_flags_to_previous_grid
const unsigned int min_cells_for_correction
void insert_timestep(const TimeStepBase *position, TimeStepBase *new_timestep)
std::vector< std::vector< bool > > refine_flags
const double refinement_threshold
virtual void wake_up(const unsigned int wakeup_level)
std::vector< std::vector< bool > > coarsen_flags
DeclException0(ExcCantComputeTimestep)
virtual void wake_up(const unsigned int)
std::vector< std::vector< std::pair< unsigned int, double > > > CorrectionRelaxations
const bool delete_and_rebuild_tria
virtual void start_sweep(const unsigned int sweep_no)
TimeStepBase(const double time)
TimeSteppingData(const unsigned int look_ahead, const unsigned int look_back)
SmartPointer< Triangulation< dim, dim >, TimeStepBase_Tria< dim > > tria