Loading [MathJax]/extensions/TeX/newcommand.js
 deal.II version GIT relicensing-3104-g8c3bc2e695 2025-04-21 18:40:00+00:00
\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}} \newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=} \newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]} \newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages Concepts
parsed_convergence_table.cc
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2019 - 2023 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
19
20#include <fstream>
21#include <set>
22
23
25
26namespace
27{
28 std::vector<std::string>
29 get_unique_component_names(const std::vector<std::string> &component_names)
30 {
31 auto elements = component_names;
32 elements.erase(std::unique(elements.begin(), elements.end()),
33 elements.end());
34 return elements;
35 }
36
37
38
39 std::vector<ComponentMask>
40 get_unique_component_masks(const std::vector<std::string> &component_names)
41 {
42 const auto unique_component_names =
43 get_unique_component_names(component_names);
44
45 std::vector<ComponentMask> masks;
46
47 std::vector<std::vector<bool>> bools(
48 unique_component_names.size(),
49 std::vector<bool>(component_names.size(), false));
50
51 unsigned int j = 0;
52 for (unsigned int i = 0; i < component_names.size(); ++i)
53 {
54 if (unique_component_names[j] != component_names[i])
55 masks.emplace_back(bools[j++]);
56 bools[j][i] = true;
57 }
58 masks.emplace_back(bools[j++]);
59 AssertDimension(j, unique_component_names.size());
60 return masks;
61 }
62} // namespace
63
64
65
67 const std::vector<std::string> &solution_names,
68 const std::vector<std::set<VectorTools::NormType>> &list_of_error_norms)
69 : ParsedConvergenceTable(solution_names,
70 list_of_error_norms,
71 2.0,
72 {"cells", "dofs"},
73 "dofs",
74 "reduction_rate_log2",
75 "",
76 3,
77 true)
78{}
79
80
81
83 const std::vector<std::string> &component_names,
84 const std::vector<std::set<VectorTools::NormType>> &list_of_error_norms,
85 const double exponent,
86 const std::set<std::string> &extra_columns,
87 const std::string &rate_key,
88 const std::string &rate_mode,
89 const std::string &error_file_name,
90 const unsigned int precision,
91 const bool compute_error)
92 : component_names(component_names)
93 , unique_component_names(get_unique_component_names(component_names))
94 , unique_component_masks(get_unique_component_masks(component_names))
95 , norms_per_unique_component(list_of_error_norms)
96 , exponent(exponent)
97 , extra_columns(extra_columns)
98 , rate_key(rate_key)
99 , rate_mode(rate_mode)
100 , precision(precision)
101 , error_file_name(error_file_name)
102 , compute_error(compute_error)
103{
106}
107
108
109
110void
112{
113 prm.add_parameter("Enable computation of the errors",
115 "When set to false, no computations are performed.");
116
117 prm.add_parameter("Error precision",
118 precision,
119 "Number of digits to use when printing the error.",
121
122 prm.add_parameter("Error file name",
124 "Set this to a filename with extension .txt, .gpl, .org, "
125 "or .tex to enable writing the convergence table to a "
126 "file.");
127
128 prm.add_parameter(
129 "List of error norms to compute",
131 "Each component is separated by a semicolon "
132 "and each norm by a comma. See the documentation of VectorTools::NormType "
133 "for a list of implemented norms. If you want to skip a component, leave "
134 "its entry empty.");
135
136
137 prm.add_parameter("Exponent for p-norms",
138 exponent,
139 "The exponent to use when computing p-norms.",
141
142 prm.add_parameter("Extra columns",
144 "Extra columns to add to the table. Available options "
145 "are dofs and cells.",
146 Patterns::List(Patterns::Selection("dofs|cells")));
147
148 prm.add_parameter("Rate key",
149 rate_key,
150 "Key to use when computing convergence rates. If "
151 "this is set to a column that is not present, or to the "
152 "empty string, then no error rates are computed.");
153
154 prm.add_parameter("Rate mode",
155 rate_mode,
156 "What type of error rate to compute. Available options are "
157 "reduction_rate_log2, reduction_rate, and none.",
159 "reduction_rate|reduction_rate_log2|none"));
160}
161
162
163
164void
166{
167 if (compute_error)
168 {
169 // Add convergence rates if the rate_key is not empty
170 if (rate_key != "")
171 {
172 bool has_key = false;
173 for (const auto &col : extra_columns)
174 {
175 if (rate_key == col)
176 has_key = true;
177
178 if (col != "")
180 }
181
182 for (const auto &extra_col : extra_column_functions)
183 if (extra_col.second.second == false)
184 {
185 if (rate_key == extra_col.first)
186 has_key = true;
188 extra_col.first);
189 }
190
191 if (has_key)
192 {
193 if (rate_mode == "reduction_rate_log2")
196 else if (rate_mode == "reduction_rate")
199 else
200 {
201 Assert(rate_mode != "none", ExcInternalError());
202 }
203 }
204 else
205 {
206 AssertThrow(rate_key != "",
208 "You specified the key <" + rate_key +
209 "> to compute convergence rates, but that key does "
210 "not exist in the current table."));
211 }
212 }
213 }
214}
215
216
217
218void
220{
221 if (compute_error)
222 {
224 table.write_text(out);
225 output_table();
226 }
227}
228
229
230
231void
233{
234 if (compute_error && error_file_name != "")
235 {
237
238 const std::string error_file_format =
239 error_file_name.substr(error_file_name.find_last_of('.') + 1);
240
241 std::ofstream table_file(error_file_name);
242
243 if (error_file_format == "tex")
244 table.write_tex(table_file);
245 else if (error_file_format == "txt")
246 table.write_text(table_file);
247 else if (error_file_format == "gpl")
248 table.write_text(table_file,
250 else if (error_file_format == "org")
252 else
253 {
254 AssertThrow(false,
255 ExcInternalError("Unrecognized file format: " +
256 error_file_format));
257 }
258 table_file.close();
259 }
260}
261
262
263
264void
266 const std::string &column_name,
267 const std::function<double()> &custom_function,
268 const bool compute_rate)
269{
270 extra_column_functions[column_name] = {custom_function, compute_rate};
271}
272
void evaluate_all_convergence_rates(const std::string &reference_column_key, const RateMode rate_mode)
void omit_column_from_convergence_rate_evaluation(const std::string &key)
void add_parameter(const std::string &entry, ParameterType &parameter, const std::string &documentation="", const Patterns::PatternBase &pattern= *Patterns::Tools::Convert< ParameterType >::to_pattern(), const bool has_to_be_set=false)
The ParsedConvergenceTable class.
void add_extra_column(const std::string &column_name, const std::function< double()> &custom_function, const bool compute_rate=true)
void add_parameters(ParameterHandler &prm)
std::set< std::string > extra_columns
ParsedConvergenceTable(const std::vector< std::string > &component_names={"u"}, const std::vector< std::set< VectorTools::NormType > > &list_of_error_norms={ {VectorTools::H1_norm, VectorTools::L2_norm, VectorTools::Linfty_norm}})
const std::vector< std::string > unique_component_names
std::vector< std::set< VectorTools::NormType > > norms_per_unique_component
std::map< std::string, std::pair< std::function< double()>, bool > > extra_column_functions
void write_text(std::ostream &out, const TextOutputFormat format=table_with_headers) const
@ table_with_separate_column_description
void write_tex(std::ostream &file, const bool with_header=true) const
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:35
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:36
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)