From 67deb9d37d48f6a769f2f36c12732e73689ad02e Mon Sep 17 00:00:00 2001 From: Luca Heltai Date: Wed, 10 Apr 2019 00:58:34 +0200 Subject: [PATCH] ParsedConvergenceTable --- .../deal.II/base/parsed_convergence_table.h | 766 ++++++++++++++++++ .../base/parsed_convergence_table_flags.h | 232 ++++++ source/base/CMakeLists.txt | 1 + source/base/parsed_convergence_table.cc | 201 +++++ 4 files changed, 1200 insertions(+) create mode 100644 include/deal.II/base/parsed_convergence_table.h create mode 100644 include/deal.II/base/parsed_convergence_table_flags.h create mode 100644 source/base/parsed_convergence_table.cc diff --git a/include/deal.II/base/parsed_convergence_table.h b/include/deal.II/base/parsed_convergence_table.h new file mode 100644 index 0000000000..669fd6ebb6 --- /dev/null +++ b/include/deal.II/base/parsed_convergence_table.h @@ -0,0 +1,766 @@ +//----------------------------------------------------------- +// +// Copyright (C) 2019 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal2lkit library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal2lkit distribution. +// +//----------------------------------------------------------- + +#ifndef dealii_base_parsed_convergence_table_h +#define dealii_base_parsed_convergence_table_h + +#include + +#include +#include +#include +#include +#include +#include +#include + +#include + +#include + +#include +#include + +#include + +#include + +DEAL_II_NAMESPACE_OPEN + +/** + * @brief The ParsedConvergenceTable class + * + * This class simplifies the construction of convergence tables, reading the + * options for the generation of the table from a parameter file. It provides a + * series of methods that can be used to compute the error given a reference + * exact solution, or the difference between two numerical solutions, or any + * other custom computation of the error, given via std::function objects. + * + * An example usage of this class is given by + * @code + * ParsedConvergenceTable table; + * + * ParameterHandler prm; + * table.add_parameters(prm); + * + * for(unsigned int i=0; i &component_names = {"u"}, + const std::vector> + &list_of_error_norms = {{ParsedConvergenceTableFlags::H1, + ParsedConvergenceTableFlags::L2, + ParsedConvergenceTableFlags::Linfty}}); + + /** + * Full constructor for ParsedConvergenceTable. + * + * @param component_names Names of the components. Repeated consecutive + * names are interpreted as components of a vector valued field; + * @param list_of_error_norms Specify what error norms to compute for each + * unique component name; + * @param extra_columns Extra columns to add; + * @param rate_key Specify the extra column by which we will compute the + * error rates; + * @param error_file_name Name of error output file (with extension txt, + * gpl, tex, or org). If different from the empty string, than + * output_table() also writes in this file in the format deduced from its + * extension; + * @param precision How many digits to use when writing the error; + * @param compute_error If set to `false`, then no error computation is + * performed, and no output is produced; + * @param output_error If set to `false`, then the call to output_table() + * will produce only an output file (if @p error_file_name is not empty), + * but no output will be written on the given stream. + * + * The parameters you specify with this constructor can be written to a + * ParameterHandler object by calling the add_parameters() method. Once you + * call the add_parameters() method, the following options will be defined in + * the given ParameterHandler object, and the parameters of the instance of + * this class will follow the modification you make to the ParameterHandler + * object at run time: + * @code + * # Listing of Parameters + * # --------------------- + * # When set to false, no computations are performed. + * set Enable computation of the errors = true + * + * # When set to false, printing of the convergence table to the stream + * # specified as input to output_table() is disabled. + * set Enable output to streams = true + * + * # Set this to a filename with extension .txt, .gpl, .org, or .tex to enable + * # writing the convergence table to a file. + * set Error file name = + * + * # Number of digits to use when printing the error. + * set Error precision = 3 + * + * # Extra columns to add to the table. Availabl options are dofs, cells, and + * # dt. + * set Extra columns = dofs, cells + * + * # Each component is separated by a semicolon, and each norm by a comma. + * # Implemented norms are Linfty, L2, W1infty, H1, and custom. If you want to + * # skip a component, use none. + * set List of error norms to compute = Linfty, L2, H1 + * + * # Key to use when computing convergence rates. If this is set to a + * # column that is not present, or to none, then no error rates are computed. + * set Rate key = dofs + * @endcode + */ + ParsedConvergenceTable( + const std::vector &component_names, + const std::vector> + &list_of_error_norms, + const std::set &extra_columns, + const ParsedConvergenceTableFlags::ExtraColumns & rate_key, + const std::string & error_file_name, + const unsigned int & precision, + const bool & compute_error, + const bool & output_error); + + /** + * Attach all the parameters in this class to entries of the parameter + * handler @p prm. Whenever the content of @p prm changes, the parameters + * of this class will be updated. + */ + void + add_parameters(ParameterHandler &prm); + + /** + * Add a row to the error table, containing the error between @p solution and + * the @p exact function, in the norm(s) specified in the parameter file. + * + * If you specify a @p weight expression during this call, then this is used + * to compute weighted errors. If the weight is different from "1.0", then + * you have to make sure that muparser is installed and enabled, otherwise an + * exception will be thrown. + */ + template + void + error_from_exact(const DoFHandlerType & vspace, + const VectorType & solution, + const Function &exact, + double dt = 0., + const std::string &weight = "1.0"); + + /** + * Same as above, with a different mapping. + */ + template + void + error_from_exact(const Mapping & mapping, + const DoFHandlerType & vspace, + const VectorType & solution, + const Function &exact, + double dt = 0., + const std::string &weight = "1.0"); + + + /** + * Call the given custom function to compute a custom error for any component + * for which you specify the ConvergenceTableFlags::custom norm. This function + * may be called several times to add different columns with different names, + * provided you use an appropriate function every time, and + * that you update @p column_name everytime. + * + * Should you wish to do so, then make sure the parameter + * @p add_table_extras is set to false for each call, except one, so + * that you only add extra columns informations once. + * + * An example usage of this class with a custom error function is the + * following: + * @code + * using namespace ParsedConvergenceTableFlags; + * ParsedConvergenceTable table({"u"}, {{"custom"}}); + * + * for(unsigned int i=0; i + void + custom_error(const std::function + & custom_error_function, + const DoFHandlerType &dh, + const std::string & column_name = "custom", + const bool add_table_extras = false, + const double dt = 0.); + + /** + * Difference between two solutions in the same vector space. + */ + template + void + difference(const DoFHandlerType &, + const VectorType &, + const VectorType &, + double dt = 0.); + + /** + * Same as above, with a non default mapping. + */ + template + void + difference(const Mapping &mapping, + const DoFHandlerType &, + const VectorType &, + const VectorType &, + double dt = 0.); + + /** + * Write the error table to the @p out stream, and (possibly) to the file + * stream specified in the parameters. + */ + void + output_table(std::ostream &out = std::cout); + +private: + /** + * Names of the solution components. + */ + const std::vector component_names; + + /** + * Same as above, but containing repeated component names only once. + */ + const std::vector unique_component_names; + + /** + * Masks for each unique component name. + */ + const std::vector unique_component_masks; + + /** + * Type of error to compute per components. + */ + std::vector> + norms_per_unique_component; + + /** + * The actual table + */ + ConvergenceTable table; + /** + * Extra columns to add to the table. + */ + std::set extra_columns; + + /** + * Wether or not to calculate the rates according to the given keys. + */ + ParsedConvergenceTableFlags::ExtraColumns rate_key; + + /** + * The precision used to output the table. + */ + unsigned int precision; + + /** + * Filename to use when writing to file. + */ + std::string error_file_name; + + /** + * Compute the error. If this is false, all functions regarding + * errors are disabled and don't do anything. + */ + bool compute_error; + + /** + * Output the error file also on screen. + */ + bool output_error; +}; + + + +#ifndef DOXYGEN +// ============================================================ +// Template instantiations +// ============================================================ +template +void +ParsedConvergenceTable::difference(const DoFHandlerType &dh, + const VectorType & solution1, + const VectorType & solution2, + double dt) +{ + AssertThrow(solution1.size() == solution2.size(), + ExcDimensionMismatch(solution1.size(), solution2.size())); + VectorType solution(solution1); + solution -= solution2; + error_from_exact(dh, + solution, + ConstantFunction( + 0, component_names.size()), + dt); +} + + + +template +void +ParsedConvergenceTable::difference( + const Mapping + & mapping, + const DoFHandlerType &dh, + const VectorType & solution1, + const VectorType & solution2, + double dt) +{ + AssertThrow(solution1.size() == solution2.size(), + ExcDimensionMismatch(solution1.size(), solution2.size())); + VectorType solution(solution1); + solution -= solution2; + error_from_exact(mapping, + dh, + solution, + ConstantFunction( + 0, component_names.size()), + dt); +} + + + +template +void +ParsedConvergenceTable::error_from_exact( + const DoFHandlerType & dh, + const VectorType & solution, + const Function &exact, + double dt, + const std::string & weight) +{ + error_from_exact(StaticMappingQ1::mapping, + dh, + solution, + exact, + dt, + weight); +} + + + +template +void +ParsedConvergenceTable::error_from_exact( + const Mapping + & mapping, + const DoFHandlerType & dh, + const VectorType & solution, + const Function &exact, + double dt, + const std::string & weight) +{ + const int dim = DoFHandlerType::dimension; + const int spacedim = DoFHandlerType::space_dimension; + if (compute_error) + { + AssertDimension(exact.n_components, component_names.size()); + + const unsigned int n_active_cells = + dh.get_triangulation().n_global_active_cells(); + const unsigned int n_dofs = dh.n_dofs(); + + for (const auto &col : extra_columns) + if (col == ParsedConvergenceTableFlags::cells) + { + table.add_value("cells", n_active_cells); + table.set_tex_caption("cells", "\\# cells"); + table.set_tex_format("cells", "r"); + } + else if (col == ParsedConvergenceTableFlags::dofs) + { + table.add_value("dofs", n_dofs); + table.set_tex_caption("dofs", "\\# dofs"); + table.set_tex_format("dofs", "r"); + } + else if (col == ParsedConvergenceTableFlags::dt) + { + table.add_value("dt", dt); + table.set_tex_caption("dt", "\\Delta t"); + table.set_tex_format("dt", "r"); + } + + + std::map + norm_to_index; + norm_to_index[ParsedConvergenceTableFlags::L2] = 0; + norm_to_index[ParsedConvergenceTableFlags::H1] = 1; + norm_to_index[ParsedConvergenceTableFlags::Linfty] = 2; + norm_to_index[ParsedConvergenceTableFlags::W1infty] = 3; + norm_to_index[ParsedConvergenceTableFlags::custom] = 4; + + const auto n_components = component_names.size(); + for (unsigned int i = 0; i < unique_component_names.size(); ++i) + { + std::vector error(5, 0.0); + const auto & norms = norms_per_unique_component[i]; + const auto & mask = unique_component_masks[i]; + +# ifndef DEAL_II_WITH_MUPARSER + Assert(weight == "1.0", + ExcMessage("If you want to use a weight which is " + "different from 1.0, you have to " + "install deal.II with muparser enabled.")); + ComponentSelectFunction select_component( + {mask.first_selected_component(), + mask.first_selected_component() + mask.n_selected_components()}, + n_components); +# else + std::string comp_sel; + { + // Select all required components + std::vector expr(n_components, "0"); + + for (unsigned int i = 0; i < n_components; ++i) + if (mask[i] == true) + expr[i] = weight; + + comp_sel = expr[0]; + for (unsigned int i = 1; i < expr.size(); ++i) + comp_sel += "; " + expr[i]; + } + FunctionParser select_component(comp_sel); +# endif + + Vector difference_per_cell( + dh.get_triangulation().n_global_active_cells()); + + QGauss q_gauss((dh.get_fe().degree + 1) * 2); + + for (const auto &norm : norms) + { + if (norm == ParsedConvergenceTableFlags::L2) + { + VectorTools::integrate_difference(mapping, + dh, + solution, + exact, + difference_per_cell, + q_gauss, + VectorTools::L2_norm, + &select_component); + } + + error[norm_to_index[ParsedConvergenceTableFlags::L2]] = + VectorTools::compute_global_error(dh.get_triangulation(), + difference_per_cell, + VectorTools::L2_norm); + difference_per_cell = 0; + + if (norm == ParsedConvergenceTableFlags::H1) + { + VectorTools::integrate_difference(mapping, + dh, + solution, + exact, + difference_per_cell, + q_gauss, + VectorTools::H1_norm, + &select_component); + } + error[norm_to_index[ParsedConvergenceTableFlags::H1]] = + VectorTools::compute_global_error(dh.get_triangulation(), + difference_per_cell, + VectorTools::H1_norm); + difference_per_cell = 0; + + if (norm == ParsedConvergenceTableFlags::W1infty) + { + VectorTools::integrate_difference(mapping, + dh, + solution, + exact, + difference_per_cell, + q_gauss, + VectorTools::W1infty_norm, + &select_component); + } + + error[norm_to_index[ParsedConvergenceTableFlags::W1infty]] = + VectorTools::compute_global_error( + dh.get_triangulation(), + difference_per_cell, + VectorTools::W1infty_seminorm); + difference_per_cell = 0; + + if (norm == ParsedConvergenceTableFlags::Linfty) + { + VectorTools::integrate_difference(mapping, + dh, + solution, + exact, + difference_per_cell, + q_gauss, + VectorTools::Linfty_norm, + &select_component); + } + + error[norm_to_index[ParsedConvergenceTableFlags::Linfty]] = + VectorTools::compute_global_error(dh.get_triangulation(), + difference_per_cell, + VectorTools::Linfty_norm); + difference_per_cell = 0; + + // Now add what we just computed. + if (norm != ParsedConvergenceTableFlags::none) + { + std::string name = unique_component_names[i] + "_" + + Patterns::Tools::to_string(norm); + std::string latex_name = + "$\\| " + unique_component_names[i] + " - " + + unique_component_names[i] + "_h \\|_{" + + Patterns::Tools::to_string(norm) + "}$"; + + table.add_value(name, error[norm_to_index[norm]]); + table.set_precision(name, precision); + table.set_scientific(name, true); + table.set_tex_caption(name, latex_name); + } + } + } + } +} + + +template +void +ParsedConvergenceTable::custom_error( + const std::function + & custom_error_function, + const DoFHandlerType &dh, + const std::string & column_name, + const bool add_table_extras, + const double dt) +{ + if (compute_error) + { + const unsigned int n_components = norms_per_unique_component.size(); + std::vector c_error(norms_per_unique_component.size()); + const unsigned int n_active_cells = + dh.get_triangulation().n_global_active_cells(); + const unsigned int n_dofs = dh.n_dofs(); + if (add_table_extras) + { + for (const auto &col : extra_columns) + { + if (col == ParsedConvergenceTableFlags::cells) + { + table.add_value("cells", n_active_cells); + table.set_tex_caption("cells", "\\# cells"); + table.set_tex_format("cells", "r"); + } + if (col == ParsedConvergenceTableFlags::dofs) + { + table.add_value("dofs", n_dofs); + table.set_tex_caption("dofs", "\\# dofs"); + table.set_tex_format("dofs", "r"); + } + if (col == ParsedConvergenceTableFlags::dt) + { + table.add_value("dt", dt); + table.set_tex_caption("dt", "\\Delta t"); + table.set_tex_format("dt", "r"); + } + } + } + + for (unsigned int j = 0; j < n_components; ++j) + { + const auto &norms = norms_per_unique_component[j]; + + for (const auto &norm : norms) + if (norm == ParsedConvergenceTableFlags::custom) + { + const double custom_error = custom_error_function(j); + + std::string name = + unique_component_names[j] + "_" + column_name; + std::string latex_name = "$\\| " + unique_component_names[j] + + " - " + unique_component_names[j] + + "_h \\|_{\\text{" + column_name + + "} $"; + + table.add_value(name, custom_error); + table.set_precision(name, precision); + table.set_scientific(name, true); + table.set_tex_caption(name, latex_name); + } + } + } +} +#endif + +DEAL_II_NAMESPACE_CLOSE + +#endif diff --git a/include/deal.II/base/parsed_convergence_table_flags.h b/include/deal.II/base/parsed_convergence_table_flags.h new file mode 100644 index 0000000000..7ff300d104 --- /dev/null +++ b/include/deal.II/base/parsed_convergence_table_flags.h @@ -0,0 +1,232 @@ +//----------------------------------------------------------- +// +// Copyright (C) 2019 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal2lkit library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal2lkit distribution. +// +//----------------------------------------------------------- + +#ifndef dealii_base_parsed_convergence_table_flags_h +#define dealii_base_parsed_convergence_table_flags_h + +#include + +#include + +DEAL_II_NAMESPACE_OPEN + +namespace ParsedConvergenceTableFlags +{ + /** + * Implemented Error norms. + */ + enum NormType + { + none = 0x00, + Linfty = 0x01, + L2 = 0x02, + W1infty = 0x04, + H1 = 0x08, + custom = 0x10, + }; + + /** + * The ExtraColumns enum. + */ + enum ExtraColumns + { + no_columns = 0x0, + dofs = 0x1, + cells = 0x2, + dt = 0x4 + }; +} // namespace ParsedConvergenceTableFlags + +namespace Patterns +{ + namespace Tools + { + template <> + struct Convert + { + /** + * Return the Correct pattern for NormType. + */ + static std::unique_ptr + to_pattern() + { + return std::make_unique( + "none|Linfty|L2|W1infty|H1|custom"); + } + + + + /** + * Convert a NormType to a string. + */ + static std::string + to_string(const ParsedConvergenceTableFlags::NormType & s, + const std::unique_ptr &p = + Convert::to_pattern()) + { + (void)p; + if (s == ParsedConvergenceTableFlags::Linfty) + return "Linfty"; + else if (s == ParsedConvergenceTableFlags::L2) + return "L2"; + else if (s == ParsedConvergenceTableFlags::W1infty) + return "W1infty"; + else if (s == ParsedConvergenceTableFlags::H1) + return "H1"; + else if (s == ParsedConvergenceTableFlags::custom) + return "custom"; + else if (s == ParsedConvergenceTableFlags::none) + return "none"; + else + { + AssertThrow(false, ExcMessage("Didn't recognize a norm type.")); + return ""; + } + } + + + + /** + * Convert a string to a NormType. + */ + static ParsedConvergenceTableFlags::NormType + to_value(const std::string & norm_string, + const std::unique_ptr &p = + Convert::to_pattern()) + { + (void)p; + ParsedConvergenceTableFlags::NormType norm = + ParsedConvergenceTableFlags::none; + + if (norm_string == "Linfty") + { + norm = ParsedConvergenceTableFlags::Linfty; + } + else if (norm_string == "L2") + { + norm = ParsedConvergenceTableFlags::L2; + } + else if (norm_string == "W1infty") + { + norm = ParsedConvergenceTableFlags::W1infty; + } + else if (norm_string == "H1") + { + norm = ParsedConvergenceTableFlags::H1; + } + else if (norm_string == "custom") + { + norm = ParsedConvergenceTableFlags::custom; + } + else if (norm_string == "none") + { + norm = ParsedConvergenceTableFlags::none; + } + else + { + AssertThrow(false, ExcMessage("Didn't recognize a norm type.")); + } + + return norm; + } + }; + + + + template <> + struct Convert + { + /** + * Return the Correct pattern for ExtraColumns. + */ + static std::unique_ptr + to_pattern() + { + return std::make_unique("none|cells|dofs|dt"); + } + + + + /** + * Convert an ExtraColumn object to a string. + */ + static std::string + to_string( + const ParsedConvergenceTableFlags::ExtraColumns &s, + const std::unique_ptr & p = + Convert::to_pattern()) + { + (void)p; + if (s == ParsedConvergenceTableFlags::cells) + return "cells"; + else if (s == ParsedConvergenceTableFlags::dofs) + return "dofs"; + else if (s == ParsedConvergenceTableFlags::dt) + return "dt"; + else if (s == ParsedConvergenceTableFlags::no_columns) + return "none"; + else + { + AssertThrow(false, ExcMessage("Didn't recognize a column type.")); + return ""; + } + } + + + + /** + * Convert a string to an ExtraColumn object. + */ + static ParsedConvergenceTableFlags::ExtraColumns + to_value( + const std::string & column_string, + const std::unique_ptr &p = + Convert::to_pattern()) + { + (void)p; + ParsedConvergenceTableFlags::ExtraColumns column = + ParsedConvergenceTableFlags::no_columns; + + if (column_string == "cells") + { + column = ParsedConvergenceTableFlags::cells; + } + else if (column_string == "dofs") + { + column = ParsedConvergenceTableFlags::dofs; + } + else if (column_string == "dt") + { + column = ParsedConvergenceTableFlags::dt; + } + else if (column_string == "none") + { + column = ParsedConvergenceTableFlags::no_columns; + } + else + { + AssertThrow(false, + ExcMessage("Didn't recognize an extra column type.")); + } + + return column; + } + }; + } // namespace Tools +} // namespace Patterns + +DEAL_II_NAMESPACE_CLOSE + +#endif diff --git a/source/base/CMakeLists.txt b/source/base/CMakeLists.txt index eacda918d5..e3c3f8ee0a 100644 --- a/source/base/CMakeLists.txt +++ b/source/base/CMakeLists.txt @@ -51,6 +51,7 @@ SET(_unity_include_src parallel.cc parameter_handler.cc parameter_acceptor.cc + parsed_convergence_table.cc parsed_function.cc partitioner.cc patterns.cc diff --git a/source/base/parsed_convergence_table.cc b/source/base/parsed_convergence_table.cc new file mode 100644 index 0000000000..0139d976d2 --- /dev/null +++ b/source/base/parsed_convergence_table.cc @@ -0,0 +1,201 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2019 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + +#include +#include +#include + +DEAL_II_NAMESPACE_OPEN + +namespace +{ + std::vector + get_unique_component_names(const std::vector &component_names) + { + auto elements = component_names; + auto last = std::unique(elements.begin(), elements.end()); + elements.resize(last - elements.begin()); + return elements; + } + + + + std::vector + get_unique_component_masks(const std::vector &component_names) + { + const auto unique_component_names = + get_unique_component_names(component_names); + + std::vector masks; + + std::vector> bools( + unique_component_names.size(), + std::vector(component_names.size(), false)); + + unsigned int j = 0; + for (unsigned int i = 0; i < component_names.size(); ++i) + { + if (unique_component_names[j] != component_names[i]) + masks.push_back(ComponentMask(bools[j++])); + bools[j][i] = true; + } + masks.push_back(ComponentMask(bools[j++])); + AssertDimension(j, unique_component_names.size()); + return masks; + } +} // namespace + +ParsedConvergenceTable::ParsedConvergenceTable( + const std::vector &solution_names, + const std::vector> + &list_of_error_norms) + : ParsedConvergenceTable(solution_names, + list_of_error_norms, + {ParsedConvergenceTableFlags::cells, + ParsedConvergenceTableFlags::dofs}, + ParsedConvergenceTableFlags::dofs, + "", + 3, + true, + true) +{} + + + +ParsedConvergenceTable::ParsedConvergenceTable( + const std::vector &component_names, + const std::vector> + &list_of_error_norms, + const std::set &extra_columns, + const ParsedConvergenceTableFlags::ExtraColumns & rate_key, + const std::string & error_file_name, + const unsigned int & precision, + const bool & compute_error, + const bool & output_error) + : component_names(component_names) + , unique_component_names(get_unique_component_names(component_names)) + , unique_component_masks(get_unique_component_masks(component_names)) + , norms_per_unique_component(list_of_error_norms) + , extra_columns(extra_columns) + , rate_key(rate_key) + , precision(precision) + , error_file_name(error_file_name) + , compute_error(compute_error) + , output_error(output_error) +{ + AssertDimension(unique_component_names.size(), list_of_error_norms.size()); +} + + + +void +ParsedConvergenceTable::add_parameters(ParameterHandler &prm) +{ + prm.add_parameter("Enable output to streams", + output_error, + "When set to false, printing of the convergence table " + "to the stream specified as input to output_table() " + "is disabled."); + + prm.add_parameter("Enable computation of the errors", + compute_error, + "When set to false, no computations are performed."); + prm.add_parameter("Error precision", + precision, + "Number of digits to use when printing the error.", + Patterns::Integer(0)); + + prm.add_parameter("Error file name", + error_file_name, + "Set this to a filename with extension .txt, .gpl, .org, " + "or .tex to enable writing the convergence table to a " + "file."); + + prm.add_parameter( + "List of error norms to compute", + norms_per_unique_component, + "Each component is separated by a semicolon, " + "and each norm by a comma. Implemented norms are Linfty, L2, " + "W1infty, H1, and custom. If you want to skip a component, use none."); + + prm.add_parameter("Extra columns", + extra_columns, + "Extra columns to add to the table. Availabl options " + "are dofs, cells, and dt."); + + prm.add_parameter("Rate key", + rate_key, + "Key to use when computing convergence rates. If " + "this is set to a column that is not present, or to none, " + "then no error rates are computed."); +} + + + +void +ParsedConvergenceTable::output_table(std::ostream &out) +{ + if (compute_error) + { + // Add convergence rates if the rate_key is not empty + if (rate_key != ParsedConvergenceTableFlags::no_columns) + { + bool has_key = false; + for (const auto &col : extra_columns) + { + if (rate_key == col) + has_key = true; + + if (col != ParsedConvergenceTableFlags::no_columns) + table.omit_column_from_convergence_rate_evaluation( + Patterns::Tools::to_string(col)); + } + if (has_key) + table.evaluate_all_convergence_rates( + Patterns::Tools::to_string(rate_key), + ConvergenceTable::reduction_rate_log2); + } + + if (output_error) + table.write_text(out); + + if (error_file_name != "") + { + const std::string error_file_format = + error_file_name.substr(error_file_name.find_last_of(".") + 1); + + std::ofstream table_file(error_file_name); + + if (error_file_format == "tex") + table.write_tex(table_file); + else if (error_file_format == "txt") + table.write_text(table_file); + else if (error_file_format == "gpl") + table.write_text( + table_file, TableHandler::table_with_separate_column_description); + else if (error_file_format == "org") + table.write_text(table_file, TableHandler::org_mode_table); + else + { + AssertThrow(false, + ExcInternalError( + std::string("Unrecognized file format: ") + + error_file_format)); + } + table_file.close(); + } + } +} +DEAL_II_NAMESPACE_CLOSE -- 2.39.5