Reference documentation for deal.II version 8.1.0
theta_timestepping.templates.h
1 // ---------------------------------------------------------------------
2 // @f$Id: theta_timestepping.templates.h 30036 2013-07-18 16:55:32Z maier @f$
3 //
4 // Copyright (C) 2006 - 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 
18 #include <deal.II/algorithms/theta_timestepping.h>
19 
20 #include <deal.II/base/parameter_handler.h>
21 #include <deal.II/lac/vector_memory.h>
22 
24 
25 namespace Algorithms
26 {
27  template <class VECTOR>
29  : vtheta(0.5), adaptive(false), op_explicit(&e), op_implicit(&i)
30  {}
31 
32 
33  template <class VECTOR>
34  void
36  {
37  op_explicit->notify(e);
38  op_implicit->notify(e);
39  }
40 
41  template <class VECTOR>
42  void
44  {
45  param.enter_subsection("ThetaTimestepping");
47  param.declare_entry("Theta", ".5", Patterns::Double());
48  param.declare_entry("Adaptive", "false", Patterns::Bool());
49  param.leave_subsection();
50  }
51 
52  template <class VECTOR>
53  void
55  {
56  param.enter_subsection("ThetaTimestepping");
57  control.parse_parameters (param);
58  vtheta = param.get_double("Theta");
59  adaptive = param.get_bool("Adaptive");
60  param.leave_subsection ();
61  }
62 
63 
64  template <class VECTOR>
65  void
67  {
69 
70  deallog.push ("Theta");
72  typename VectorMemory<VECTOR>::Pointer aux(mem);
73  aux->reinit(*out(0));
74 
75  control.restart();
76 
78 
79  // The data used to compute the
80  // vector associated with the old
81  // timestep
82  VECTOR *p = out(0);
84  src1.add(p, "Previous iterate");
85  src1.merge(in);
86 
88  src2.add(p, "Previous iterate");
89 
90  p = aux;
92  out1.add(p, "Result");
93  // The data provided to the inner
94  // solver
95  src2.add(p, "Previous time");
96  src2.merge(in);
97 
98  if (output != 0)
99  (*output) << 0U << out;
100 
101  for (unsigned int count = 1; d_explicit.time < control.final(); ++count)
102  {
103  const bool step_change = control.advance();
105  d_explicit.step = (1.-vtheta)*control.step();
107  deallog << "Time step:" << d_implicit.time << std::endl;
108 
109  op_explicit->notify(Events::new_time);
110  op_implicit->notify(Events::new_time);
111  if (step_change)
112  {
115  }
116 
117  // Compute
118  // (I + (1-theta)dt A) u
119  (*op_explicit)(out1, src1);
120  (*op_implicit)(out, src2);
121 
122  if (output != 0 && control.print())
123  (*output) << count << out;
124 
126  }
127  deallog.pop();
128  }
129 }
130 
131 DEAL_II_NAMESPACE_CLOSE
void pop()
void add(DATA &v, const std::string &name)
Definition: named_data.h:253
void parse_parameters(ParameterHandler &param)
ThetaTimestepping(Operator< VECTOR > &op_explicit, Operator< VECTOR > &op_implicit)
void enter_subsection(const std::string &subsection)
virtual void operator()(NamedData< VECTOR * > &out, const NamedData< VECTOR * > &in)
#define Assert(cond, exc)
Definition: exceptions.h:299
double time
The current time.
SmartPointer< Operator< VECTOR >, ThetaTimestepping< VECTOR > > op_implicit
const Event new_timestep_size
SmartPointer< OutputOperator< VECTOR >, ThetaTimestepping< VECTOR > > output
void merge(NamedData< DATA2 > &)
Definition: named_data.h:278
bool get_bool(const std::string &entry_name) const
double step
The current step size times something.
static void declare_parameters(ParameterHandler &param)
void push(const std::string &text)
bool leave_subsection()
double get_double(const std::string &entry_name) const
void declare_entry(const std::string &entry, const std::string &default_value, const Patterns::PatternBase &pattern=Patterns::Anything(), const std::string &documentation=std::string())
::ExceptionBase & ExcNotImplemented()
const Event new_time
SmartPointer< Operator< VECTOR >, ThetaTimestepping< VECTOR > > op_explicit