]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
introduce AnyData, OperatorBase and test that we do not break compatibility
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 10 Mar 2014 11:07:13 +0000 (11:07 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 10 Mar 2014 11:07:13 +0000 (11:07 +0000)
git-svn-id: https://svn.dealii.org/trunk@32633 0785d39b-7218-0410-832d-ea1e28bc413d

12 files changed:
deal.II/include/deal.II/algorithms/any_data.h [new file with mode: 0644]
deal.II/include/deal.II/algorithms/newton.h
deal.II/include/deal.II/algorithms/newton.templates.h
deal.II/include/deal.II/algorithms/operator.h
deal.II/include/deal.II/algorithms/operator.templates.h
deal.II/include/deal.II/base/named_data.h
deal.II/source/algorithms/operator.cc
tests/algorithms/CMakeLists.txt [new file with mode: 0644]
tests/algorithms/any_data_01.cc [new file with mode: 0644]
tests/algorithms/any_data_01.output [new file with mode: 0644]
tests/algorithms/operator_compatibility.cc [new file with mode: 0644]
tests/algorithms/operator_compatibility.output [new file with mode: 0644]

diff --git a/deal.II/include/deal.II/algorithms/any_data.h b/deal.II/include/deal.II/algorithms/any_data.h
new file mode 100644 (file)
index 0000000..5b90063
--- /dev/null
@@ -0,0 +1,248 @@
+// ---------------------------------------------------------------------
+// $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
index 23ac4519a45f6eff5932310dc1db2be84f7a4a13..97ac90b4db228f9dcfcd76715f5f0b4697b759b7 100644 (file)
@@ -21,6 +21,7 @@
 #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
 
@@ -46,7 +47,7 @@ namespace Algorithms
    * #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
@@ -105,7 +106,12 @@ namespace Algorithms
      * 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 &);
 
index 16c6b74483064e8d6ec7dc145632f496fbd01d63..acb61f4712c1c43a719eb231a2a568bd0381c7ba 100644 (file)
@@ -94,11 +94,18 @@ namespace Algorithms
   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;
@@ -108,19 +115,16 @@ namespace Algorithms
     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)
index de940de86f07cd1f3074209beae9fd065148f7f5..0bc6575e77592d640a7138df09631360668aa394 100644 (file)
@@ -19,6 +19,7 @@
 #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>
 
@@ -40,6 +41,8 @@ DEAL_II_NAMESPACE_OPEN
 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.
@@ -58,29 +61,29 @@ namespace Algorithms
    * 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
@@ -93,14 +96,73 @@ namespace Algorithms
     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;
   };
 
   /**
index 173d10b9b733de1a6ad25d2a818338018988d5c8..e2a7458569b3f32cd0e0cd5e36c5c57fa25aaa95 100644 (file)
@@ -22,23 +22,84 @@ DEAL_II_NAMESPACE_OPEN
 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;
   }
 
 
index d177e6a1a1d99b1e195c19c9a9d36420987e6f41..20b549080c3fc1930c79d9ec923178467c4af7ea 100644 (file)
@@ -28,6 +28,11 @@ DEAL_II_NAMESPACE_OPEN
 
 
 /**
+ * @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.
@@ -252,7 +257,9 @@ inline
 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);
 }
@@ -263,7 +270,9 @@ inline
 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);
@@ -277,7 +286,9 @@ inline
 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)
     {
@@ -294,7 +305,9 @@ inline
 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));
@@ -342,7 +355,11 @@ inline
 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];
 }
index 8f0bcc05169732f2e4200d1b69ac31375c96aada..75fbbaa24f1b4e8eecf4cf7d1317097fc9186b39 100644 (file)
@@ -37,6 +37,21 @@ DEAL_II_NAMESPACE_OPEN
 
 namespace Algorithms
 {
+  OperatorBase::~OperatorBase()
+  {}
+  
+  void OperatorBase::notify(const Event &e)
+  {
+    notifications += e;
+  }
+  
+  void
+  OperatorBase::clear_events ()
+  {
+    notifications.clear();
+  }
+  
+
 #include "operator.inst"
 }
 
diff --git a/tests/algorithms/CMakeLists.txt b/tests/algorithms/CMakeLists.txt
new file mode 100644 (file)
index 0000000..3b96045
--- /dev/null
@@ -0,0 +1,5 @@
+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()
diff --git a/tests/algorithms/any_data_01.cc b/tests/algorithms/any_data_01.cc
new file mode 100644 (file)
index 0000000..0422741
--- /dev/null
@@ -0,0 +1,92 @@
+// ---------------------------------------------------------------------
+// $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);
+}
diff --git a/tests/algorithms/any_data_01.output b/tests/algorithms/any_data_01.output
new file mode 100644 (file)
index 0000000..24a366e
--- /dev/null
@@ -0,0 +1,15 @@
+
+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
diff --git a/tests/algorithms/operator_compatibility.cc b/tests/algorithms/operator_compatibility.cc
new file mode 100644 (file)
index 0000000..422f936
--- /dev/null
@@ -0,0 +1,95 @@
+// ---------------------------------------------------------------------
+// $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);
+}
diff --git a/tests/algorithms/operator_compatibility.output b/tests/algorithms/operator_compatibility.output
new file mode 100644 (file)
index 0000000..97807d9
--- /dev/null
@@ -0,0 +1,15 @@
+
+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

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.