]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add expand_instantiations in a separate folder
authormaier <maier@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 8 Aug 2012 19:35:52 +0000 (19:35 +0000)
committermaier <maier@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 8 Aug 2012 19:35:52 +0000 (19:35 +0000)
git-svn-id: https://svn.dealii.org/branches/branch_cmake@25786 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/common/expand_instantiations/CMakeLists.txt [new file with mode: 0644]
deal.II/common/expand_instantiations/expand_instantiations.cc [new file with mode: 0644]
deal.II/common/expand_instantiations/template-arguments [new file with mode: 0644]
deal.II/common/expand_instantiations/template-arguments.in [new file with mode: 0644]

diff --git a/deal.II/common/expand_instantiations/CMakeLists.txt b/deal.II/common/expand_instantiations/CMakeLists.txt
new file mode 100644 (file)
index 0000000..bbed267
--- /dev/null
@@ -0,0 +1,25 @@
+#TODO template-arguments.in
+
+ADD_EXECUTABLE(expand_instantiations expand_instantiations.cc)
+
+GET_TARGET_PROPERTY (
+  expand_instantiations_exe
+  expand_instantiations
+  LOCATION
+  )
+
+MACRO(macro_expand_instantiations inst_in_files)
+  FOREACH (inst_in_file ${inst_in_files})
+    STRING(REGEX REPLACE "\\.in$" "" inst_file "${inst_in_file}" )
+    ADD_CUSTOM_COMMAND (
+      OUTPUT ${CMAKE_CURRENT_BINARY_DIR}/${inst_file}
+      DEPENDS
+        expand_instantiations
+      COMMAND ${expand_instantiations_exe}
+      ARGS
+        ~/template-arguments
+        < ${CMAKE_CURRENT_SOURCE_DIR}/${inst_in_file}
+        > ${CMAKE_CURRENT_BINARY_DIR}/${inst_file}
+      )
+  ENDFOREACH()
+ENDMACRO()
diff --git a/deal.II/common/expand_instantiations/expand_instantiations.cc b/deal.II/common/expand_instantiations/expand_instantiations.cc
new file mode 100644 (file)
index 0000000..8f62d4e
--- /dev/null
@@ -0,0 +1,505 @@
+//----------------------------  expand_instantiations.cc  -------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2007, 2008, 2009, 2010, 2011 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//----------------------------  expand_instantiations.cc  -------------------
+
+// This is the program that we use to generate explicit instantiations for a
+// variety of template arguments. It takes two kinds of input files. The first
+// is given as arguments on the command line and contains entries of the
+// following form:
+// --------------------
+// REAL_SCALARS    := { double; float; long double }
+// COMPLEX_SCALARS := { std::complex<double>;
+//                      std::complex<float>;
+//                      std::complex<long double> }
+// VECTORS := { Vector<double>; Vector<float>; Vector<long double> }
+// --------------------
+//
+// The input file of this form is typically located in the common/ directory
+// and is generated by ./configure to contain the list of vectors etc that
+// make sense for the current configuration. For example, the list of VECTORS
+// is going to contain PETSc vectors if so configured.
+//
+// The second intput is read from the command line and consists of a sequence
+// of statements of the following form:
+// --------------------
+// for (u,v:VECTORS; z:SCALARS) { f(u, z, const v &); }
+// --------------------
+// Here, everything between {...} will be copied as many times as there are
+// combinations of arguments u,v in the list of substitutions given by
+// VECTORS. For each copy, the arguments u,v will be replaced by one of these
+// combinations.
+
+// Author: Wolfgang Bangerth, 2007
+
+
+#include <iostream>
+#include <fstream>
+#include <map>
+#include <list>
+#include <cstdlib>
+#include <string>
+
+// a map from the keys in the expansion lists to the list itself. For
+// instance, the example above will lead to the entry
+//      expansion_lists[REAL_SCALARS] = (double, float, long double)
+// in this map, among others
+std::map<std::string, std::list<std::string> >  expansion_lists;
+
+
+
+/* ======================== auxiliary functions ================= */
+
+
+// replace all occurrences of 'pattern' by 'substitute' in 'in', and
+// return the result
+std::string
+replace_all (const std::string &in,
+            const std::string &pattern,
+            const std::string &substitute)
+{
+  std::string x = in;
+  while (x.find(pattern) != std::string::npos)
+    x.replace (x.find(pattern),
+              pattern.size(),
+              substitute);
+  return x;
+}
+
+
+// extract from the start of #in the part of the string that ends with one of
+// the characters in #delim_list. The extracted part is deleted from #in
+// and returned. We skip characters in delim_list if they are preceded
+// by a backslash
+std::string
+get_substring_with_delim (std::string       &in,
+                         const std::string &delim_list)
+{
+  std::string x;
+  while (in.size() != 0)
+    {
+                                      // stop copying to the result
+                                      // if the current character is
+                                      // a delimiter, but only if the
+                                      // previous character was not a
+                                      // backslash
+      if ((delim_list.find (in[0]) != std::string::npos)
+         &&
+         !((x.size() > 0)
+           &&
+           (x[x.size()-1] == '\\')))
+       break;
+
+      x += in[0];
+      in.erase (0, 1);
+    }
+
+  return x;
+}
+
+
+// delete all whitespace at the beginning of the given argument
+void
+skip_space (std::string &in)
+{
+  while ((in.size() != 0)
+        &&
+        ((in[0] == ' ') || (in[0] == '\t') || (in[0] == '\n')))
+    in.erase (0, 1);
+}
+
+
+std::string remove_comments (std::string line)
+{
+  const std::string::size_type double_slash_comment = line.find ("//");
+  if (double_slash_comment != std::string::npos)
+    line.erase (double_slash_comment, std::string::npos);
+
+  const std::string::size_type slash_star_comment_begin = line.find ("/*");
+  if (slash_star_comment_begin != std::string::npos)
+    {
+      const std::string::size_type slash_star_comment_end = line.find ("*/");
+      if (slash_star_comment_end == std::string::npos)
+       {
+         std::cerr << "The program can currently only handle /* block */"
+                   << "comments that start and end within the same line."
+                   << std::endl;
+         std::exit (1);
+       }
+      line.erase (slash_star_comment_begin,
+                 slash_star_comment_end - slash_star_comment_begin + 2);
+    }
+
+  return line;
+}
+
+
+// read the whole file specified by the stream given as argument into a string
+// for simpler parsing, and return it
+std::string read_whole_file (std::istream &in)
+{
+  std::string whole_file;
+  while (in)
+    {
+      std::string line;
+      getline (in, line);
+
+      whole_file += remove_comments (line);
+      whole_file += '\n';
+    }
+                                  // substitute tabs by spaces, multiple
+                                  // spaces by single ones
+  for (unsigned int i=0; i<whole_file.size(); ++i)
+    if (whole_file[i] == '\t')
+      whole_file[i] = ' ';
+  while (whole_file.find("  ") != std::string::npos)
+    whole_file.replace (whole_file.find("  "), 2, " ");
+
+  return whole_file;
+}
+
+
+
+// split a given string assumed to consist of a list of substrings delimited
+// by a particular character into its components
+std::list<std::string>
+split_string_list (const std::string &s,
+                  const char         delimiter)
+{
+  std::string tmp = s;
+  std::list<std::string> split_list;
+
+                                  // split the input list
+  while (tmp.length() != 0)
+    {
+      std::string name;
+      name = tmp;
+
+      if (name.find(delimiter) != std::string::npos)
+       {
+         name.erase (name.find(delimiter), std::string::npos);
+         tmp.erase (0, tmp.find(delimiter)+1);
+       }
+      else
+       tmp = "";
+
+      skip_space (name);
+
+      while ((name.size() != 0) && (name[name.length()-1] == ' '))
+       name.erase (name.length()-1, 1);
+
+      split_list.push_back (name);
+    }
+
+  return split_list;
+}
+
+
+
+// return the given list but without empty entries
+std::list<std::string>
+delete_empty_entries (const std::list<std::string> &list)
+{
+  std::list<std::string> return_list;
+  for (std::list<std::string>::const_iterator i = list.begin();
+       i != list.end(); ++i)
+    if (*i != "")
+      return_list.push_back (*i);
+
+  return return_list;
+}
+
+
+
+// determine whether a given substring at position #pos and length #length in
+// the string #text is a real token, i.e. not just part of another word
+bool is_real_token (const std::string &text,
+                   const std::string::size_type pos,
+                   const std::string::size_type length)
+{
+  static const std::string token_chars ("abcdefghijklmnopqrstuvwxyz"
+                                       "ABCDEFGHIJKLMNOPQRSTUVWXYZ"
+                                       "0123456789"
+                                       "_");
+  if ((pos != 0) && (token_chars.find (text[pos-1]) != std::string::npos))
+    return false;
+
+  if ((pos+length < text.size()) &&
+      (token_chars.find (text[pos+length]) != std::string::npos))
+    return false;
+
+  return true;
+}
+
+
+// substitute all occurrences of #token in #text by #substitute. because a
+// replacement token could be a templated class like std::complex<double> and
+// because the token to the substituted may be a template argument itself, we
+// surround the substitution by a space which shouldn't matter in C++
+std::string substitute_tokens (const std::string &text,
+                              const std::string &token,
+                              const std::string &substitute)
+{
+  std::string x_text = text;
+  std::string::size_type pos = 0;
+  while ((pos = x_text.find(token, pos)) != std::string::npos)
+    {
+      if (is_real_token (x_text, pos, token.size()))
+       {
+         x_text.replace (pos, token.size(),
+                         std::string(" ")+substitute+std::string(" "));
+         pos += substitute.size()+2;
+       }
+      else
+       ++pos;
+    }
+
+  return x_text;
+}
+
+
+
+/* ======================== the main functions ================= */
+
+
+// read and parse the expansion lists like
+//     REAL_SCALARS    := { double; float; long double }
+// as specified at the top of the file and store them in
+// the global expansion_lists variable
+void read_expansion_lists (const std::string &filename)
+{
+  std::ifstream in (filename.c_str());
+
+  if (! in)
+    {
+      std::cerr << "Instantiation list file can not be read!"
+               << std::endl;
+      std::exit (1);
+    }
+
+                                  // read the entire file into a string for
+                                  // simpler processing. replace end-of-line
+                                  // characters by spaces
+  std::string whole_file = read_whole_file (in);
+
+                                  // now process entries of the form
+                                  // NAME := { class1; class2; ...}.
+  while (whole_file.size() != 0)
+    {
+      const std::string
+       name = get_substring_with_delim (whole_file, " :");
+
+      skip_space (whole_file);
+      if (whole_file.find (":=") != 0)
+       {
+         std::cerr << "Invalid entry <" << name << '>' << std::endl;
+         std::exit (1);
+       }
+      whole_file.erase (0, 2);
+      skip_space (whole_file);
+      if (whole_file.find ("{") != 0)
+       {
+         std::cerr << "Invalid entry <" << name << '>' << std::endl;
+         std::exit (1);
+       }
+      whole_file.erase (0, 1);
+      skip_space (whole_file);
+
+      std::string
+       expansion = get_substring_with_delim (whole_file, "}");
+
+      if (whole_file.find ("}") != 0)
+       {
+         std::cerr << "Invalid entry <" << name << '>' << std::endl;
+         std::exit (1);
+       }
+      whole_file.erase (0, 1);
+      skip_space (whole_file);
+
+                                      // assign but remove empty entries;
+                                      // this may happen if an expansion list
+                                      // ends in a semicolon (then we get an
+                                      // empty entry at the end), or if there
+                                      // are multiple semicolons after each
+                                      // other (this may happen if, for
+                                      // example, we have "Vector<double>;
+                                      // TRILINOS_VECTOR;" and if
+                                      // TRILINOS_VECTOR is an empty
+                                      // expansion after running ./configure)
+      expansion_lists[name]
+       = delete_empty_entries (split_string_list (expansion, ';'));
+    }
+}
+
+
+
+// produce all combinations of substitutions of the tokens given in the
+// #substitutions list in #text and output it to std::cout
+void substitute (const std::string &text,
+                const std::list<std::pair<std::string, std::string> > &substitutions)
+{
+                                  // do things recursively: if the list of
+                                  // substitutions has a single entry, then
+                                  // process all of them. otherwise, process
+                                  // the first in the list and call the
+                                  // function recursively with the rest of
+                                  // the substitutions
+  if (substitutions.size() > 1)
+    {
+                                      // do the first substitution, then call
+                                      // function recursively
+      const std::string name    = substitutions.front().first,
+                       pattern = substitutions.front().second;
+
+      const std::list<std::pair<std::string, std::string> >
+       rest_of_substitutions (++substitutions.begin(),
+                              substitutions.end());
+
+      for (std::list<std::string>::const_iterator
+            expansion = expansion_lists[pattern].begin();
+          expansion != expansion_lists[pattern].end();
+          ++expansion)
+       {
+         std::string new_text
+           = substitute_tokens (text, name, *expansion);
+
+         substitute (new_text, rest_of_substitutions);
+       }
+    }
+  else if (substitutions.size() == 1)
+    {
+                                      // do the substitutions
+      const std::string name    = substitutions.front().first,
+                       pattern = substitutions.front().second;
+
+      for (std::list<std::string>::const_iterator
+            expansion = expansion_lists[pattern].begin();
+          expansion != expansion_lists[pattern].end();
+          ++expansion)
+       std::cout << substitute_tokens (text, name, *expansion)
+                 << std::endl;
+    }
+  else
+    {
+      std::cout << text
+               << std::endl;
+    }
+}
+
+
+
+// process the list of instantiations given in the form
+// --------------------
+// for (u,v:VECTORS; z:SCALARS) { f(u, z, const v &); }
+// --------------------
+void process_instantiations ()
+{
+  std::string whole_file = read_whole_file (std::cin);
+
+                                  // process entries of the form
+                                  //  for (X:Y; A:B) { INST }
+  while (whole_file.size() != 0)
+    {
+      skip_space (whole_file);
+      if (whole_file.find ("for") != 0)
+       {
+         std::cerr << "Invalid instantiation list: missing 'for'" << std::endl;
+         std::exit (1);
+       }
+      whole_file.erase (0, 3);
+      skip_space (whole_file);
+      if (whole_file.find ("(") != 0)
+       {
+         std::cerr << "Invalid instantiation list: missing '('" << std::endl;
+         std::exit (1);
+       }
+      whole_file.erase (0, 1);
+      skip_space (whole_file);
+
+      const std::list<std::string>
+       substitutions_list
+       = split_string_list (get_substring_with_delim (whole_file,
+                                                      ")"),
+                            ';');
+      if (whole_file.find (")") != 0)
+       {
+         std::cerr << "Invalid instantiation list: missing ')'" << std::endl;
+         std::exit (1);
+       }
+      whole_file.erase (0, 1);
+      skip_space (whole_file);
+
+                                      // process the header
+      std::list<std::pair<std::string, std::string> >
+       substitutions;
+
+      for (std::list<std::string>::const_iterator
+            s = substitutions_list.begin();
+          s != substitutions_list.end(); ++s)
+       {
+         const std::list<std::string>
+           names_and_type = split_string_list (*s, ':');
+         if (names_and_type.size() != 2)
+           {
+             std::cerr << "Invalid instantiation header" << std::endl;
+             std::exit (1);
+           }
+
+         const std::list<std::string>
+           names = split_string_list (names_and_type.front(), ',');
+
+         for (std::list<std::string>::const_iterator
+                x = names.begin(); x != names.end(); ++x)
+           substitutions.push_back (std::make_pair (*x,
+                                                    names_and_type.back()));
+       }
+
+                                      // now read the part in {...}
+      skip_space (whole_file);
+      if (whole_file.find ("{") != 0)
+       {
+         std::cerr << "Invalid substitution text" << std::endl;
+         std::exit (1);
+       }
+      whole_file.erase (0, 1);
+      skip_space (whole_file);
+      const std::string text_to_substitute
+       = get_substring_with_delim (whole_file, "}");
+      whole_file.erase (0,1);
+      skip_space (whole_file);
+
+                                      // now produce the
+                                      // substitutions. first replace
+                                      // all occurrences of "\{" by
+                                      // "{"
+      substitute (replace_all(replace_all(text_to_substitute, "\\{", "{"),
+                             "\\}", "}"),
+                 substitutions);
+    }
+}
+
+
+
+int main (int argc, char **argv)
+{
+  if (argc < 2)
+    {
+      std::cerr << "Usage: " << std::endl
+               << "  expand_instantiations class_list_files < in_file > out_file"
+               << std::endl;
+      std::exit (1);
+    }
+
+  for (int i=1; i<argc; ++i)
+    read_expansion_lists (argv[i]);
+
+  process_instantiations ();
+}
diff --git a/deal.II/common/expand_instantiations/template-arguments b/deal.II/common/expand_instantiations/template-arguments
new file mode 100644 (file)
index 0000000..4940694
--- /dev/null
@@ -0,0 +1,75 @@
+REAL_SCALARS    := { double; float; long double }
+COMPLEX_SCALARS := { std::complex<double>;
+                     std::complex<float>;
+                     std::complex<long double> }
+
+DERIVATIVE_TENSORS := { double;
+                       Tensor<1,deal_II_dimension>;
+                       Tensor<2,deal_II_dimension> }
+
+DEAL_II_VEC_TEMPLATES := { Vector; BlockVector }
+
+SERIAL_VECTORS := { Vector<double>;
+                    Vector<float> ;
+                   Vector<long double>;
+
+                   BlockVector<double>;
+                    BlockVector<float>;
+                    BlockVector<long double>;
+
+                   parallel::distributed::Vector<double>;
+                   parallel::distributed::Vector<float> ;
+                   parallel::distributed::Vector<long double>;
+
+                    parallel::distributed::BlockVector<double>;
+                    parallel::distributed::BlockVector<float> ;
+                    parallel::distributed::BlockVector<long double>;
+
+                   ;
+                   ;
+                   ;
+                   ;
+
+                   ;
+                   ;
+                   ;
+                   ;
+                  }
+
+DOFHANDLERS := { DoFHandler<deal_II_dimension>;
+                 hp::DoFHandler<deal_II_dimension>;
+                MGDoFHandler<deal_II_dimension> }
+
+DOFHANDLER_TEMPLATES := { DoFHandler;
+                          hp::DoFHandler;
+                         MGDoFHandler }
+
+TRIANGULATION_AND_DOFHANDLER_TEMPLATES := { Triangulation;
+                                           DoFHandler;
+                                           hp::DoFHandler;
+                                           MGDoFHandler}
+
+TRIANGULATION_AND_DOFHANDLERS := { Triangulation<deal_II_dimension, deal_II_space_dimension>;
+                                   DoFHandler<deal_II_dimension, deal_II_space_dimension>;
+                                   hp::DoFHandler<deal_II_dimension, deal_II_space_dimension>;
+                                  MGDoFHandler<deal_II_dimension, deal_II_space_dimension> }
+
+
+FEVALUES_BASES := { FEValuesBase<deal_II_dimension>;
+                   FEFaceValuesBase<deal_II_dimension> }
+
+SPARSITY_PATTERNS := { SparsityPattern;
+                       CompressedSparsityPattern;
+                       CompressedSetSparsityPattern;
+                       CompressedSimpleSparsityPattern;
+                       ;
+
+                       BlockSparsityPattern;
+                       BlockCompressedSparsityPattern;
+                       BlockCompressedSetSparsityPattern;
+                       BlockCompressedSimpleSparsityPattern;
+                       ; }
+
+DIMENSIONS := { 1; 2; 3 }
+
+SPACE_DIMENSIONS := { 1; 2; 3 }
diff --git a/deal.II/common/expand_instantiations/template-arguments.in b/deal.II/common/expand_instantiations/template-arguments.in
new file mode 100644 (file)
index 0000000..09d4d33
--- /dev/null
@@ -0,0 +1,75 @@
+REAL_SCALARS    := { double; float; long double }
+COMPLEX_SCALARS := { std::complex<double>;
+                     std::complex<float>;
+                     std::complex<long double> }
+
+DERIVATIVE_TENSORS := { double;
+                       Tensor<1,deal_II_dimension>;
+                       Tensor<2,deal_II_dimension> }
+
+DEAL_II_VEC_TEMPLATES := { Vector; BlockVector }
+
+SERIAL_VECTORS := { Vector<double>;
+                    Vector<float> ;
+                   Vector<long double>;
+
+                   BlockVector<double>;
+                    BlockVector<float>;
+                    BlockVector<long double>;
+
+                   parallel::distributed::Vector<double>;
+                   parallel::distributed::Vector<float> ;
+                   parallel::distributed::Vector<long double>;
+
+                    parallel::distributed::BlockVector<double>;
+                    parallel::distributed::BlockVector<float> ;
+                    parallel::distributed::BlockVector<long double>;
+
+                   @DEAL_II_EXPAND_TRILINOS_VECTOR@;
+                   @DEAL_II_EXPAND_TRILINOS_MPI_VECTOR@;
+                   @DEAL_II_EXPAND_PETSC_VECTOR@;
+                   @DEAL_II_EXPAND_PETSC_MPI_VECTOR@;
+
+                   @DEAL_II_EXPAND_TRILINOS_BLOCKVECTOR@;
+                   @DEAL_II_EXPAND_TRILINOS_MPI_BLOCKVECTOR@;
+                   @DEAL_II_EXPAND_PETSC_BLOCKVECTOR@;
+                   @DEAL_II_EXPAND_PETSC_MPI_BLOCKVECTOR@;
+                  }
+
+DOFHANDLERS := { DoFHandler<deal_II_dimension>;
+                 hp::DoFHandler<deal_II_dimension>;
+                MGDoFHandler<deal_II_dimension> }
+
+DOFHANDLER_TEMPLATES := { DoFHandler;
+                          hp::DoFHandler;
+                         MGDoFHandler }
+
+TRIANGULATION_AND_DOFHANDLER_TEMPLATES := { Triangulation;
+                                           DoFHandler;
+                                           hp::DoFHandler;
+                                           MGDoFHandler}
+
+TRIANGULATION_AND_DOFHANDLERS := { Triangulation<deal_II_dimension, deal_II_space_dimension>;
+                                   DoFHandler<deal_II_dimension, deal_II_space_dimension>;
+                                   hp::DoFHandler<deal_II_dimension, deal_II_space_dimension>;
+                                  MGDoFHandler<deal_II_dimension, deal_II_space_dimension> }
+
+
+FEVALUES_BASES := { FEValuesBase<deal_II_dimension>;
+                   FEFaceValuesBase<deal_II_dimension> }
+
+SPARSITY_PATTERNS := { SparsityPattern;
+                       CompressedSparsityPattern;
+                       CompressedSetSparsityPattern;
+                       CompressedSimpleSparsityPattern;
+                       @DEAL_II_EXPAND_TRILINOS_SPARSITY_PATTERN@;
+
+                       BlockSparsityPattern;
+                       BlockCompressedSparsityPattern;
+                       BlockCompressedSetSparsityPattern;
+                       BlockCompressedSimpleSparsityPattern;
+                       @DEAL_II_EXPAND_TRILINOS_BLOCK_SPARSITY_PATTERN@; }
+
+DIMENSIONS := { 1; 2; 3 }
+
+SPACE_DIMENSIONS := { 1; 2; 3 }

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.