Reference documentation for deal.II version Git 7026f387cc 2019-10-15 14:19:01 -0400
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
theta_timestepping.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2010 - 2018 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 
17 #ifndef dealii_theta_timestepping_h
18 #define dealii_theta_timestepping_h
19 
20 #include <deal.II/algorithms/operator.h>
21 #include <deal.II/algorithms/timestep_control.h>
22 
23 #include <deal.II/base/smartpointer.h>
24 
25 DEAL_II_NAMESPACE_OPEN
26 
27 // Forward declaration
28 #ifndef DOXYGEN
29 class ParameterHandler;
30 #endif
31 
32 namespace Algorithms
33 {
43  struct TimestepData
44  {
46  double time;
48  double step;
49  };
50 
202  template <typename VectorType>
204  {
205  public:
211  ThetaTimestepping(OperatorBase &op_explicit, OperatorBase &op_implicit);
212 
224  virtual void
225  operator()(AnyData &out, const AnyData &in) override;
226 
230  virtual void
231  notify(const Event &) override;
232 
237  void
238  set_output(OutputOperator<VectorType> &output);
239 
243  static void
244  declare_parameters(ParameterHandler &param);
245 
249  void
250  parse_parameters(ParameterHandler &param);
251 
255  double
256  current_time() const;
257 
261  double
262  theta() const;
263 
267  double
268  theta(double new_theta);
269 
276  const TimestepData &
277  explicit_data() const;
278 
285  const TimestepData &
286  implicit_data() const;
287 
292  timestep_control();
293 
294  private:
300 
305  double vtheta;
309  bool adaptive;
310 
315 
320 
321 
333 
345 
351  };
352 
353 
354  template <typename VectorType>
355  inline const TimestepData &
357  {
358  return d_explicit;
359  }
360 
361 
362  template <typename VectorType>
363  inline const TimestepData &
365  {
366  return d_implicit;
367  }
368 
369 
370  template <typename VectorType>
371  inline TimestepControl &
373  {
374  return control;
375  }
376 
377  template <typename VectorType>
378  inline void
380  {
381  output = &out;
382  }
383 
384 
385  template <typename VectorType>
386  inline double
388  {
389  return vtheta;
390  }
391 
392 
393  template <typename VectorType>
394  inline double
396  {
397  const double tmp = vtheta;
398  vtheta = new_theta;
399  return tmp;
400  }
401 
402 
403  template <typename VectorType>
404  inline double
406  {
407  return control.now();
408  }
409 } // namespace Algorithms
410 
411 DEAL_II_NAMESPACE_CLOSE
412 
413 #endif
TimestepControl & timestep_control()
const TimestepData & explicit_data() const
double time
The current time.
SmartPointer< OperatorBase, ThetaTimestepping< VectorType > > op_implicit
SmartPointer< OutputOperator< VectorType >, ThetaTimestepping< VectorType > > output
void set_output(OutputOperator< VectorType > &output)
double step
The current step size times something.
SmartPointer< OperatorBase, ThetaTimestepping< VectorType > > op_explicit
const TimestepData & implicit_data() const