--- /dev/null
+/*---------------------------- 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 ---------------------------*/
--- /dev/null
+/*---------------------------- 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 ---------------------------*/
--- /dev/null
+/*---------------------------- 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 ---------------------------*/
--- /dev/null
+/*---------------------------- 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);
+
+