1 // --------------------------------------------------------------------- 2 // 3 // Copyright (C) 1999 - 2018 by the deal.II authors 4 // 5 // This file is part of the deal.II library. 6 // 7 // The deal.II library is free software; you can use it, redistribute 8 // it, and/or modify it under the terms of the GNU Lesser General 9 // Public License as published by the Free Software Foundation; either 10 // version 2.1 of the License, or (at your option) any later version. 11 // The full text of the license can be found in the file LICENSE.md at 12 // the top level directory of deal.II. 13 // 14 // --------------------------------------------------------------------- 15 16 #ifndef dealii_convergence_table_h 17 #define dealii_convergence_table_h 18 19 20 #include <deal.II/base/config.h> 21 22 #include <deal.II/base/table_handler.h> 23 24 DEAL_II_NAMESPACE_OPEN 25 26 27 /** 28 * The ConvergenceTable class is an application to the TableHandler class and 29 * stores some convergence data, such as residuals of the cg-method, or some 30 * evaluated <i>L<sup>2</sup></i>-errors of discrete solutions, etc, and 31 * evaluates convergence rates or orders. 32 * 33 * The already implemented #RateMode's are #reduction_rate, where the 34 * convergence rate is the quotient of two following rows, and 35 * #reduction_rate_log2, that evaluates the order of convergence. These 36 * standard evaluations are useful for global refinement, for local refinement 37 * this may not be an appropriate method, as the convergence rates should be 38 * set in relation to the number of cells or the number of DoFs. The 39 * implementations of these non-standard methods is left to a user. 40 * 41 * The number of cells and the number of DoFs may be added to the table by 42 * calling e.g. <tt>add_value("n cells", n_cells)</tt>. The table data is 43 * also added by calling add_value(). Before the output of the table the 44 * functions evaluate_convergence_rates() and evaluate_all_convergence_rates() 45 * may be called. 46 * 47 * There are two possibilities of how to evaluate the convergence rates of 48 * multiple columns in the same RateMode. 49 * <ol> 50 * <li> call evaluate_convergence_rates() for all wanted columns 51 * <li> call omit_column_from_convergence_rate_evaluation() for all columns 52 * for which this evaluation is not desired and then 53 * evaluate_all_convergence_rates() to evaluate the convergence rates of all 54 * columns that have not been flagged for omission. 55 * </ol> 56 * 57 * A detailed discussion of this class can also be found in the step-7 and 58 * step-13 example programs. 59 * 60 * @ingroup textoutput 61 */ 62 class ConvergenceTable : public TableHandler 63 { 64 public: 65 /** 66 * Constructor. 67 */ 68 ConvergenceTable() = default; 69 70 /** 71 * Rate in relation to the rows. 72 */ 73 enum RateMode 74 { 75 /** 76 * Do not do anything. 77 */ 78 none, 79 /** 80 * Quotient of values in the previous row and in this row. 81 */ 82 reduction_rate, 83 /** 84 * Logarithm of #reduction_rate to the base 2 representing the order of 85 * convergence when halving the grid size, e.g. from h to h/2. 86 */ 87 reduction_rate_log2 88 }; 89 90 /** 91 * Evaluate the convergence rates of the data column 92 * <tt>data_column_key</tt> due to the #RateMode in relation to the 93 * reference column <tt>reference_column_key</tt>. Be sure that the value 94 * types of the table entries of the data column and the reference data 95 * column is a number, i.e. double, float, (unsigned) int, and so on. 96 * 97 * As this class has no information on the space dimension upon which the 98 * reference column vs. the value column is based upon, it needs to be 99 * passed as last argument to this method. The <i>default dimension for the 100 * reference column</i> is 2, which is appropriate for the number of cells 101 * in 2D. If you work in 3D, set the number to 3. If the reference column is 102 * $1/h$, remember to set the dimension to 1 also when working in 3D to get 103 * correct rates. 104 * 105 * The new rate column and the data column will be merged to a supercolumn. 106 * The tex caption of the supercolumn will be (by default) the same as the 107 * one of the data column. This may be changed by using the 108 * <tt>set_tex_supercaption (...)</tt> function of the base class 109 * TableHandler. 110 * 111 * This method behaves in the following way: 112 * 113 * If RateMode is reduction_rate, then the computed output is $ 114 * \frac{e_{n-1}/k_{n-1}}{e_n/k_n}, $ where $k$ is the reference column (no 115 * dimension dependence!). 116 * 117 * If RateMode is reduction_rate_log2, then the computed output is $ dim 118 * \frac{\log |e_{n-1}/e_{n}|}{\log |k_n/k_{n-1}|} $. 119 * 120 * This is useful, for example, if we use as reference key the number of 121 * degrees of freedom or better, the number of cells. Assuming that the 122 * error is proportional to $ C (1/\sqrt{k})^r $ in 2D, then this method 123 * will produce the rate $r$ as a result. For general dimension, as 124 * described by the last parameter of this function, the formula needs to be 125 * $ C (1/\sqrt[dim]{k})^r $. 126 * 127 * @note Since this function adds columns to the table after several rows 128 * have already been filled, it switches off the auto fill mode of the 129 * TableHandler base class. If you intend to add further data with auto 130 * fill, you will have to re-enable it after calling this function. 131 */ 132 void 133 evaluate_convergence_rates(const std::string &data_column_key, 134 const std::string &reference_column_key, 135 const RateMode rate_mode, 136 const unsigned int dim = 2); 137 138 139 /** 140 * Evaluate the convergence rates of the data column 141 * <tt>data_column_key</tt> due to the #RateMode. Be sure that the value 142 * types of the table entries of the data column is a number, i.e. double, 143 * float, (unsigned) int, and so on. 144 * 145 * The new rate column and the data column will be merged to a supercolumn. 146 * The tex caption of the supercolumn will be (by default) the same as the 147 * one of the data column. This may be changed by using the 148 * set_tex_supercaption() function of the base class TableHandler. 149 * 150 * @note Since this function adds columns to the table after several rows 151 * have already been filled, it switches off the auto fill mode of the 152 * TableHandler base class. If you intend to add further data with auto 153 * fill, you will have to re-enable it after calling this function. 154 */ 155 void 156 evaluate_convergence_rates(const std::string &data_column_key, 157 const RateMode rate_mode); 158 159 /** 160 * Omit this column <tt>key</tt> (not supercolumn!) from the evaluation of 161 * the convergence rates of `all' columns (see the following two functions). 162 * 163 * The Column::flag==1 is reserved for omitting the column from convergence 164 * rate evaluation. 165 */ 166 void 167 omit_column_from_convergence_rate_evaluation(const std::string &key); 168 169 /** 170 * Evaluate convergence rates due to the <tt>rate_mode</tt> in relation to 171 * the reference column <tt>reference_column_key</tt>. This function 172 * evaluates the rates of ALL columns except of the columns that are to be 173 * omitted (see previous function) and except of the columns that are 174 * previously evaluated rate columns. This function allows to evaluate the 175 * convergence rate for almost all columns of a table without calling 176 * evaluate_convergence_rates() for each column separately. 177 * 178 * Example: Columns like <tt>n cells</tt> or <tt>n dofs</tt> columns may be 179 * wanted to be omitted in the evaluation of the convergence rates. Hence 180 * they should omitted by calling the 181 * omit_column_from_convergence_rate_evaluation(). 182 */ 183 void 184 evaluate_all_convergence_rates(const std::string &reference_column_key, 185 const RateMode rate_mode); 186 187 /** 188 * Evaluate convergence rates due to the <tt>rate_mode</tt>. This function 189 * evaluates the rates of ALL columns except of the columns that are to be 190 * omitted (see previous function) and except of the columns that are 191 * previously evaluated rate columns. This function allows to evaluate the 192 * convergence rate for almost all columns of a table without calling 193 * evaluate_convergence_rates() for each column separately. 194 * 195 * Example: Columns like <tt>n cells</tt> or <tt>n dofs</tt> columns may be 196 * wanted to be omitted in the evaluation of the convergence rates. Hence 197 * they should omitted by calling the 198 * omit_column_from_convergence_rate_evaluation(). 199 */ 200 void 201 evaluate_all_convergence_rates(const RateMode rate_mode); 202 203 /** 204 * @addtogroup Exceptions 205 * @{ 206 */ 207 208 /** 209 * Exception 210 */ 211 DeclException1(ExcRateColumnAlreadyExists, 212 std::string, 213 << "Rate column <" << arg1 << "> does already exist."); 214 //@} 215 }; 216 217 218 DEAL_II_NAMESPACE_CLOSE 219 220 #endif 221