]> https://gitweb.dealii.org/ - dealii.git/commitdiff
exclude_from_rates -> compute_rate
authorLuca Heltai <luca.heltai@sissa.it>
Wed, 1 May 2019 13:11:52 +0000 (15:11 +0200)
committerLuca Heltai <luca.heltai@sissa.it>
Wed, 1 May 2019 13:22:07 +0000 (15:22 +0200)
include/deal.II/base/parsed_convergence_table.h
source/base/parsed_convergence_table.cc
tests/base/parsed_convergence_table_10.cc

index 2132019ab632e71eb2b2396532b4a44a8667de67..32b14fe04a4dd73bb59ae287dca06740313027a7 100644 (file)
@@ -222,7 +222,7 @@ public:
    * # Number of digits to use when printing the error.
    * set Error precision                  = 3
    *
-   * # Extra columns to add to the table. Available options are dofs and cells
+   * # Extra columns to add to the table. Available options are dofs and cells.
    * set Extra columns                    = dofs, cells
    *
    * # The exponent to use when computing p-norms.
@@ -306,7 +306,7 @@ public:
    * error_from_exact() or difference().
    *
    * Make sure you add all extra columns before the first call to
-   * error_from_exact() or difference(), since adding additional columns to the
+   * error_from_exact() or difference(). Adding additional columns to the
    * convergence table after you already started filling the table will trigger
    * an exception.
    *
@@ -321,7 +321,7 @@ public:
    *        return dt;
    * };
    *
-   * table.add_extra_column("dt", dt_function,  true);
+   * table.add_extra_column("dt", dt_function, false);
    *
    * for (unsigned int i = 0; i < n_cycles; ++i)
    *   {
@@ -343,9 +343,9 @@ public:
    * provided that you use the following parameter file (only non default
    * entries are shown here):
    * @code
-   * set Extra columns =
-   * set List of error norms to compute   = L2_norm
-   * set Rate key                         = dt
+   * set Extra columns                  =
+   * set List of error norms to compute = L2_norm
+   * set Rate key                       = dt
    * @endcode
    *
    *
@@ -353,13 +353,16 @@ public:
    * @param custom_function Function that will be called to fill the given
    * entry. You need to make sure that the scope of this function is valid
    * up to the call to error_from_exact() or difference();
-   * @param exclude_from_rates If set to true, then no error rates will be
-   * computed for this column.
+   * @param compute_rate If set to true, then this column will be included in
+   * the list of columns for which error rates are computed. You may want to set
+   * this to false if you want to compute error rates with respect to this
+   * column. In this case, you should also specify @p column_name as the rate
+   * key in the parameter file.
    */
   void
   add_extra_column(const std::string &            column_name,
                    const std::function<double()> &custom_function,
-                   const bool &                   exclude_from_rates = false);
+                   const bool &                   compute_rate = true);
 
   /**
    * Difference between two solutions in the same vector space.
@@ -434,7 +437,7 @@ private:
   std::vector<std::set<VectorTools::NormType>> norms_per_unique_component;
 
   /**
-   * Exponent to use in the p-norm types.
+   * Exponent to use in p-norm types.
    */
   double exponent;
 
@@ -469,8 +472,8 @@ private:
   std::string error_file_name;
 
   /**
-   * Compute the error. If this is false, all functions regarding
-   * errors are disabled and don't do anything.
+   * Compute the error. If this is false, all methods that perform the
+   * computation of the error are disabled and don't do anything.
    */
   bool compute_error;
 };
index c4e1445788d8f5177c8fd777acc87222599af626..72761b30655f1a0b00c494f01cf6c19b54e841ab 100644 (file)
@@ -171,12 +171,11 @@ ParsedConvergenceTable::prepare_table_for_output()
                 has_key = true;
 
               if (col != "")
-                table.omit_column_from_convergence_rate_evaluation(
-                  Patterns::Tools::to_string(col));
+                table.omit_column_from_convergence_rate_evaluation(col);
             }
 
           for (const auto &extra_col : extra_column_functions)
-            if (extra_col.second.second)
+            if (extra_col.second.second == false)
               {
                 if (rate_key == extra_col.first)
                   has_key = true;
@@ -188,12 +187,10 @@ ParsedConvergenceTable::prepare_table_for_output()
             {
               if (rate_mode == "reduction_rate_log2")
                 table.evaluate_all_convergence_rates(
-                  Patterns::Tools::to_string(rate_key),
-                  ConvergenceTable::reduction_rate_log2);
+                  rate_key, ConvergenceTable::reduction_rate_log2);
               else if (rate_mode == "reduction_rate")
                 table.evaluate_all_convergence_rates(
-                  Patterns::Tools::to_string(rate_key),
-                  ConvergenceTable::reduction_rate);
+                  rate_key, ConvergenceTable::reduction_rate);
               else
                 {
                   Assert(rate_mode != "none", ExcInternalError());
@@ -264,9 +261,9 @@ void
 ParsedConvergenceTable::add_extra_column(
   const std::string &            column_name,
   const std::function<double()> &custom_function,
-  const bool &                   exclude_from_rates)
+  const bool &                   compute_rate)
 {
-  extra_column_functions[column_name] = {custom_function, exclude_from_rates};
+  extra_column_functions[column_name] = {custom_function, compute_rate};
 }
 
 DEAL_II_NAMESPACE_CLOSE
index 834f99592330beb8b06aef09a7d56785b0503525..d217a138df51e45ee0cfdea610085950b2d47c26 100644 (file)
@@ -66,7 +66,7 @@ main()
       auto dt    = [&]() { return i + 1.0; };
 
       table.add_extra_column("cycle", cycle);
-      table.add_extra_column("dt", dt, true);
+      table.add_extra_column("dt", dt, false);
       table.error_from_exact(dh, sol, exact);
     }
 

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.