--- /dev/null
+// ---------------------------------------------------------------------
+// $Id: named_data.h 30036 2013-07-18 16:55:32Z maier $
+//
+// Copyright (C) 2000 - 2013 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+#ifndef __deal2__any_data_h
+#define __deal2__any_data_h
+
+#include <deal.II/base/config.h>
+#include <deal.II/base/exceptions.h>
+#include <deal.II/base/subscriptor.h>
+
+#include <boost/any.hpp>
+#include <vector>
+#include <algorithm>
+#include <typeinfo>
+
+DEAL_II_NAMESPACE_OPEN
+
+/**
+ * Store any amount of any type of data accessible by an identifier string.
+ *
+ * @todo Deprecate access by index after NamedData has been deprecated
+ * for long enough, then change to a map.
+ */
+class AnyData :
+ public Subscriptor
+{
+ public:
+ /// Number of stored data objects.
+ unsigned int size() const;
+
+ /// Add a new data object
+ template <typename type>
+ void add(type entry, const std::string& name);
+
+ /**
+ * @brief Merge the data of another NamedData to the end of this object.
+ */
+ void merge(const AnyData& other);
+
+ /**
+ * @brief Access to stored data object by name.
+ */
+ template <typename type>
+ type entry (const std::string& name);
+
+ /// Read-only access to stored data object by name.
+ template <typename type>
+ const type entry (const std::string& name) const;
+
+ /// Dedicated read only access by name.
+ template <typename type>
+ const type read (const std::string& name) const;
+
+ /**
+ * @brief Access to stored data object by index.
+ */
+ template <typename type>
+ type entry (const unsigned int i);
+
+ /// Read-only access to stored data object by index.
+ template <typename type>
+ const type entry (const unsigned int i) const;
+
+ /// Dedicated read only access.
+ template <typename type>
+ const type read (const unsigned int i) const;
+
+ /// Name of object at index.
+ const std::string &name(const unsigned int i) const;
+
+ /// Find index of a named object
+ unsigned int find(const std::string &name) const;
+
+ /// Find out if object is of a certain type
+ template <typename type>
+ bool is_type(const unsigned int i) const;
+
+ /// The requested type and the stored type are different
+ DeclException2(ExcTypeMismatch,
+ char*, char*,
+ << "The requested type " << arg1
+ << " and the stored type " << arg2
+ << " must coincide");
+
+ private:
+ /// The stored data
+ std::vector<boost::any> data;
+ /// The names of the stored data
+ std::vector<std::string> names;
+};
+
+
+unsigned int
+inline
+AnyData::size () const
+{
+ AssertDimension(data.size(), names.size());
+ return data.size();
+}
+
+
+template <typename type>
+inline
+type
+AnyData::entry (const unsigned int i)
+{
+ AssertIndexRange(i, size());
+ Assert(data[i].type() == typeid(type),
+ ExcTypeMismatch(typeid(type).name(),data[i].type().name()));
+ return boost::any_cast<type>(data[i]);
+}
+
+
+template <typename type>
+inline
+const type
+AnyData::entry (const unsigned int i) const
+{
+ AssertIndexRange(i, size());
+ const type* p = boost::any_cast<type>(&data[i]);
+ Assert(p != 0,
+ ExcTypeMismatch(typeid(type).name(),data[i].type().name()));
+ return *p;
+}
+
+
+template <typename type>
+inline
+const type
+AnyData::read(const unsigned int i) const
+{
+ AssertIndexRange(i, size());
+ Assert(data[i].type() == typeid(type),
+ ExcTypeMismatch(typeid(type).name(),data[i].type().name()));
+ return boost::any_cast<const type>(data[i]);
+}
+
+
+inline
+const std::string&
+AnyData::name(const unsigned int i) const
+{
+ AssertIndexRange(i, size());
+ return names[i];
+}
+
+
+inline
+unsigned int
+AnyData::find(const std::string& n) const
+{
+ typename std::vector<std::string>::const_iterator it =
+ std::find(names.begin(), names.end(), n);
+
+ Assert(it != names.end(), ExcMessage("An entry with this name does not exist"));
+
+ return it - names.begin();
+}
+
+
+template <typename type>
+inline
+bool
+AnyData::is_type(const unsigned int i) const
+{
+ return data[i].type() == typeid(type);
+}
+
+
+template <typename type>
+inline
+type
+AnyData::entry (const std::string& n)
+{
+ const unsigned int i = find(n);
+ /* Assert(dynamic_cast<type*>(&data[i]) != 0, */
+ /* ExcTypeMismatch(typeid(type).name(),data[i].type().name())); */
+ return boost::any_cast<type>(data[i]);
+}
+
+
+template <typename type>
+inline
+const type
+AnyData::entry (const std::string& n) const
+{
+ const unsigned int i = find(n);
+ Assert(data[i].type() == typeid(type),
+ ExcTypeMismatch(typeid(type).name(),data[i].type().name()));
+ return boost::any_cast<const type>(data[i]);
+}
+
+
+template <typename type>
+inline
+const type
+AnyData::read(const std::string& n) const
+{
+ const unsigned int i = find(n);
+ Assert(data[i].type() == typeid(type),
+ ExcTypeMismatch(typeid(type).name(),data[i].type().name()));
+ return boost::any_cast<const type>(data[i]);
+}
+
+
+template <typename type>
+inline
+void
+AnyData::add(type ent, const std::string& n)
+{
+ boost::any e = ent;
+ data.push_back(e);
+ names.push_back(n);
+}
+
+
+inline
+void
+AnyData::merge(const AnyData& other)
+{
+ for (unsigned int i=0; i<other.size(); ++i)
+ {
+ names.push_back(other.names[i]);
+ data.push_back(other.data[i]);
+ }
+}
+
+
+//----------------------------------------------------------------------//
+
+
+
+DEAL_II_NAMESPACE_CLOSE
+
+#endif
#include <deal.II/base/smartpointer.h>
#include <deal.II/lac/solver_control.h>
#include <deal.II/algorithms/operator.h>
+#include <deal.II/algorithms/any_data.h>
DEAL_II_NAMESPACE_OPEN
* #inverse_derivative. It is up to this object to implement
* reassembling accordingly.
*
- * <h3>Contents of the NamedData objects</h3>
+ * <h3>Contents of the AnyData objects</h3>
*
* The only value used by the Newton method is the first vector in the
* parameter <tt>out</tt> of operator()(). It serves as the start
* down to the objects
* #residual and #inverse_derivative.
*/
- virtual void operator() (NamedData<VECTOR *> &out, const NamedData<VECTOR *> &in);
+ virtual void operator() (AnyData &out, const AnyData &in);
+
+ /**
+ * @deprecated Use the function using AnyData
+ */
+ virtual void operator() (NamedData<VECTOR*> &out, const NamedData<VECTOR*> &in);
virtual void notify(const Event &);
template <class VECTOR>
void
Newton<VECTOR>::operator() (NamedData<VECTOR *> &out, const NamedData<VECTOR *> &in)
+ {
+ Operator<VECTOR>::operator() (out, in);
+ }
+
+ template <class VECTOR>
+ void
+ Newton<VECTOR>::operator() (AnyData &out, const AnyData &in)
{
Assert (out.size() == 1, ExcNotImplemented());
deallog.push ("Newton");
- VECTOR &u = *out(0);
+ VECTOR &u = *out.entry<VECTOR*>(0);
if (debug>2)
deallog << "u: " << u.l2_norm() << std::endl;
typename VectorMemory<VECTOR>::Pointer res(mem);
res->reinit(u);
- NamedData<VECTOR *> src1;
- NamedData<VECTOR *> src2;
- VECTOR *p = &u;
- src1.add(p, "Newton iterate");
+ AnyData src1;
+ AnyData src2;
+ src1.add<const VECTOR*>(&u, "Newton iterate");
src1.merge(in);
- p = res;
- src2.add(p, "Newton residual");
+ src2.add<const VECTOR*>(res, "Newton residual");
src2.merge(src1);
- NamedData<VECTOR *> out1;
- out1.add(p, "Residual");
- p = Du;
- NamedData<VECTOR *> out2;
- out2.add(p, "Update");
+ AnyData out1;
+ out1.add<VECTOR*>(res, "Residual");
+ AnyData out2;
+ out2.add<VECTOR*>(Du, "Update");
unsigned int step = 0;
// fill res with (f(u), v)
#define __deal2__operator_h
#include <deal.II/base/config.h>
+#include <deal.II/algorithms/any_data.h>
#include <deal.II/base/named_data.h>
#include <deal.II/base/event.h>
namespace Algorithms
{
/**
+ * @todo Update this documentation and the one of Operator
+ *
* The abstract base class of all algorithms in this library. An
* operator is an object with an operator(), which transforms a set
* of named vectors into another set of named vectors.
* the innermost method in such a nested system will have to compute a
* residual using values from all outer iterations. Since the depth
* and order of such a nesting is hardly predictable when designing a
- * general tool, we use NamedData to access these vectors. Typically,
+ * general tool, we use AnyData to access these vectors. Typically,
* the first vector in <tt>out</tt> contains the start vector when
* operator()() is called, and the solution when the function
* returns. The object <tt>in</tt> is providing additional information
* and forwarded to the inner Operator objects of the nested
* iteration.
*
- * @author Guido Kanschat, 2010
+ * @author Guido Kanschat
+ * @date 2014
*/
- template <class VECTOR>
- class Operator : public Subscriptor
+ class OperatorBase : public Subscriptor
{
public:
/**
* The virtual destructor.
*/
- ~Operator();
+ ~OperatorBase();
/**
* The actual operation, which
* is implemented in a derived class.
*/
- virtual void operator() (NamedData<VECTOR *> &out, const NamedData<VECTOR *> &in) = 0;
+ virtual void operator() (AnyData &out, const AnyData &in) = 0;
/**
* Register an event triggered
void clear_events();
protected:
/**
- * Accumulate reasons for
- * reassembling here. If any of
+ * Accumulate events here. If any of
* those is set, the function
* solve() of a terminal
* application must take care
* of reassembling the matrix.
*/
Event notifications;
+
+ };
+
+ /**
+ * @deprecated This class has been replaced by OperatorBase.
+ *
+ * The abstract base class of all algorithms in this library. An
+ * operator is an object with an operator(), which transforms a set
+ * of named vectors into another set of named vectors.
+ *
+ * Furthermore, an operator can be notified of parameter changes by
+ * the calling routine. The outer iteration can notify() the Operator
+ * of an Event, which could be for instance a change of mesh, a
+ * different time step size or too slow convergence of Newton's
+ * method, which would then trigger reassembling of a matrix or
+ * similar things.
+ *
+ * <h3>Usage for nested iterations</h3>
+ *
+ * This is probably the most prominent use for Operator, where an
+ * outer iterative method calls an inner solver and so on. Typically,
+ * the innermost method in such a nested system will have to compute a
+ * residual using values from all outer iterations. Since the depth
+ * and order of such a nesting is hardly predictable when designing a
+ * general tool, we use NamedData to access these vectors. Typically,
+ * the first vector in <tt>out</tt> contains the start vector when
+ * operator()() is called, and the solution when the function
+ * returns. The object <tt>in</tt> is providing additional information
+ * and forwarded to the inner Operator objects of the nested
+ * iteration.
+ *
+ * @author Guido Kanschat, 2010
+ */
+ template <class VECTOR>
+ class Operator : public OperatorBase
+ {
+ public:
+ Operator();
+
+ /**
+ * Implementation of the function in the base class in order to do
+ * compatibility conversions between the old and the new
+ * interface.
+ */
+ virtual void operator() (AnyData &out, const AnyData &in);
+
+ /**
+ * @deprecated It is in particular this function which should not be used anymore.
+ *
+ * The actual operation, which
+ * is implemented in a derived class.
+ */
+ virtual void operator() (NamedData<VECTOR *> &out, const NamedData<VECTOR *> &in);
+
+ private:
+ /**
+ * While we are providing compatibility functions to the old
+ * interface, this variable will ensure there is no endless loop.
+ */
+ bool compatibility_flag;
};
/**
namespace Algorithms
{
template <class VECTOR>
- Operator<VECTOR>::~Operator()
+ Operator<VECTOR>::Operator()
+ : compatibility_flag(false)
{}
-
template <class VECTOR>
- void Operator<VECTOR>::notify(const Event &e)
+ void
+ Operator<VECTOR>::operator() (AnyData& out, const AnyData& in)
{
- notifications += e;
+ // Had this function been overloaded in a derived clas, it would
+ // not have been called. Therefore, we have to start the
+ // compatibility engine. But before, we have to avoid an endless loop.
+ Assert(!compatibility_flag, ExcMessage("Compatibility resolution of Operator generates and endless loop\n"
+ "Please provide an operator() in a derived class"));
+ compatibility_flag = true;
+
+ NamedData<VECTOR*> new_out;
+ for (unsigned int i=0;i<out.size();++i)
+ {
+ if (out.is_type<VECTOR*>(i))
+ new_out.add(out.entry<VECTOR*>(i), out.name(i));
+ else
+ deallog << "Cannot convert AnyData argument " << out.name(i) << " to NamedData"
+ << std::endl;
+ }
+
+ NamedData<VECTOR*> new_in;
+ for (unsigned int i=0;i<in.size();++i)
+ {
+ // deallog << "Convert " << in.name(i) << std::endl;
+ if (in.is_type<VECTOR*>(i))
+ {
+ // This const cast is due to the wrong constness handling
+ // in NamedData. As soon as NamedData is gone, this code
+ // will not be necessary anymore. And deprecating begins
+ // now.
+ VECTOR* p = const_cast<VECTOR*> (in.entry<VECTOR*>(i));
+ new_in.add(p, in.name(i));
+ }
+ else if (in.is_type<const VECTOR*>(i))
+ {
+ // This const cast is due to the wrong constness handling
+ // in NamedData. As soon as NamedData is gone, this code
+ // will not be necessary anymore. And deprecating begins
+ // now.
+ VECTOR* p = const_cast<VECTOR*> (in.entry<const VECTOR*>(i));
+ new_in.add(p, in.name(i));
+ }
+ else
+ deallog << "Cannot convert AnyData argument " << in.name(i)
+ << " to NamedData" << std::endl;
+ }
+ this->operator() (new_out, new_in);
+ compatibility_flag = false;
}
template <class VECTOR>
void
- Operator<VECTOR>::clear_events ()
+ Operator<VECTOR>::operator() (NamedData<VECTOR*>& out, const NamedData<VECTOR*>& in)
{
- notifications.clear();
+ // Had this function been overloaded in a derived clas, it would
+ // not have been called. Therefore, we have to start the
+ // compatibility engine. But before, we have to avoid an endless loop.
+ Assert(!compatibility_flag, ExcMessage("Compatibility resolution of Operator generates and endless loop\n"
+ "Please provide an operator() in a derived class"));
+ compatibility_flag = true;
+
+ AnyData new_out;
+ for (unsigned int i=0;i<out.size();++i)
+ new_out.add(out(i), out.name(i));
+
+ AnyData new_in;
+ for (unsigned int i=0;i<in.size();++i)
+ new_in.add(in(i), in.name(i));
+
+ this->operator() (new_out, new_in);
+ compatibility_flag = false;
}
/**
+ * @deprecated The use of this class is deprecated and AnyData should
+ * be used instead. It is not only more flexible, the way constness
+ * was assured in this class never worked the way it was
+ * intended. Therefore, the according checks have been disabled.
+ *
* This class is a collection of DATA objects. Each of the pointers
* has a name associated, enabling identification by this name rather
* than the index in the vector only.
void
NamedData<DATA>::add(DATA &v, const std::string &n)
{
- Assert(!is_constant, ExcConstantObject());
+ // see operator() below
+
+ // Assert(!is_constant, ExcConstantObject());
names.push_back(n);
data.push_back(v);
}
void
NamedData<DATA>::add(const DATA &v, const std::string &n)
{
- Assert(!is_constant, ExcConstantObject());
+ // see operator() below
+
+ // Assert(!is_constant, ExcConstantObject());
DATA &aux = const_cast<DATA &>(v);
data.push_back(aux);
names.push_back(n);
void
NamedData<DATA>::merge(NamedData<DATA2> &other)
{
- Assert(!is_constant, ExcConstantObject());
+ // see operator() below
+
+ // Assert(!is_constant, ExcConstantObject());
for (unsigned int i=0; i<other.size(); ++i)
{
void
NamedData<DATA>::merge(const NamedData<DATA2> &other)
{
- Assert(!is_constant, ExcConstantObject());
+ // see operator() below
+
+ // Assert(!is_constant, ExcConstantObject());
for (unsigned int i=0; i<other.size(); ++i)
{
names.push_back(other.name(i));
DATA &
NamedData<DATA>::operator() (unsigned int i)
{
- Assert(!is_constant, ExcConstantObject());
+ // Remove this assertion, since is_const() does not really assert
+ // what we would like to. This can only be fixed by using constness
+ // in AnyData.
+
+ // Assert(!is_constant, ExcConstantObject());
AssertIndexRange(i, size());
return data[i];
}
namespace Algorithms
{
+ OperatorBase::~OperatorBase()
+ {}
+
+ void OperatorBase::notify(const Event &e)
+ {
+ notifications += e;
+ }
+
+ void
+ OperatorBase::clear_events ()
+ {
+ notifications.clear();
+ }
+
+
#include "operator.inst"
}
--- /dev/null
+CMAKE_MINIMUM_REQUIRED(VERSION 2.8.9)
+INCLUDE(${DEAL_II_SOURCE_DIR}/cmake/setup_testsubproject.cmake)
+PROJECT(testsuite CXX)
+INCLUDE(${DEAL_II_TARGET_CONFIG})
+DEAL_II_PICKUP_TESTS()
--- /dev/null
+// ---------------------------------------------------------------------
+// $Id: named_data.cc 31349 2013-10-20 19:07:06Z maier $
+//
+// Copyright (C) 2000 - 2013 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+
+
+#include "../tests.h"
+#include <deal.II/base/any_data.h>
+#include <deal.II/base/logstream.h>
+
+#include <fstream>
+
+double d1 = 17.;
+
+void fill(AnyData& data)
+{
+ int i = 7;
+ data.add(i, " i 7");
+ data.add<double>(d1, " d 17.");
+ data.add<double*>(&d1, " d* 17.");
+ data.add<const double*>(&d1, "cd* 17.");
+ d1 = 18.;
+}
+
+
+void extract(const AnyData& data)
+{
+ // This set of tests with old functionality. Remove when deprecating
+ // index access
+ for (unsigned int i=0;i<data.size();++i)
+ deallog << i << '\t' << data.name(i) << std::endl;
+ deallog << data.name(0)
+ << '\t' << data.entry<int>(0) << std::endl;
+ deallog << data.name(1)
+ << '\t' << data.entry<double>(1) << std::endl;
+ deallog << data.name(2)
+ << '\t' << *data.entry<double*>(2) << std::endl;
+ deallog << data.name(3)
+ << '\t' << *data.entry<const double*>(3) << std::endl;
+
+ // From here on keep after deprecation
+ int i1 = data.entry<int>(" i 7");
+ const int i2 = data.entry<const int>(" i 7");
+ double d = data.entry<double>(" d 17.");
+ double* p2 = data.entry<double*>(" d* 17.");
+ const double * p3 = data.entry<const double*>("cd* 17.");
+
+ deallog << i1 << std::endl
+ << i2 << std::endl
+ << d << std::endl
+ << *p2 << std::endl
+ << *p3 << std::endl;
+ try
+ {
+ double* p3a = data.entry<double*>("cd* 17.");
+ deallog << p3a;
+ }
+ // catch(ExceptionBase e)
+ // {
+ // deallog << e.what() << std::endl;
+ // }
+ catch(...)
+ {
+ deallog << "Exception duly thrown" << std::endl;
+ }
+
+}
+
+int main()
+{
+ deal_II_exceptions::disable_abort_on_exception();
+
+ std::ofstream logfile("output");
+ deallog.attach(logfile);
+ deallog.depth_console(0);
+
+ AnyData data;
+ fill(data);
+ extract(data);
+}
--- /dev/null
+
+DEAL::0 i 7
+DEAL::1 d 17.
+DEAL::2 d* 17.
+DEAL::3 cd* 17.
+DEAL:: i 7 7
+DEAL:: d 17. 17.0000
+DEAL:: d* 17. 18.0000
+DEAL::cd* 17. 18.0000
+DEAL::7
+DEAL::7
+DEAL::17.0000
+DEAL::18.0000
+DEAL::18.0000
+DEAL::Exception duly thrown
--- /dev/null
+// ---------------------------------------------------------------------
+// $Id: newton_01.cc 31349 2013-10-20 19:07:06Z maier $
+//
+// Copyright (C) 2008 - 2013 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+// Test the compatibility functionality for Operator::operator()
+
+#include "../tests.h"
+#include <deal.II/algorithms/operator.h>
+#include <deal.II/lac/vector.h>
+
+using namespace Algorithms;
+
+class Op1 : public Operator<Vector<double> >
+{
+ public:
+ virtual void operator() (AnyData& out, const AnyData& in)
+ {
+ Vector<double>* u = out.entry<Vector<double>*>("u");
+ deallog << "u " << (*u)(0) << std::endl;
+
+ const Vector<double>* v = in.entry<Vector<double>*>("v");
+ deallog << "v " << (*v)(0) << std::endl;
+ }
+};
+
+class Op2 : public Operator<Vector<double> >
+{
+ public:
+ virtual void operator() (NamedData<Vector<double>*>& out,
+ const NamedData<Vector<double>*>& in)
+ {
+ Vector<double>* u = out(out.find("u"));
+ deallog << "u " << (*u)(0) << std::endl;
+
+ const Vector<double> * const v = in(in.find("v"));
+ deallog << "v " << (*v)(0) << std::endl;
+ }
+};
+
+void test(Operator<Vector<double> >& op)
+{
+ Vector<double> u(1);
+ Vector<double>* p = &u;
+ Vector<double> v(1);
+ u(0) = 3.;
+ v(0) = 7.;
+
+ AnyData o1;
+ o1.add (p, "u");
+ AnyData i1;
+ i1.add (&v, "v");
+
+ NamedData<Vector<double>*> o2;
+ o2.add (p, "u");
+ NamedData<Vector<double>*> i2;
+ i2.add (&v, "v");
+
+ deallog << "Call with AnyData" << std::endl;
+ op(o1, i1);
+ deallog << "Call with NamedData" << std::endl;
+ op(o2, i2);
+}
+
+int main()
+{
+ // deal_II_exceptions::disable_abort_on_exception();
+
+ std::ofstream logfile("output");
+ deallog.attach(logfile);
+ deallog.depth_console(0);
+
+ deallog << "Operator with AnyData" << std::endl;
+ Op1 op1;
+ test(op1);
+
+ deallog << "Operator with NamedData" << std::endl;
+ Op2 op2;
+ test(op2);
+
+ // deallog << "Operator with NamedData" << std::endl;
+ // Op2 op2;
+ // test(op1);
+}
--- /dev/null
+
+DEAL::Operator with AnyData
+DEAL::Call with AnyData
+DEAL::u 3.00000
+DEAL::v 7.00000
+DEAL::Call with NamedData
+DEAL::u 3.00000
+DEAL::v 7.00000
+DEAL::Operator with NamedData
+DEAL::Call with AnyData
+DEAL::u 3.00000
+DEAL::v 7.00000
+DEAL::Call with NamedData
+DEAL::u 3.00000
+DEAL::v 7.00000