]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
first version of an implementation of a TableHandler and a ConvergenceTable class.
authorhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 9 Jul 1999 16:42:52 +0000 (16:42 +0000)
committerhartmann <hartmann@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 9 Jul 1999 16:42:52 +0000 (16:42 +0000)
git-svn-id: https://svn.dealii.org/trunk@1566 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/base/include/base/convergence_table.h [new file with mode: 0644]
deal.II/base/include/base/table_handler.h [new file with mode: 0644]
deal.II/base/source/convergence_table.cc [new file with mode: 0644]
deal.II/base/source/table_handler.cc [new file with mode: 0644]

diff --git a/deal.II/base/include/base/convergence_table.h b/deal.II/base/include/base/convergence_table.h
new file mode 100644 (file)
index 0000000..c01373b
--- /dev/null
@@ -0,0 +1,93 @@
+/*----------------------------   convergence_table.h     ---------------------------*/
+/*      $Id$                 */
+#ifndef __convergence_table_H
+#define __convergence_table_H
+/*----------------------------   convergence_table.h     ---------------------------*/
+
+
+#include <base/table_handler.h>
+
+
+/**
+ * The #ConvergenceTable# class is an application to the #TableHandler# class and
+ * stores some convergence data as residuals of the cg-method, or some evaluated
+ * $L^2$-errors of discrete solutions, etc,
+ * and evaluates convergence rates or orders.
+ * The already implemented #RateMode# is #ConvergenceTable::standard#, 
+ * where the convergence rate is the quotient of two following rows, and 
+ * #ConvergenceTable::standard_order#, that evaluates the order of convergence.
+ * These standard evaluations are useful for global refinement, for local refinement
+ * this may not be a appropriate method, as the convergence rates should be
+ * set in relation to the number of cells or the number of DoFs. The
+ * implementations of these not standard methods is left to the fantasy of a
+ * fanatic user.
+ *
+ * \subsection{Usage}
+ * The number of cells and the number of DoFs are added to the table by
+ * calling #add_run(unsigned int ncells, unsigned int ndofs#. The
+ * data is added by #add_value(...)# of the base class #TableHandler#.
+ * Before the output of the table the function #evaluate_convergence_rates#
+ * may be (also multiply) called. 
+ *  
+ *
+ * @author Ralf Hartmann, 1999
+ */
+class ConvergenceTable: public TableHandler
+{
+  public:
+                                    /**
+                                     * Constructor.
+                                     */
+    ConvergenceTable();
+    
+    enum RateMode
+    {
+         none, standard, standard_order, n_cells, n_dofs
+    };
+                                    /**
+                                     * Adds the basic information
+                                     * (number of cells, number of DoFs)
+                                     * of a new run.
+                                     */
+    void add_run(unsigned int ncells, unsigned int ndofs);
+
+                                    /**
+                                     * Evaluates the convergence rates of the
+                                     * column #key# due to a #RateMode# and
+                                     * merges the rate column and this column
+                                     * to a supercolumn.
+                                     */
+    void evaluate_convergence_rates(const string &key, const RateMode conv_rate);
+
+                                    /**
+                                     * Evaluate the convergence rates of
+                                     * all (not "n cells" or "n dofs") 
+                                     * columns due to a #RateMode#.
+                                     */
+    void evaluate_convergence_rates(const RateMode conv_rate);
+
+                                    /**
+                                     * Exception
+                                     */
+    DeclException0 (ExcWrongValueType)
+
+  private:
+                                    /**
+                                     * Preset string "n cells".
+                                     */
+    string n_cells_string;
+
+                                    /**
+                                     * Preset string "n dofs".
+                                     */
+    string n_dofs_string;
+};
+
+    
+
+
+
+/*----------------------------   convergence_table.h     ---------------------------*/
+/* end of #ifndef __convergence_table_H */
+#endif
+/*----------------------------   convergence_table.h     ---------------------------*/
diff --git a/deal.II/base/include/base/table_handler.h b/deal.II/base/include/base/table_handler.h
new file mode 100644 (file)
index 0000000..27ccb2e
--- /dev/null
@@ -0,0 +1,382 @@
+/*----------------------------   table_handler.h     ---------------------------*/
+/*      $Id$                 */
+#ifndef __table_handler_H
+#define __table_handler_H
+/*----------------------------   table_handler.h     ---------------------------*/
+
+
+
+#include <base/exceptions.h>
+
+#include <map>
+#include <vector>
+#include <string>
+#include <iostream>
+
+#include <fstream.h>
+
+
+
+/**
+ * Abstract base class for the #TableEntry# class. See there.
+ * This class is not to be used by the user.
+ *
+ * @author Ralf Hartmann, 1999
+ */
+class TableEntryBase 
+{
+  public:
+                                    /**
+                                     * Constructor.
+                                     */
+    TableEntryBase();
+                                    /**
+                                     * Virtual destructor.
+                                     */    
+    virtual ~TableEntryBase();
+    
+                                    /**
+                                     * Write the table entry as text.
+                                     */
+    virtual void write_text (ostream &) const =0;
+
+                                    /**
+                                     * Write the table entry in tex format.
+                                     */
+    virtual void write_tex (ostream &) const =0;
+};
+
+
+/**
+ * A #TableEntry# stores the value of an table entry.
+ * The value type of this table entry is arbitrary. For
+ * a #TableEntry<typename value_type># with not common value
+ * type overload the output function.
+ * This class is not to be used by the user.
+ *
+ * For more detail see the #TableHandler# class.
+ *
+ * @author Ralf Hartmann, 1999
+ */
+template <typename value_type>
+class TableEntry : public TableEntryBase 
+{
+  public:
+                                    /**
+                                     * Constructor.
+                                     */
+    TableEntry(const value_type value);
+
+                                    /**
+                                     * Destructor.
+                                     */
+    virtual ~TableEntry();
+
+                                    /**
+                                     * Returns the value of this
+                                     * table entry.
+                                     */
+    value_type value() const;
+
+                                    /**
+                                     * Write the table entry as text.
+                                     */
+    virtual void write_text (ostream &out) const;
+
+                                    /**
+                                     * Write the table entry in tex format.
+                                     */
+    virtual void write_tex (ostream &out) const;
+
+  private:
+                                    /**
+                                     * Stored value.
+                                     */
+    const value_type val;
+};
+
+/**
+ * The #TableHandler# stores #TableEntries# of arbitrary value type and 
+ * writes the table as text or in tex format to a tex file. The value type
+ * actually may vary from column to column as from row to row. 
+ *
+ * \subsection{Usage}
+ *
+ * The most important function is the template function 
+ * #add_value(const string &key, const value_type value)#, that adds a column
+ * with the name #key# to the table if this column does not yet exist and adds the 
+ * value of #value_type# (e.g. unsigned int, double, string, ...) to this column.
+ * After the table is complete there are different possibilities of output, e.g.
+ * into a tex file with #write_tex(ofstream &file)# or as text with 
+ * #write_text (ostream &out)#. Two (or more) columns may be merged into a
+ * 'supercolumn' by twice (or multiple) calling #add_column_to_supercolumn(...)#,
+ * see there. Additionally there is a function to set for each column 
+ * the precision of the output of numbers, and there are several functions to
+ * prescribe the format and the captions the columns are written with in tex mode.
+ *
+ * \subsection{Example}
+ * This is a simple example demonstrating the usage of this class. The first column
+ * includes the numbers i=1..n, the second 1^2...n^2, the third sqrt(1)...sqrt(n),
+ * where the second and third columns are merged into one supercolumn 
+ * with the superkey `squares and roots'. Additionally the first column is
+ * aligned to the right (the default was `centered') and the precision of
+ * the square roots are given to be 6 (instead of 4 as default).
+ *
+ * \begin{verbatim}
+ * TableHandler table();
+ *
+ * for (unsigned int i=1; i<=n; ++i)
+ *   {
+ *     table.add_value("numbers", i);
+ *     table.add_value("squares", i*i);
+ *     table.add_value("square roots", sqrt(i));
+ *   }
+ *                                  // merge the second and third column
+ * table.add_column_to_supercolumn("squares", "squares and roots");
+ * table.add_column_to_supercolumn("square roots", "squares and roots");
+ *
+ *                                  // additional settings
+ * table.set_tex_format("numbers", "r");
+ * table.set_precision("square roots", 6);
+ *
+ *                                  // output
+ * ofstream out_file("number_table.tex");
+ * table.write_tex(out_file);
+ * out_file.close();
+ * \end{verbatim}
+ *
+ * @author Ralf Hartmann, 1999
+ */
+class TableHandler
+{
+  public:
+
+                                    /**
+                                     * Constructor.
+                                     */
+    TableHandler();
+    
+                                    /**
+                                     * Adds a column (if not yet existent) with
+                                     * the key #key#
+                                     * and adds the value of #value_type#
+                                     * to the column.
+                                     */
+    template <typename value_type>
+    void add_value(const string &key, const value_type value);
+    
+                                    /**
+                                     * Creates a sypercolumn (if not yet
+                                     * existent) and includes column to it.
+                                     * The keys of the column and the supercolumn
+                                     * are key and superkey respectively.
+                                     * To merge two columns c1 and c2 to
+                                     * a supercolumn sc, hence call
+                                     * #add_column_to_supercolumn(c1,sc)# and
+                                     * #add_column_to_supercolumn(c2,sc)#.
+                                     *
+                                     * Concerning the order of the columns,
+                                     * the supercolumn replaces the first
+                                     * column that is added to the supercolumn.
+                                     * Within the supercolumn the order
+                                     * the order of output follows the order
+                                     * the columns are added to the supercolumn.
+                                     */
+    void add_column_to_supercolumn(const string &key, const string &superkey);
+
+                                    /**
+                                     * Change the order of columns and
+                                     * supercolumns in the table.
+                                     */
+    void set_column_order(const vector<string> &new_order);
+    
+                                    /**
+                                     * Sets the output format of a column.
+                                     */
+    void set_precision(const string &key, const unsigned int precision);
+
+                                    /**
+                                     * Sets the caption of the column #key#
+                                     * for tex output.
+                                     */
+    void set_tex_caption(const string &key, const string &tex_caption);
+
+                                    /**
+                                     * Sets the caption the the supercolumn
+                                     * #superkey# for tex output.
+                                     */
+    void set_tex_supercaption(const string &superkey, const string &tex_supercaption);
+
+                                    /**
+                                     * Sets the tex output format of a column.
+                                     * e.g. "c", "r", "l".
+                                     */
+    void set_tex_format(const string &key, const string &format);
+
+                                    /**
+                                     * Write table as formatted text, e.g.
+                                     * to the standard output.
+                                     */
+    void write_text (ostream &out) const;
+
+                                    /**
+                                     * Write table as a tex file.
+                                     */
+    void write_tex (ofstream &file) const;
+   
+                                    /**
+                                     * Exception
+                                     */
+    DeclException1 (ExcColumnNotExistent,
+                   string,
+                   << "Column <" << arg1 << "> does not exist.");
+
+                                    /**
+                                     * Exception
+                                     */
+    DeclException1 (ExcSuperColumnNotExistent,
+                   string,
+                   << "Supercolumn <" << arg1 << "> does not exist.");
+
+                                    /**
+                                     * Exception
+                                     */
+    DeclException1 (ExcColumnOrSuperColumnNotExistent,
+                   string,
+                   << "Column or supercolumn <" << arg1 << "> does not exist.");
+    
+                                    /**
+                                     * Exception
+                                     */
+    DeclException4 (ExcWrongNumberOfDataEntries,
+                   string, int, string, int,
+                   << "Column <" << arg1 << "> has got " << arg2
+                   << "rows, but Column <" << arg3 << "> has got " << arg4 << ".");
+
+                                    /**
+                                     * Exception
+                                     */
+    DeclException1 (ExcUndefinedTexFormat,
+                   string,
+                   << "<" << arg1 << "> is not a tex column format. Use l,c,r.");
+
+  protected:
+
+
+    struct Column
+    {
+                                        /**
+                                         * Constructor needed by stl_map.h
+                                         */
+       Column();
+
+                                        /**
+                                         * Constructor.
+                                         */
+       Column(const string &tex_caption);
+       
+                                        /**
+                                         * Destructor.
+                                         */
+       ~Column();
+       
+       vector<TableEntryBase *> entries;
+       
+                                        /**
+                                         * The caption of the column in tex output.
+                                         * By default, this is the key string that
+                                         * is given to the #TableHandler# by
+                                         * #TableHandler::add_value(...)#. This may
+                                         * be changed by calling
+                                         * #TableHandler::set_tex_caption(...)#.
+                                         */
+       string tex_caption;
+
+                                        /**
+                                         * The column format in tex output.
+                                         * By default, this is #"c"#, meaning
+                                         * `centered'. This may
+                                         * be changed by calling
+                                         * #TableHandler::set_tex_format(...)#
+                                         * with #"c", "r", "l"# for centered,
+                                         * right or left.
+                                         */
+       
+       string tex_format;
+
+                                        /**
+                                         * Double or float entries are written with
+                                         * this (by the user set) precision.
+                                         * The default is 4.
+                                         */
+       unsigned int precision;
+
+                                        /**
+                                         * Flag that may be used by derived classes.
+                                         */
+       unsigned int flag;
+    };
+
+                                    /**
+                                     * Help function that gives a vector of
+                                     * the keys of all columns that are mentioned
+                                     * in column_order, where each supercolumn key
+                                     * is replaced by its subcolumn keys.
+                                     *
+                                     * This function implicitly checks the
+                                     * consistency of the data.
+                                     */
+    void get_selected_columns(vector<string> &sel_columns) const;
+    
+                                    /**
+                                     * Builtin function, that checks if
+                                     * the number of rows is equal in every
+                                     * row. This function is e.g. called before
+                                     * writing output.
+                                     */
+    unsigned int check_n_rows() const;
+
+                                    /**
+                                     * Stores the column and 
+                                     * supercolumn keys in the user wanted
+                                     * order.
+                                     * By default this is the order of
+                                     * adding the columns. This order may be
+                                     * changed by #set_column_order(...)#.
+                                     */
+    vector<string> column_order;
+
+                                    /**
+                                     * Maps the column keys to the columns
+                                     * (not supercolumns).
+                                     */
+    map<string,Column> columns;
+
+                                    /**
+                                     * Maps each supercolumn key to the
+                                     * the keys of its subcolumns in the right order.
+                                     * It is allowed that a supercolumn has got
+                                     * the same key as a column.
+                                     */
+    map<string, vector<string> > supercolumns;
+
+                                    /**
+                                     * Maps the supercolumn keys to the
+                                     * captions of the supercolumns that
+                                     * are used in tex output.
+                                     *
+                                     * By default these are just the
+                                     * supercolumn keys but they may be changed
+                                     * by #set_tex_supercaptions(...)#.
+                                     */
+    map<string, string> tex_supercaptions;
+};
+
+
+
+
+
+/*----------------------------   table_handler.h     ---------------------------*/
+/* end of #ifndef __table_handler_H */
+#endif
+/*----------------------------   table_handler.h     ---------------------------*/
diff --git a/deal.II/base/source/convergence_table.cc b/deal.II/base/source/convergence_table.cc
new file mode 100644 (file)
index 0000000..fdcca35
--- /dev/null
@@ -0,0 +1,99 @@
+/*----------------------------   convergence_table.cc     ---------------------------*/
+/*      $Id$                 */
+
+
+#include <base/convergence_table.h>
+#include <math.h>
+
+
+ConvergenceTable::ConvergenceTable():
+               n_cells_string("n cells"),
+               n_dofs_string("n dofs")   {}
+
+
+void ConvergenceTable::add_run(unsigned int ncells, unsigned int ndofs)
+{
+  add_value(n_cells_string, ncells);
+  add_value(n_dofs_string, ndofs);
+}
+
+
+void ConvergenceTable::evaluate_convergence_rates(const string &key,
+                                                 const RateMode rate_mode)
+{
+  Assert(columns.count(key), ExcColumnNotExistent(key));
+  
+  vector<TableEntryBase *> &entries=columns[key].entries;
+  string rate_key=key;
+
+  const unsigned int n=entries.size();
+  
+  vector<double> values(n);
+  for (unsigned int i=0; i<n; ++i)
+    {
+      if (dynamic_cast<TableEntry<double>*>(entries[i]) != 0)
+       values[i]=dynamic_cast<TableEntry<double>*>(entries[i])->value();
+      else if (dynamic_cast<TableEntry<float>*>(entries[i]) != 0)
+       values[i]=dynamic_cast<TableEntry<float>*>(entries[i])->value();
+      else
+       Assert(false, ExcWrongValueType());
+    }
+  
+  switch (rate_mode)
+    {
+      case standard:
+           rate_key+="s";
+           add_value(rate_key, string("-"));
+           for (unsigned int i=1; i<n; ++i)
+             add_value(rate_key, values[i-1]/values[i]);
+           break;
+      case standard_order:
+           rate_key+="so";
+           add_value(rate_key, string("-"));
+           for (unsigned int i=1; i<n; ++i)
+             add_value(rate_key, log(values[i-1]/values[i])/log(2));
+           break;
+      case n_cells:
+      case n_dofs:
+      case none:
+           break;
+      default:
+           ExcNotImplemented();  
+    }
+
+  Assert(columns.count(rate_key), ExcInternalError());  
+  columns[rate_key].flag=1;
+  set_precision(rate_key, 3);
+
+  string superkey="s"+key;
+  if (!supercolumns.count(superkey))
+    {
+      add_column_to_supercolumn(key, superkey);
+      set_tex_supercaption(superkey, columns[key].tex_caption);
+    }
+
+  add_column_to_supercolumn(rate_key, superkey);
+}
+
+
+void ConvergenceTable::evaluate_convergence_rates(const RateMode rate_mode)
+{
+                                  // make sure that no convergence rates
+                                  // are evaluated for the n_cells and the
+                                  // n_dofs column.
+  Assert(columns.count(n_cells_string), ExcInternalError());  
+  Assert(columns.count(n_dofs_string), ExcInternalError());
+  
+  columns[n_cells_string].flag=1;
+  columns[n_dofs_string].flag=1;
+  
+  check_n_rows();
+
+  for (map<string, Column>::const_iterator col_iter=columns.begin();
+       col_iter!=columns.end(); ++col_iter)
+    if (!col_iter->second.flag)
+      evaluate_convergence_rates(col_iter->first, rate_mode);
+}
+
+
+/*----------------------------   convergence_table.cc     ---------------------------*/
diff --git a/deal.II/base/source/table_handler.cc b/deal.II/base/source/table_handler.cc
new file mode 100644 (file)
index 0000000..d379349
--- /dev/null
@@ -0,0 +1,410 @@
+/*----------------------------   table_handler.cc     ---------------------------*/
+/*      $Id$                 */
+
+
+
+#include <base/table_handler.h>
+
+#include <iomanip>
+
+
+
+
+TableEntryBase::TableEntryBase ()
+{}
+
+
+TableEntryBase::~TableEntryBase ()
+{}
+
+/*---------------------------------------------------------------------*/
+
+template <typename number>
+TableEntry<number>::TableEntry(const number value):
+               val(value)  {}
+
+template <typename number>
+TableEntry<number>::~TableEntry()
+{}
+
+template <typename number>
+number TableEntry<number>::value() const
+{
+  return val;
+}
+
+
+template <typename number>
+void TableEntry<number>::write_tex (ostream &out) const
+{
+  out << val;
+}
+
+template <typename number>
+void TableEntry<number>::write_text (ostream &out) const
+{
+  out << val;
+}
+
+
+/*---------------------------------------------------------------------*/
+
+TableHandler::Column::Column(const string &tex_caption):
+               tex_caption(tex_caption),
+               tex_format("c"),
+               precision(4),
+               flag(0)  {}
+
+
+TableHandler::Column::Column():
+               tex_caption(),
+               tex_format("c"),
+               precision(4),
+               flag(0)  {}
+
+
+TableHandler::Column::~Column()
+{
+  for (unsigned int i=0; i<entries.size(); ++i)
+    delete entries[i];
+}
+
+/*---------------------------------------------------------------------*/
+
+
+TableHandler::TableHandler()
+{}
+
+
+template <typename number>
+void TableHandler::add_value(const string &key, number value)
+{
+  if (!columns.count(key))
+    {
+      pair<string, Column> new_column(key, Column(key));
+      columns.insert(new_column);
+      column_order.push_back(key);
+    }
+
+  TableEntry<number> *entry_ptr=new TableEntry<number>(value);
+  const map<string, Column>::iterator col_iter=columns.find(key);
+  Assert(col_iter!=columns.end(), ExcInternalError());
+  
+  col_iter->second.entries.push_back(entry_ptr);
+}
+
+
+
+void TableHandler::add_column_to_supercolumn(const string &key, const string &superkey)
+{
+  Assert(columns.count(key), ExcColumnNotExistent(key));
+
+  if (!supercolumns.count(superkey))
+    {
+      pair<string, vector<string> > new_column(superkey, vector<string>());
+      supercolumns.insert(new_column);
+                                      // replace key in column_order
+                                      // by superkey
+      for (unsigned int j=0; j<column_order.size(); ++j)
+       if (column_order[j]==key)
+         {
+           column_order[j]=superkey;
+           break;
+         }
+    } 
+  else 
+    {
+                                      // remove key from column_order
+                                      // for erase we need an iterator
+      for (vector<string>::iterator order_iter=column_order.begin(); 
+          order_iter!=column_order.end(); ++order_iter)
+       if (*order_iter==key)
+         {
+           column_order.erase(order_iter);
+           break;
+         }
+    }
+
+  if (supercolumns.count(superkey))
+    {
+      supercolumns[superkey].push_back(key);
+                                      // By default set the
+                                      // tex_supercaption to superkey
+      pair<string, string> new_tex_supercaption(superkey, superkey);
+      tex_supercaptions.insert(new_tex_supercaption);
+    }
+  else
+    Assert(false, ExcInternalError());
+}
+
+
+
+void TableHandler::set_column_order(const vector<string> &new_order)
+{
+  for (unsigned int j=0; j<new_order.size(); ++j)
+    Assert(supercolumns.count(new_order[j]) || columns.count(new_order[j]),
+          ExcColumnOrSuperColumnNotExistent(new_order[j]));
+
+  column_order=new_order;
+}
+
+
+void TableHandler::set_tex_caption(const string &key, const string &tex_caption)
+{
+  Assert(columns.count(key), ExcColumnNotExistent(key));
+  columns[key].tex_caption=tex_caption;
+}
+
+
+void TableHandler::set_tex_supercaption(const string &superkey, const string &tex_supercaption)
+{
+  Assert(supercolumns.count(superkey), ExcSuperColumnNotExistent(superkey));
+  Assert(tex_supercaptions.count(superkey), ExcInternalError());
+  tex_supercaptions[superkey]=tex_supercaption;
+}
+
+
+void TableHandler::set_tex_format(const string &key, const string &tex_format)
+{
+  Assert(columns.count(key), ExcColumnNotExistent(key));
+  Assert(tex_format=="l" || tex_format=="c" || tex_format=="r",
+        ExcUndefinedTexFormat(tex_format));
+  columns[key].tex_format=tex_format;
+}
+
+void TableHandler::set_precision(const string &key, unsigned int precision)
+{
+  Assert(columns.count(key), ExcColumnNotExistent(key));
+  columns[key].precision=precision;
+}
+
+void TableHandler::write_text(ostream &out) const
+{
+  vector<string> sel_columns;
+  get_selected_columns(sel_columns);
+
+
+                                  // write the caption line
+  for (unsigned int j=0; j<column_order.size(); ++j)
+    {
+      string key=column_order[j];
+      const map<string, vector<string> >::const_iterator 
+       super_iter=supercolumns.find(key);
+
+      unsigned int max_string_size=6;
+      if (super_iter!=supercolumns.end())
+       max_string_size+=8*(super_iter->second.size()-1);
+       
+      if (key.size()>max_string_size)
+       key.erase(max_string_size);
+      
+      out.setf(ios::left);
+      out << setw(max_string_size)
+//       << setfill('*') 
+         << key.c_str() << "\t";
+    }
+  out << endl;
+
+  const unsigned int n_rows=check_n_rows();
+  for (unsigned int i=0; i<n_rows; ++i)
+    {
+      const unsigned int n_cols=sel_columns.size();
+      
+      for (unsigned int j=0; j<n_cols; ++j)
+       {
+         string key=sel_columns[j];
+                                          // avoid `column[key]'
+         const map<string, Column>::const_iterator
+           col_iter=columns.find(key);
+         Assert(col_iter!=columns.end(), ExcInternalError());
+         
+         const Column &column=col_iter->second;
+         
+         out << setprecision(column.precision);
+         column.entries[i]->write_tex(out);
+         
+         if (j<n_cols-1)
+           out << "\t";
+       }
+      out << endl;
+    }
+}
+
+
+
+void TableHandler::write_tex(ofstream &out) const
+{
+  out << "\\documentclass[10pt]{report}" << endl
+      << "\\usepackage{float}" << endl << endl << endl
+      << "\\begin{document}" << endl
+      << "\\begin{table}[H]" << endl
+      << "\\begin{center}" << endl
+      << "\\begin{tabular}{|";
+
+  vector<string> sel_columns;
+  get_selected_columns(sel_columns);
+
+                                  // write the column formats
+  for (unsigned int j=0; j<column_order.size(); ++j)
+    {
+      string key=column_order[j];
+                                      // avoid `supercolumns[key]'
+      const map<string, vector<string> >::const_iterator 
+       super_iter=supercolumns.find(key);
+
+      if (super_iter!=supercolumns.end())
+       {
+         const unsigned int n_subcolumns=super_iter->second.size();
+         for (unsigned int k=0; k<n_subcolumns; ++k)
+           {
+                                              // avoid `columns[supercolumns[key]]'
+             const map<string, Column>::const_iterator
+               col_iter=columns.find(super_iter->second[k]);
+             Assert(col_iter!=columns.end(), ExcInternalError());
+             
+             out << col_iter->second.tex_format << "|";
+           }
+       }
+      else
+       {
+                                          // avoid `columns[key]';
+         const map<string, Column>::const_iterator
+           col_iter=columns.find(key);
+         Assert(col_iter!=columns.end(), ExcInternalError());
+         out << col_iter->second.tex_format << "|";
+       }
+    }
+  out << "} \\hline" << endl;
+
+                                  // write the caption line of the table
+  
+  for (unsigned int j=0; j<column_order.size(); ++j)
+    {
+      string key=column_order[j];
+      const map<string, vector<string> >::const_iterator 
+       super_iter=supercolumns.find(key);
+
+      if (super_iter!=supercolumns.end())
+       {
+         const unsigned int n_subcolumns=super_iter->second.size();
+                                          // avoid use of `tex_supercaptions[key]'
+         map<string,string>::const_iterator 
+           tex_super_cap_iter=tex_supercaptions.find(key);
+         out << endl << "\\multicolumn{" << n_subcolumns << "}{|c|}{" 
+             << tex_super_cap_iter->second << "}";
+       }
+      else
+       {
+                                          // col_iter->second=columns[col];
+         const map<string, Column>::const_iterator
+           col_iter=columns.find(key);
+         Assert(col_iter!=columns.end(), ExcInternalError());
+         out << col_iter->second.tex_caption;
+       }
+      if (j<column_order.size()-1)
+       out << " & ";
+    }
+  out << "\\\\ \\hline" << endl;
+
+                                  // write the n rows
+  const unsigned int n_rows=check_n_rows();
+  for (unsigned int i=0; i<n_rows; ++i)
+    {
+      const unsigned int n_cols=sel_columns.size();
+      
+      for (unsigned int j=0; j<n_cols; ++j)
+       {
+         string key=sel_columns[j];
+                                          // avoid `column[key]'
+         const map<string, Column>::const_iterator
+           col_iter=columns.find(key);
+         Assert(col_iter!=columns.end(), ExcInternalError());
+         
+         const Column &column=col_iter->second;
+         
+         out << setprecision(column.precision);
+         column.entries[i]->write_tex(out);
+         
+         if (j<n_cols-1)
+           out << " & ";
+       }
+      out << "\\\\ \\hline" << endl;
+    }
+
+  string caption="table";
+  
+  out << "\\end{tabular}" << endl
+      << "\\end{center}" << endl
+//      << "\\caption{"  << caption << "}" << endl
+      << "\\end{table}" << endl
+      << "\\end{document}" << endl;
+}
+
+
+unsigned int TableHandler::check_n_rows() const
+{
+  map<string, Column>::const_iterator col_iter=columns.begin();
+  unsigned int n=col_iter->second.entries.size();
+  string first_name=col_iter->first;
+
+  for (++col_iter; col_iter!=columns.end(); ++col_iter)
+    Assert(col_iter->second.entries.size()==n, ExcWrongNumberOfDataEntries(
+      col_iter->first, col_iter->second.entries.size(), first_name, n));
+
+  return n;
+}
+
+
+
+void TableHandler::get_selected_columns(vector<string> &sel_columns) const
+{
+  sel_columns.clear();
+  
+  for (unsigned int j=0; j<column_order.size(); ++j)
+    {
+      string key=column_order[j];
+      const map<string, vector<string> >::const_iterator 
+       super_iter=supercolumns.find(key);
+
+      if (super_iter!=supercolumns.end())
+       {
+                                          // i.e. key is a supercolumn key
+         const unsigned int n_subcolumns=super_iter->second.size();
+         for (unsigned int j=0; j<n_subcolumns; ++j)
+           {
+             const string sub_key=super_iter->second[j];
+             Assert(columns.count(sub_key), ExcInternalError());
+             sel_columns.push_back(super_iter->second[j]);
+           }
+       }
+      else
+       {
+         Assert(columns.count(key), ExcInternalError());
+                                          // i.e. key is a column key
+         sel_columns.push_back(key);
+       }
+    }
+}
+
+
+template class TableEntry<double>;
+template class TableEntry<float>;
+template class TableEntry<int>;
+template class TableEntry<unsigned int>;
+template class TableEntry<string>;
+
+
+
+template
+void TableHandler::add_value(const string &, double);
+
+template
+void TableHandler::add_value(const string &, int);
+
+template
+void TableHandler::add_value(const string &, unsigned int);
+
+template
+void TableHandler::add_value(const string &, string);
+
+

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.