]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added tensor and point types.
authorLuca Heltai <luca.heltai@sissa.it>
Wed, 19 Jul 2017 18:56:00 +0000 (20:56 +0200)
committerLuca Heltai <luca.heltai@sissa.it>
Sun, 23 Jul 2017 14:57:14 +0000 (16:57 +0200)
Moved some functions to internal namespaces.

include/deal.II/base/patterns_tools.h
source/base/parameter_handler.cc

index 23b51d1a49cb1f1ccb34868760afbb402e604f49..b15f5569d91ed380f050a21238bc4ebfcfb1d2f0 100644 (file)
@@ -29,6 +29,7 @@
 
 #include <map>
 #include <vector>
+#include <array>
 #include <string>
 #include <memory>
 #include <sstream>
 
 DEAL_II_NAMESPACE_OPEN
 
+using namespace Patterns;
 
 /**
- * Namespace for a few class and functions that act on values and patterns.
+ * Namespace for a few classes and functions that act on values and patterns, and
+ * allow to convert from non elementary types to strings and viceversa.
+ *
+ * A typical usage of these tools is in the following example:
+ *
+ * @code
+ * typedef std::vector<unsigned int> T;
+ *
+ * T vec(3);
+ * vec[0] = 1;
+ * vec[1] = 3;
+ * vec[2] = 5;
+ *
+ * auto pattern = PatternsTools::Convert<T>::to_pattern();
+ *
+ * std::cout << pattern->description() << std::endl;
+ * // [List of <[Integer]> of length 0...4294967295 (inclusive)]
+ *
+ * auto s = PatternsTools::Convert<T>::to_string(vec);
+ * std::cout << s << std::endl;
+ * // 1, 2, 3
+ *
+ * auto vec = PatternsTools::Convert<T>::to_value("2,3,4,5");
+ * // now vec has size 4, and contains the elements 2,3,4,5
+ *
+ * std::cout << internal::RankInfo<T>::list_rank << std::endl; // Outputs 1
+ * std::cout << internal::RankInfo<T>::map_rank  << std::endl; // Outputs 0
+ * @endcode
+ *
+ * Convert<T> is used by the function PatternsTools::add_parameter() in this
+ * namespace. Internally it uses the internal::RankInfo<T> class to decide how
+ * many different separators are required to convert the given type to a string.
+ *
+ * For example, to write vectors of vectors, the default is to use "," for the first
+ * (inner) separator, and ";" for the second (outer) separator, i.e.
+ *
+ * @code
+ * std::vector<std::vector<unsigned int>> vec;
+ * vec = Convert<decltype(vec)>::to_value("1,2,3 ; 4,5,6");
+ *
+ * s = convert<decltype(vec[0])>::to_string(vec[0]);
+ * // s now contains the string "1,2,3"
+ * @endcode
+ *
+ * Separators for List and Map compatible types are selected according to the
+ * rank of the list and map objects, using the arrays
+ * PatternsTools::internal::default_list_separator and
+ * PatternsTools::internal::default_map_separator.
+ *
+ * They are currently set to:
+ *
+ * @code
+ * default_list_separator{{" ",  ","  ,  ";"  ,  "|"  ,   "%"}};
+ * default_map_separator {{" ",  ":"  ,  "="  ,  "@"  ,   "#"}};
+ * @endcode
+ *
+ * When one needs a mixture of Lists and Map types, their RankInfo is computed by taking
+ * the maximum of the vector_rank of the Key and of the Value type, so that, for example,
+ * it is possible to have the following
+ * @code
+ * ... // Build compare class
+ * std::map<std::vector<unsigned int>, std::vector<double>, compare> map;
+ *
+ * map = convert<decltype(map)>::to_value("1,2,3 : 5.0,6.0,7.0  ; 8,9,10 : 11.0,12.0,13.0");
+ *
+ * @endcode
+ *
+ * Some non elementary types are supported, like Point<dim>(), or std::complex<double>.
+ * If you wish to support more types, you have to specialize the Convert struct as well
+ * as the RankInfo struct.
  *
  * @ingroup input
+ * @author Luca Heltai, 2017
  */
 namespace PatternsTools
 {
-  using namespace Patterns;
-
-
-  /**
-   * Store information about the rank types of the given class.
-   *
-   * A class has Rank equal to the number of different separators
-   * that are required to specify its element(s) in a string.
-   *
-   * This class is used to detect wether the class T is compatible
-   * with a Patterns::List pattern or with a Patterns::Map pattern, and
-   * to choose among a default list of separators for complex types,
-   * like vectors of vectors.
-   *
-   * Objects like Point<dim> or std::complex<double> are vector-likes, and
-   * have vector_rank 1. Elementary types, like `int`, `unsigned int`, `double`, etc.
-   * have vector_rank 0. `std::vector`, `std::list` and in general containers have rank
-   * equal to 1 + vector_rank of the contained type. Similarly for map types.
-   *
-   * A class with vector_rank_type::value = 0 is either elementary or a
-   * map. A class with map_rank_type::value = 0 is either a List compatible
-   * class, or an elementary type.
-   *
-   * Elementary types are not compatible with Patterns::List, but
-   * non elementary types, like Point<dim>(), or std::complex<double>, are
-   * compatible with the List type. Adding more compatible types is a matter
-   * of adding a specialization of this struct for the given type.
-   */
-  template<class T, class Enable=void>
-  struct RankInfo
-  {
-    typedef std::integral_constant<int, 0>::type vector_rank_type;
-    typedef std::integral_constant<int, 0>::type map_rank_type;
-  };
-
-
-  /**
-   * Return the default list separator for an object with the given vector rank.
-   *
-   * Objects like Point<dim> or std::complex<double> are vector-likes, and
-   * have rank 1. Elementary types, like `int`, `unsigned int`, `double`, etc.
-   * have rank 0. `std::vector`, `std::list` and in general containers have rank
-   * equal to 1 + rank of the contained type.
-   *
-   * This function helps in constructing patterns for non elementary types,
-   * like for example std::vector<std::vector<std::complex<double>>>.
-   */
-  std::string default_list_separator(unsigned int);
-
-
-
-  /**
-   * Return the default map separator for an object with the given map rank.
-   *
-   * This function helps in constructing patterns for non elementary types,
-   * like for example std::map<unsigned int, std::map<unsigned int, double>>
-   */
-  std::string default_map_separator(unsigned int);
-
-
   /**
    * Converter class. This class is used to generate strings and Patterns associated to
    * the given type, and to convert from a string to the given type and viceversa.
@@ -136,13 +148,16 @@ namespace PatternsTools
     /**
      * Return a string containing a textual version of the variable s. Use the pattern
      * passed to perform the conversion, or create and use a default one.
+     *
+     * The created string is
      */
     static std::string to_string(const T &s,
                                  const std::unique_ptr<PatternBase> &p = Convert<T>::to_pattern()) = delete;
 
 
     /**
-     * Convert a string to a value, using the given pattern, or a default one.
+     * Convert a string to a value, using the given pattern. Use the pattern
+     * passed to perform the conversion, or create and use a default one.
      */
     static T to_value(const std::string &s,
                       const std::unique_ptr<PatternBase> &p = Convert<T>::to_pattern()) = delete;
@@ -177,7 +192,40 @@ namespace PatternsTools
                   << "The string " << arg1 << " does not match the pattern \""
                   << arg2.description() << "\"");
   //@}
-
+  namespace internal
+  {
+    /**
+     * Store information about the rank types of the given class.
+     *
+     * A class has Rank equal to the number of different separators
+     * that are required to uniquely identify its element(s) in a string.
+     *
+     * This class is used to detect wether the class T is compatible
+     * with a Patterns::List pattern or with a Patterns::Map pattern.
+     *
+     * Objects like Point<dim> or std::complex<double> are vector-likes, and
+     * have vector_rank 1. Elementary types, like `int`, `unsigned int`, `double`, etc.
+     * have vector_rank 0. `std::vector`, `std::list` and in general containers have rank
+     * equal to 1 + vector_rank of the contained type. Similarly for map types.
+     *
+     * A class with list_rank::value = 0 is either elementary or a
+     * map. A class with map_rank::value = 0 is either a List compatible
+     * class, or an elementary type.
+     *
+     * Elementary types are not compatible with Patterns::List, but
+     * non elementary types, like Point<dim>(), or std::complex<double>, are
+     * compatible with the List type. Adding more compatible types is a matter
+     * of adding a specialization of this struct for the given type.
+     *
+     * @author Luca Heltai, 2017
+     */
+    template<class T, class Enable=void>
+    struct RankInfo
+    {
+      static constexpr int list_rank = 0;
+      static constexpr int map_rank = 0;
+    };
+  }
 }
 
 // ---------------------- inline and template functions --------------------
@@ -193,7 +241,7 @@ namespace PatternsTools
       if (std::is_integral<T>::value)
         return std_cxx14::make_unique<Integer>(std::numeric_limits<T>::min(), std::numeric_limits<T>::max());
       else if (std::is_floating_point<T>::value)
-        return std_cxx14::make_unique<Double>(std::numeric_limits<T>::min(), std::numeric_limits<T>::max());
+        return std_cxx14::make_unique<Double>(-std::numeric_limits<T>::max(), std::numeric_limits<T>::max());
       else
         {
           AssertThrow(false, ExcNotImplemented());
@@ -203,10 +251,13 @@ namespace PatternsTools
 
     static std::string to_string(const T &value, const std::unique_ptr<PatternBase> &p = Convert<T>::to_pattern())
     {
-      std::stringstream str;
-      str << value;
-      AssertThrow(p->match(str.str()), ExcNoMatch(str.str(), *p));
-      return str.str();
+      std::string str;
+      if (std::is_same<T, unsigned char>())
+        str = std::to_string((int)value);
+      else
+        str = std::to_string(value);
+      AssertThrow(p->match(str), ExcNoMatch(str, *p));
+      return str;
     }
 
     static T to_value(const std::string &s,
@@ -214,17 +265,28 @@ namespace PatternsTools
     {
       AssertThrow(p->match(s), ExcNoMatch(s, *p));
       std::istringstream is(s);
-      T i;
-      is >> i;
+      T value;
+      if (std::is_same<T,unsigned char>::value)
+        {
+          int i;
+          is >> i;
+          value=i;
+        }
+      else
+        is >> value;
       AssertThrow(!is.fail(), ExcMessage("Failed to convert from \"" + s + "\" to the type \""
                                          +boost::core::demangle(typeid(T).name()) + "\""));
-      return i;
+      return value;
     }
   };
 
-  //specialize a type for all of the STL containers and maps
+
   namespace internal
   {
+    const std::array<std::string, 5> default_list_separator {{" ",  ","  ,  ";"  ,  "|"  ,   "%"}};
+    const std::array<std::string, 5> default_map_separator {{" ",  ":"  ,  "="  ,  "@"  ,   "#"}};
+
+    //specialize a type for all of the STL containers and maps
     template <typename T>       struct is_stl_container:std::false_type {};
     template <typename T, std::size_t N> struct is_stl_container<std::array    <T,N>>    :std::true_type {};
     template <typename... Args> struct is_stl_container<std::vector            <Args...>>:std::true_type {};
@@ -258,26 +320,39 @@ namespace PatternsTools
     static constexpr bool const value = internal::is_stl_map<std::decay_t<T>>::value;
   };
 
-
-  // Rank of vector types
-  template<class T>
-  struct RankInfo<T, typename std::enable_if<is_stl_container<T>::value>::type >
+  namespace internal
   {
-    typedef typename std::integral_constant<int, RankInfo<typename T::value_type>::vector_rank_type::value+1>::type vector_rank_type;
-    typedef typename std::integral_constant<int, 0>::type map_rank_type;
-  };
+    // Rank of vector types
+    template<class T>
+    struct RankInfo<T, typename std::enable_if<is_stl_container<T>::value>::type >
+    {
+      static constexpr int list_rank = RankInfo<typename T::value_type>::list_rank + 1;
+      static constexpr int map_rank = RankInfo<typename T::value_type>::map_rank;
+    };
 
 
-  // Rank of vector types
-  template<class T>
-  struct RankInfo<T, typename std::enable_if<is_stl_map<T>::value>::type >
-  {
-    typedef typename std::integral_constant<int, std::max(RankInfo<typename T::key_type>::vector_rank_type::value,
-                                                          RankInfo<typename T::value_type>::vector_rank_type::value)+1>::type vector_rank_type;
-    typedef typename std::integral_constant<int, std::max(RankInfo<typename T::key_type>::map_rank_type::value,
-                                                          RankInfo<typename T::value_type>::map_rank_type::value)+1>::type map_rank_type;
-  };
+    // Rank of map types
+    template<class T>
+    struct RankInfo<T, typename std::enable_if<is_stl_map<T>::value>::type >
+    {
+      static constexpr int list_rank = std::max(RankInfo<typename T::key_type>::list_rank,
+                                                RankInfo<typename T::mapped_type>::list_rank) + 1;
+      static constexpr int map_rank = std::max(RankInfo<typename T::key_type>::map_rank,
+                                               RankInfo<typename T::mapped_type>::map_rank) + 1;
+    };
+
+    // Rank of Tensor types
+    template <int rank, int dim, class Number>
+    struct RankInfo<Tensor<rank, dim, Number>>
+    {
+      static constexpr int list_rank = rank + RankInfo<Number>::list_rank;
+      static constexpr int map_rank = RankInfo<Number>::map_rank;
+    };
 
+    template <int dim, class Number>
+    struct RankInfo<Point<dim,Number>>:RankInfo<Tensor<1,dim,Number>> {};
+
+  }
 
 
   template<class T>
@@ -286,8 +361,9 @@ namespace PatternsTools
     static std::unique_ptr<PatternBase> to_pattern()
     {
       return std_cxx14::make_unique<List>(*Convert<typename T::value_type>::to_pattern(),
-                                          0, std::numeric_limits<unsigned int>::max(),
-                                          default_list_separator(RankInfo<T>::vector_rank_type::value));
+                                          std::numeric_limits<typename T::size_type>::min(),
+                                          std::numeric_limits<typename T::size_type>::max(),
+                                          internal::default_list_separator[internal::RankInfo<T>::list_rank]);
     }
 
     static std::string to_string(const T &t,
@@ -309,8 +385,7 @@ namespace PatternsTools
       for (unsigned int i=1; i<vec.size(); ++i)
         s += p->get_separator() + " " + vec[i];
 
-      Assert(p->match(s), ExcMessage("No match for " + s +
-                                     " with pattern " + p->description()));
+      AssertThrow(pattern->match(s), ExcNoMatch(s, *p));
       return s;
     }
 
@@ -321,8 +396,7 @@ namespace PatternsTools
                       const std::unique_ptr<PatternBase> &pattern = Convert<T>::to_pattern())
     {
 
-      AssertThrow(pattern->match(s), ExcMessage("No match for " + s +
-                                                " using pattern " + pattern->description()));
+      AssertThrow(pattern->match(s), ExcNoMatch(s, *pattern));
 
       auto p = dynamic_cast<const Patterns::List *>(pattern.get());
       AssertThrow(p, ExcMessage("I need a List pattern to convert a string to a List type."));
@@ -347,8 +421,8 @@ namespace PatternsTools
       return std_cxx14::make_unique<Map>(*Convert<typename T::key_type>::to_pattern(),
                                          *Convert<typename T::mapped_type>::to_pattern(),
                                          0, std::numeric_limits<unsigned int>::max(),
-                                         default_list_separator(RankInfo<T>::vector_rank_type::value),
-                                         default_map_separator(RankInfo<T>::map_rank_type::value)
+                                         internal::default_list_separator[internal::RankInfo<T>::list_rank],
+                                         internal::default_map_separator[internal::RankInfo<T>::map_rank]
                                         );
     }
 
@@ -375,20 +449,15 @@ namespace PatternsTools
       for (unsigned int i=1; i<vec.size(); ++i)
         s += p->get_separator() + " " + vec[i];
 
-      Assert(p->match(s), ExcMessage("No match for " + s +
-                                     " with pattern " + p->description()));
+      AssertThrow(p->match(s), ExcNoMatch(s, *p));
       return s;
     }
 
-    /**
-     * Convert a string to a value, using the given pattern, or a default one.
-     */
     static T  to_value(const std::string &s,
                        const std::unique_ptr<PatternBase> &pattern = Convert<T>::to_pattern())
     {
 
-      AssertThrow(pattern->match(s), ExcMessage("No match for " + s +
-                                                " using pattern " + pattern->description()));
+      AssertThrow(pattern->match(s), ExcNoMatch(s, *pattern));
 
       auto p = dynamic_cast<const Patterns::Map *>(pattern.get());
       AssertThrow(p, ExcMessage("I need a Map pattern to convert a string to a List type."));
@@ -410,8 +479,111 @@ namespace PatternsTools
     }
   };
 
+
+  template<int rank, int dim, class Number>
+  struct Convert<Tensor<rank,dim,Number>>
+  {
+    typedef Tensor<rank,dim,Number> T;
+    static std::unique_ptr<PatternBase> to_pattern()
+    {
+      return std_cxx14::make_unique<List>(*Convert<typename T::value_type>::to_pattern(),
+                                          dim, dim,
+                                          internal::default_list_separator[internal::RankInfo<T>::list_rank]);
+    }
+
+    static std::string to_string(const T &t,
+                                 const std::unique_ptr<PatternBase> &pattern = Convert<T>::to_pattern())
+    {
+
+      auto p = dynamic_cast<const Patterns::List *>(pattern.get());
+      AssertThrow(p, ExcMessage("I need a List pattern to convert a string to a List type."));
+      auto base_p = p->get_base_pattern().clone();
+      std::vector<std::string> vec(dim);
+
+      for (unsigned int i=0; i<dim; ++i)
+        vec[i] = Convert<typename T::value_type>::to_string(t[i], base_p);
+
+      std::string s;
+      if (vec.size() > 0)
+        s = vec[0];
+      for (unsigned int i=1; i<vec.size(); ++i)
+        s += p->get_separator() + " " + vec[i];
+
+      AssertThrow(p->match(s), ExcNoMatch(s, *p));
+      return s;
+    }
+
+    static T to_value(const std::string &s,
+                      const std::unique_ptr<PatternBase> &pattern = Convert<T>::to_pattern())
+    {
+
+      AssertThrow(pattern->match(s), ExcNoMatch(s, *pattern));
+
+      auto p = dynamic_cast<const Patterns::List *>(pattern.get());
+      AssertThrow(p, ExcMessage("I need a List pattern to convert a string to a List type."));
+
+      auto base_p = p->get_base_pattern().clone();
+      T t;
+
+      auto v = Utilities::split_string_list(s,p->get_separator());
+      unsigned int i=0;
+      for (auto str : v)
+        t[i++] =  Convert<typename T::value_type>::to_value(str, base_p);
+
+      return t;
+    }
+  };
+
+
+  template<int dim, class Number>
+  struct Convert<Point<dim,Number>>
+  {
+
+    typedef Point<dim,Number> T;
+
+    static std::unique_ptr<PatternBase> to_pattern()
+    {
+      return Convert<Tensor<1,dim,Number>>::to_pattern();
+    }
+
+    static std::string to_string(const T &t,
+                                 const std::unique_ptr<PatternBase> &pattern = Convert<T>::to_pattern())
+    {
+      return Convert<Tensor<1,dim,Number>>::to_string(Tensor<1,dim,Number>(t),pattern);
+    }
+
+    static T to_value(const std::string &s,
+                      const std::unique_ptr<PatternBase> &pattern = Convert<T>::to_pattern())
+    {
+      return T(Convert<Tensor<1,dim,Number>>::to_value(s,pattern));
+    }
+  };
+
+
+  template <class ParameterType>
+  void add_parameter(const std::string           &entry,
+                     ParameterType                &parameter,
+                     ParameterHandler            &prm,
+                     const std::string           &documentation,
+                     const Patterns::PatternBase &pattern)
+  {
+
+    static_assert(std::is_const<ParameterType>::value == false,
+                  "You tried to add a parameter using a type "
+                  "that is const. Use a non-const type.");
+
+    prm.declare_entry(entry, PatternsTools::Convert<ParameterType,void>::to_string(parameter, pattern.clone()),
+                      pattern, documentation);
+
+    auto action = [&] (const std::string &val)
+    {
+      parameter = PatternsTools::Convert<ParameterType,void>::to_value(val, pattern.clone());
+    };
+    prm.add_action(entry, action);
+  }
 }
 
+
 DEAL_II_NAMESPACE_CLOSE
 
 #endif
index 5aab9be8f768dff9a4b61ae9fa15a8a80fb4713e..c8626df734bd8adec802f6401cad587c6cb8154f 100644 (file)
@@ -305,7 +305,8 @@ namespace Patterns
     std::istringstream str(test_string);
 
     double d;
-    if (!(str >> d))
+    str >> d;
+    if (str.fail())
       return false;
 
     if (!has_only_whitespace (str))

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.