Reference documentation for deal.II version 8.1.0
time_dependent.h
1 // ---------------------------------------------------------------------
2 // @f$Id: time_dependent.h 30036 2013-07-18 16:55:32Z maier @f$
3 //
4 // Copyright (C) 1999 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__time_dependent_h
18 #define __deal2__time_dependent_h
19 
20 
21 /*---------------------------- time-dependent.h ---------------------------*/
22 
23 
24 #include <deal.II/base/config.h>
26 #include <deal.II/base/subscriptor.h>
27 #include <deal.II/base/smartpointer.h>
28 
29 #include <vector>
30 #include <utility>
31 
33 
34 // forward declarations
35 class TimeStepBase;
36 template <typename number> class Vector;
37 template <int dim, int spacedim> class Triangulation;
38 
355 {
356 public:
368  {
374  TimeSteppingData (const unsigned int look_ahead,
375  const unsigned int look_back);
376 
406  const unsigned int look_ahead;
407 
428  const unsigned int look_back;
429  };
430 
437  {
438  forward, backward
439  };
440 
444  TimeDependent (const TimeSteppingData &data_primal,
445  const TimeSteppingData &data_dual,
446  const TimeSteppingData &data_postprocess);
447 
448 
457  virtual ~TimeDependent ();
458 
491  void insert_timestep (const TimeStepBase *position,
492  TimeStepBase *new_timestep);
493 
505  void add_timestep (TimeStepBase *new_timestep);
506 
531  void delete_timestep (const unsigned int position);
532 
544  void solve_primal_problem ();
545 
557  void solve_dual_problem ();
558 
570  void postprocess ();
571 
612  template <typename InitFunctionObject, typename LoopFunctionObject>
613  void do_loop (InitFunctionObject init_function,
614  LoopFunctionObject loop_function,
615  const TimeSteppingData &timestepping_data,
616  const Direction direction);
617 
618 
646  virtual void start_sweep (const unsigned int sweep_no);
647 
675  virtual void end_sweep ();
676 
681  virtual void end_sweep (const unsigned int n_threads) DEAL_II_DEPRECATED;
682 
688  std::size_t memory_consumption () const;
689 
693  DeclException0 (ExcInvalidPosition);
694 
695 protected:
704  std::vector<SmartPointer<TimeStepBase,TimeDependent> > timesteps;
705 
711  unsigned int sweep_no;
712 
720 
728 
736 
737 private:
738 
744  void end_sweep (const unsigned int begin_timestep,
745  const unsigned int end_timestep);
746 };
747 
748 
749 
759 class TimeStepBase : public Subscriptor
760 {
761 public:
767  {
768  primal_problem = 0x0,
769  dual_problem = 0x1,
770  postprocess = 0x2
771  };
772 
777  TimeStepBase (const double time);
778 
783  virtual ~TimeStepBase ();
784 
811  virtual void wake_up (const unsigned int);
812 
830  virtual void sleep (const unsigned int);
831 
856  virtual void start_sweep ();
857 
866  virtual void end_sweep ();
867 
879  virtual void init_for_primal_problem ();
880 
885  virtual void init_for_dual_problem ();
886 
891  virtual void init_for_postprocessing ();
892 
905  virtual void solve_primal_problem () = 0;
906 
923  virtual void solve_dual_problem ();
924 
943  virtual void postprocess_timestep ();
944 
949  double get_time () const;
950 
957  unsigned int get_timestep_no () const;
958 
983  double get_backward_timestep () const;
984 
992  double get_forward_timestep () const;
993 
1007  virtual std::size_t memory_consumption () const;
1008 
1012  DeclException0 (ExcCantComputeTimestep);
1013 
1014 protected:
1020 
1026 
1033  unsigned int sweep_no;
1034 
1043  unsigned int timestep_no;
1044 
1048  const double time;
1049 
1059  unsigned int next_action;
1060 
1061 private:
1071  void set_previous_timestep (const TimeStepBase *previous);
1072 
1082  void set_next_timestep (const TimeStepBase *next);
1083 
1094  void set_timestep_no (const unsigned int step_no);
1095 
1103  void set_sweep_no (const unsigned int sweep_no);
1104 
1105 
1115  TimeStepBase (const TimeStepBase &);
1116 
1127 
1128  // make the manager object a friend
1129  friend class TimeDependent;
1130 };
1131 
1132 
1133 
1134 
1144 {
1152  template <int dim>
1153  struct Flags
1154  {
1160  Flags ();
1161 
1167  Flags (const bool delete_and_rebuild_tria,
1168  const unsigned int wakeup_level_to_build_grid,
1169  const unsigned int sleep_level_to_delete_grid);
1170 
1194 
1212  const unsigned int wakeup_level_to_build_grid;
1213 
1220  const unsigned int sleep_level_to_delete_grid;
1221 
1225  DeclException1 (ExcInvalidParameter,
1226  int,
1227  << "The parameter " << arg1 << " has an invalid value.");
1228  };
1229 
1230 
1231 
1352  template <int dim>
1354  {
1363  typedef std::vector<std::vector<std::pair<unsigned int, double> > > CorrectionRelaxations;
1364 
1369  static CorrectionRelaxations default_correction_relaxations;
1370 
1377  RefinementFlags (const unsigned int max_refinement_level = 0,
1378  const unsigned int first_sweep_with_correction = 0,
1379  const unsigned int min_cells_for_correction = 0,
1380  const double cell_number_corridor_top = (1<<dim),
1381  const double cell_number_corridor_bottom = 1,
1382  const CorrectionRelaxations &correction_relaxations = CorrectionRelaxations(),
1383  const unsigned int cell_number_correction_steps = 0,
1384  const bool mirror_flags_to_previous_grid = false,
1385  const bool adapt_grids = false);
1386 
1404  const unsigned int max_refinement_level;
1405 
1414  const unsigned int first_sweep_with_correction;
1415 
1416 
1423  const unsigned int min_cells_for_correction;
1424 
1434 
1439 
1444  const std::vector<std::vector<std::pair<unsigned int,double> > > correction_relaxations;
1445 
1454  const unsigned int cell_number_correction_steps;
1455 
1490 
1495  const bool adapt_grids;
1496 
1500  DeclException1 (ExcInvalidValue,
1501  double,
1502  << "The following value does not fulfill the requirements: " << arg1);
1503  };
1504 
1505 
1506 
1513  template <int dim>
1515  {
1519  RefinementData (const double refinement_threshold,
1520  const double coarsening_threshold=0);
1521 
1533  const double refinement_threshold;
1534 
1541  const double coarsening_threshold;
1542 
1546  DeclException1 (ExcInvalidValue,
1547  double,
1548  << "The following value does not fulfill the requirements: " << arg1);
1549  };
1550 }
1551 
1552 
1553 
1554 
1572 template <int dim>
1574 {
1575 public:
1584 
1585 
1592  {
1593  grid_refinement = 0x1000
1594  };
1595 
1596 
1610  TimeStepBase_Tria ();
1611 
1634  TimeStepBase_Tria (const double time,
1636  const Flags &flags,
1637  const RefinementFlags &refinement_flags = RefinementFlags());
1638 
1645  virtual ~TimeStepBase_Tria ();
1646 
1681  virtual void wake_up (const unsigned int wakeup_level);
1682 
1708  virtual void sleep (const unsigned int);
1709 
1732  void refine_grid (const RefinementData data);
1733 
1742  virtual void init_for_refinement ();
1743 
1755  virtual void get_tria_refinement_criteria (Vector<float> &criteria) const = 0;
1756 
1764  void save_refine_flags ();
1765 
1779  virtual std::size_t memory_consumption () const;
1780 
1784  DeclException0 (ExcGridNotDeleted);
1785 
1786 protected:
1787 
1803 
1814 
1821  const Flags flags;
1822 
1829  const RefinementFlags refinement_flags;
1830 
1831 private:
1839  std::vector<std::vector<bool> > refine_flags;
1840 
1844  std::vector<std::vector<bool> > coarsen_flags;
1845 
1852  void restore_grid ();
1853 };
1854 
1855 
1856 
1857 
1858 
1859 /*----------------------------- template functions ------------------------------*/
1860 
1861 template <typename InitFunctionObject, typename LoopFunctionObject>
1862 void TimeDependent::do_loop (InitFunctionObject init_function,
1863  LoopFunctionObject loop_function,
1864  const TimeSteppingData &timestepping_data,
1865  const Direction direction)
1866 {
1867  // the following functions looks quite
1868  // disrupted due to the recurring switches
1869  // for forward and backward running loops.
1870  //
1871  // I chose to switch at every place where
1872  // it is needed, since it is so easy
1873  // to overlook something when you change
1874  // some code at one place when it needs
1875  // to be changed at a second place, here
1876  // for the other direction, also.
1877 
1878  const unsigned int n_timesteps = timesteps.size();
1879 
1880  // initialize the time steps for
1881  // a round of this loop
1882  for (unsigned int step=0; step<n_timesteps; ++step)
1883  switch (direction)
1884  {
1885  case forward:
1886  init_function (static_cast<typename InitFunctionObject::argument_type>
1887  (&*timesteps[step]));
1888  break;
1889  case backward:
1890  init_function (static_cast<typename InitFunctionObject::argument_type>
1891  (&*timesteps[n_timesteps-@ref step_1 "step-1"]));
1892  break;
1893  };
1894 
1895 
1896  // wake up the first few time levels
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)
1900  switch (direction)
1901  {
1902  case forward:
1903  if (step+look_ahead >= 0)
1904  timesteps[step+look_ahead]->wake_up(look_ahead);
1905  break;
1906  case backward:
1907  if (n_timesteps-(step+look_ahead) < n_timesteps)
1908  timesteps[n_timesteps-(step+look_ahead)]->wake_up(look_ahead);
1909  break;
1910  };
1911 
1912 
1913  for (unsigned int step=0; step<n_timesteps; ++step)
1914  {
1915  // first thing: wake up the
1916  // timesteps ahead as necessary
1917  for (unsigned int look_ahead=0;
1918  look_ahead<=timestepping_data.look_ahead; ++look_ahead)
1919  switch (direction)
1920  {
1921  case forward:
1922  if (step+look_ahead < n_timesteps)
1923  timesteps[step+look_ahead]->wake_up(look_ahead);
1924  break;
1925  case backward:
1926  if (n_timesteps > (step+look_ahead))
1927  timesteps[n_timesteps-(step+look_ahead)-1]->wake_up(look_ahead);
1928  break;
1929  };
1930 
1931 
1932  // actually do the work
1933  switch (direction)
1934  {
1935  case forward:
1936  loop_function (static_cast<typename LoopFunctionObject::argument_type>
1937  (&*timesteps[step]));
1938  break;
1939  case backward:
1940  loop_function (static_cast<typename LoopFunctionObject::argument_type>
1941  (&*timesteps[n_timesteps-@ref step_1 "step-1"]));
1942  break;
1943  };
1944 
1945  // let the timesteps behind sleep
1946  for (unsigned int look_back=0;
1947  look_back<=timestepping_data.look_back; ++look_back)
1948  switch (direction)
1949  {
1950  case forward:
1951  if (step>=look_back)
1952  timesteps[step-look_back]->sleep(look_back);
1953  break;
1954  case backward:
1955  if (n_timesteps-(step-look_back) <= n_timesteps)
1956  timesteps[n_timesteps-(step-look_back)-1]->sleep(look_back);
1957  break;
1958  };
1959  };
1960 
1961  // make the last few timesteps sleep
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)
1966  switch (direction)
1967  {
1968  case forward:
1969  if ((step-look_back >= 0)
1970  &&
1971  (step-look_back < static_cast<int>(n_timesteps)))
1972  timesteps[step-look_back]->sleep(look_back);
1973  break;
1974  case backward:
1975  if ((step-look_back >= 0)
1976  &&
1977  (step-look_back < static_cast<int>(n_timesteps)))
1978  timesteps[n_timesteps-(step-look_back)-1]->sleep(look_back);
1979  break;
1980  };
1981 }
1982 
1983 DEAL_II_NAMESPACE_CLOSE
1984 
1985 /*---------------------------- time-dependent.h ---------------------------*/
1986 #endif
1987 /*---------------------------- time-dependent.h ---------------------------*/
void set_sweep_no(const unsigned int sweep_no)
virtual void solve_primal_problem()=0
std::size_t memory_consumption() const
virtual void start_sweep()
void refine_grid(const RefinementData data)
unsigned int next_action
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)
virtual void end_sweep()
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
DeclException0(ExcInvalidPosition)
virtual ~TimeStepBase()
void set_previous_timestep(const TimeStepBase *previous)
unsigned int timestep_no
void do_loop(InitFunctionObject init_function, LoopFunctionObject loop_function, const TimeSteppingData &timestepping_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 ~TimeDependent()
double get_time() const
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
const double time
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 sleep_level_to_delete_grid
DeclException1(ExcInvalidValue, double,<< "The following value does not fulfill the requirements: "<< arg1)
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
virtual void end_sweep()
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
void save_refine_flags()
unsigned int sweep_no
void insert_timestep(const TimeStepBase *position, TimeStepBase *new_timestep)
std::vector< std::vector< bool > > refine_flags
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)
void postprocess()
std::vector< std::vector< std::pair< unsigned int, double > > > CorrectionRelaxations
unsigned int sweep_no
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