From d7bb11b8fedb02f310a90b18742888e8946878c2 Mon Sep 17 00:00:00 2001 From: kronbichler Date: Fri, 10 May 2013 16:51:44 +0000 Subject: [PATCH] Fix one more issue, augment documentation of new parameter. git-svn-id: https://svn.dealii.org/trunk@29497 0785d39b-7218-0410-832d-ea1e28bc413d --- .../include/deal.II/base/convergence_table.h | 122 +++++++----------- deal.II/source/base/convergence_table.cc | 4 +- tests/base/convergence_table_04.cc | 18 +-- 3 files changed, 55 insertions(+), 89 deletions(-) diff --git a/deal.II/include/deal.II/base/convergence_table.h b/deal.II/include/deal.II/base/convergence_table.h index 3e6e069bb9..7ab30b21da 100644 --- a/deal.II/include/deal.II/base/convergence_table.h +++ b/deal.II/include/deal.II/base/convergence_table.h @@ -92,7 +92,11 @@ public: * * As this class has no information on the space dimension upon which the * reference column vs. the value column is based upon, it needs to be - * passed as last argument to this method. The default dimension is two. + * passed as last argument to this method. The default dimension for the + * reference column is 2, which is appropriate for the number of cells + * in 2D. If you work in 3D, set the number to 3. If the reference column is + * $1/h$, remember to set the dimension to 1 also when working in 3D to get + * correct rates. * * The new rate column and the data column will be merged to a * supercolumn. The tex caption of the supercolumn will be (by default) the @@ -103,9 +107,8 @@ public: * This method behaves in the following way: * * If RateMode is reduction_rate, then the computed output is - * $ \frac{e_{n-1}/k_{n-1}}{e_n/k_n} $. - * - * Where $k$ is the reference column. + * $ \frac{e_{n-1}/k_{n-1}}{e_n/k_n}, $ + * where $k$ is the reference column (no dimension dependence!). * * If RateMode is reduction_rate_log2, then the computed output is * $ @@ -113,8 +116,9 @@ public: * $. * * This is useful, for example, if we use as reference key the number of - * degrees of freedom. Assuming that the error is proportional to $ C - * (1/\sqrt{k})^r $, then this method will produce the rate $r$ as a result. + * degrees of freedom or better, the number of cells. Assuming that the + * error is proportional to $ C (1/\sqrt{k})^r $ in 2D, then this method + * will produce the rate $r$ as a result. * * @note Since this function adds columns to the table after several rows * have already been filled, it switches off the auto fill mode of the @@ -129,103 +133,65 @@ public: /** - * Evaluates the convergence rates of the - * data column data_column_key - * due to the #RateMode. - * Be sure that the value types of the - * table entries of the data column - * is a number, i.e. double, float, - * (unsigned) int, and so on. + * Evaluates the convergence rates of the data column + * data_column_key due to the #RateMode. Be sure that the value + * types of the table entries of the data column is a number, i.e. double, + * float, (unsigned) int, and so on. * - * The new rate column and the data column - * will be merged to a supercolumn. The - * tex caption of the supercolumn will be - * (by default) the same as the one of the - * data column. This may be changed by using - * the set_tex_supercaption() function - * of the base class TableHandler. + * The new rate column and the data column will be merged to a + * supercolumn. The tex caption of the supercolumn will be (by default) the + * same as the one of the data column. This may be changed by using the + * set_tex_supercaption() function of the base class TableHandler. * - * @note Since this function adds columns - * to the table after several rows have - * already been filled, it switches off - * the auto fill mode of the TableHandler - * base class. If you intend to add - * further data with auto fill, you will - * have to re-enable it after calling - * this function. + * @note Since this function adds columns to the table after several rows + * have already been filled, it switches off the auto fill mode of the + * TableHandler base class. If you intend to add further data with auto + * fill, you will have to re-enable it after calling this function. */ void evaluate_convergence_rates (const std::string &data_column_key, const RateMode rate_mode); /** - * Omit this column key - * (not supercolumn!) from the - * evaluation of the convergence rates - * of `all' columns (see the following - * two functions). + * Omit this column key (not supercolumn!) from the evaluation of + * the convergence rates of `all' columns (see the following two functions). * - * The Column::flag==1 is reserved for - * omitting the column from convergence + * The Column::flag==1 is reserved for omitting the column from convergence * rate evalution. */ void omit_column_from_convergence_rate_evaluation(const std::string &key); /** - * Evaluates convergence rates - * due to the rate_mode - * in relation to the reference - * column - * reference_column_key. This - * function evaluates the rates - * of ALL columns except of the - * columns that are to be omitted - * (see previous function) and - * execpt of the columns that are - * previously evaluated rate - * columns. This function allows - * to evaluate the convergence - * rate for almost all columns of - * a table without calling - * evaluate_convergence_rates() - * for each column separately. + * Evaluates convergence rates due to the rate_mode in relation to + * the reference column reference_column_key. This function + * evaluates the rates of ALL columns except of the columns that are to be + * omitted (see previous function) and execpt of the columns that are + * previously evaluated rate columns. This function allows to evaluate the + * convergence rate for almost all columns of a table without calling + * evaluate_convergence_rates() for each column separately. * * Example: - * Columns like n cells or - * n dofs columns may be wanted - * to be omitted in the evaluation - * of the convergence rates. Hence they - * should omitted by calling the - * omit_column_from_convergence_rate_evaluation(). + * Columns like n cells or n dofs columns may be wanted to + * be omitted in the evaluation of the convergence rates. Hence they should + * omitted by calling the omit_column_from_convergence_rate_evaluation(). */ void evaluate_all_convergence_rates(const std::string &reference_column_key, const RateMode rate_mode); /** - * Evaluates convergence rates - * due to the rate_mode. This - * function evaluates the rates of - * ALL columns except of the - * columns that are to be omitted - * (see previous function) - * and execpt of the columns that are - * previously - * evaluated rate columns. - * This function allows to evaluate - * the convergence rate for almost all - * columns of a table without calling - * evaluate_convergence_rates() - * for each column seperately. + * Evaluates convergence rates due to the rate_mode. This function + * evaluates the rates of ALL columns except of the columns that are to be + * omitted (see previous function) and execpt of the columns that are + * previously evaluated rate columns. This function allows to evaluate the + * convergence rate for almost all columns of a table without calling + * evaluate_convergence_rates() for each column seperately. * * Example: - * Columns like n cells or - * n dofs columns may be wanted - * to be omitted in the evaluation - * of the convergence rates. Hence they - * should omitted by calling the - * omit_column_from_convergence_rate_evaluation(). + * Columns like n cells or n dofs columns may be wanted to + * be omitted in the evaluation of the convergence rates. Hence they should + * omitted by calling the omit_column_from_convergence_rate_evaluation(). */ void evaluate_all_convergence_rates(const RateMode rate_mode); diff --git a/deal.II/source/base/convergence_table.cc b/deal.II/source/base/convergence_table.cc index 0fdf15a838..df3555ff05 100644 --- a/deal.II/source/base/convergence_table.cc +++ b/deal.II/source/base/convergence_table.cc @@ -40,7 +40,7 @@ void ConvergenceTable::evaluate_convergence_rates(const std::string &data_column std::vector &entries=columns[data_column_key].entries; std::vector &ref_entries=columns[reference_column_key].entries; - std::string rate_key = data_column_key; + std::string rate_key = data_column_key+"..."; const unsigned int n = entries.size(); const unsigned int n_ref = ref_entries.size(); @@ -80,7 +80,7 @@ void ConvergenceTable::evaluate_convergence_rates(const std::string &data_column } break; case reduction_rate_log2: - rate_key += "red.rate"; + rate_key += "red.rate.log2"; no_rate_entries = columns[rate_key].entries.size(); // Calculate all missing rate values: for (unsigned int i = no_rate_entries; i