using namespace Algorithms;
-class Explicit
- : public Operator<Vector<double> >
+class Explicit : public Operator<Vector<double> >
{
public:
Explicit(const FullMatrix<double> &matrix);
- void operator() (NamedData<Vector<double>*> &out,
- const NamedData<Vector<double>*> &in);
+ void operator() (AnyData &out, const AnyData &in);
- void initialize_timestep_data(const TimestepData &);
private:
- const TimestepData *timestep_data;
SmartPointer<const FullMatrix<double>, Explicit> matrix;
FullMatrix<double> m;
};
-class Implicit
- : public Operator<Vector<double> >
+class Implicit : public Operator<Vector<double> >
{
public:
Implicit(const FullMatrix<double> &matrix);
- void operator() (NamedData<Vector<double>*> &out,
- const NamedData<Vector<double>*> &in);
-
- void initialize_timestep_data(const TimestepData &);
-private:
- const TimestepData *timestep_data;
+ void operator() (AnyData &out, const AnyData &in);
+
+ private:
SmartPointer<const FullMatrix<double>, Implicit> matrix;
FullMatrix<double> m;
};
int main()
{
FullMatrix<double> matrix(2);
- matrix(0,0) = 1.;
- matrix(1,1) = 1.;
- matrix(0,1) = 31.4;
- matrix(1,0) = -31.4;
+ matrix(0,0) = 0.;
+ matrix(1,1) = 0.;
+ matrix(0,1) = 3.14;
+ matrix(1,0) = -3.14;
OutputOperator<Vector<double> > out;
out.initialize_stream(std::cout);
Explicit op_explicit(matrix);
Implicit op_implicit(matrix);
ThetaTimestepping<Vector<double> > solver(op_explicit, op_implicit);
- op_explicit.initialize_timestep_data(solver.explicit_data());
- op_implicit.initialize_timestep_data(solver.implicit_data());
solver.set_output(out);
Vector<double> value(2);
value(0) = 1.;
- NamedData<Vector<double>*> indata;
- NamedData<Vector<double>*> outdata;
- Vector<double> *p = &value;
- outdata.add(p, "value");
+ AnyData indata;
+ AnyData outdata;
+ outdata.add(&value, "value");
solver.notify(Events::initial);
solver(outdata, indata);
void
-Explicit::initialize_timestep_data(const TimestepData &t)
-{
- timestep_data = &t;
-}
-
-
-void
-Explicit::operator() (NamedData<Vector<double>*> &out, const NamedData<Vector<double>*> &in)
+Explicit::operator() (AnyData &out, const AnyData &in)
{
+ const double timestep = *in.read_ptr<double>("Timestep");
if (this->notifications.test(Events::initial) || this->notifications.test(Events::new_timestep_size))
{
- m.equ(-timestep_data->step, *matrix);
+ m.equ(-timestep, *matrix);
for (unsigned int i=0; i<m.m(); ++i)
m(i,i) += 1.;
}
this->notifications.clear();
- unsigned int i = in.find("Previous iterate");
- m.vmult(*out(0), *in(i));
+ m.vmult(*out.entry<Vector<double>*>(0),
+ *in.read_ptr<Vector<double> >("Previous iterate"));
}
void
-Implicit::initialize_timestep_data(const TimestepData &t)
-{
- timestep_data = &t;
-}
-
-
-void
-Implicit::operator() (NamedData<Vector<double>*> &out, const NamedData<Vector<double>*> &in)
+Implicit::operator() (AnyData &out, const AnyData &in)
{
+ const double timestep = *in.read_ptr<double>("Timestep");
if (this->notifications.test(Events::initial) || this->notifications.test(Events::new_timestep_size))
{
- m.equ(timestep_data->step, *matrix);
+ m.equ(timestep, *matrix);
for (unsigned int i=0; i<m.m(); ++i)
m(i,i) += 1.;
m.gauss_jordan();
}
this->notifications.clear();
-
- unsigned int i = in.find("Previous time");
- m.vmult(*out(0), *in(i));
+ m.vmult(*out.entry<Vector<double>*>(0),
+ *in.read_ptr<Vector<double> >("Previous time"));
}
* @until End of declarations
*
* These operators will be implemented after the main program. But let
- * us look at how they get used. First, let us define a matrix to be
+ * us look first at how they get used. First, let us define a matrix to be
* used for our system and also an OutputOperator in order to write
* the data of each timestep to a file.
*
* @until out.initialize
*
* Now we create objects for the implicit and explicit parts of the
- * steps as well as the ThetaTimestepping itself. Notice how the
- * TimestepData of ThetaTimestepping gets forwarded to the inner
- * operators. There are two different data objects, because the
- * timestep size is modified by #theta.
+ * steps as well as the ThetaTimestepping itself. We initialize the
+ * timestepping with the output operator in order to be able to see
+ * the output in every step.
*
* @until set_output
*
* is filled with the initial value and is also the vector where the
* solution at each timestep will be. Because the interface of
* Operator has to be able to handle several vectors, we need to store
- * it in a NamedData object. Notice, that we need to create the
- * intermediate pointer <tt>p</tt>. If we would use
- * <code>&value</code> directly in the <code>add</code> function, the
- * resulting object would be constant.
+ * it in an AnyData object. Since our problem has no additional
+ * parameters, the input AnyData object remains empty.
*
* @until add
*
- * Finally, we are ready to tell the solver, that we are looknig at
+ * Finally, we are ready to tell the solver, that we are starting at
* the initial timestep and run it.
*
- * @until outdata
- * @skip Explicit::initialize
+ * @until }
+ *
+ * First the constructor, which simply copies the system matrix into
+ * the member pointer for later use.
+ *
+ * @skip Explicit::
+ * @until }
*
* Now we need to study the application of the implicit and explicit
* operator. We assume that the pointer <code>matrix</code> points to
- * the matrix created in the main program, and that
- * <code>timestep_data</code> points to the correct data object of
- * ThetaTimestepping.
+ * the matrix created in the main program (the constructor did this
+ * for us). Here, we first get the time step size from the AnyData
+ * object that was provided as input. Then, if we are in the first
+ * step or if the timestep has changed, we fill the local matrix
+ * <tt>m</tt>, such that with the given matrix \f$M\f$, it becomes
+ * \f[
+ * m = I - \Delta t M.
+ * \f]
+ * After we have worked off the notifications, we clear them, such
+ * that the matrix is only generated when necessary.
*
* @skipline void
- * @until vmult
+ * @until clear
+ *
+ * Now we multiply the input vector with the new matrix and store on output.
+ *
* @until }
+ * The code for the implicit operator is almost the same, except
+ * that we change the sign in front of the timestep and use the
+ * inverse of t he matrix.
*
+ * @until vmult
+ * @until }
* @author Guido Kanschat
* @date 2010
*/
projects / dealii-svn.git / commitdiff