deal.II version GIT relicensing-1989-gd7a2c90e4e 2024-10-14 01:50: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\}}\)
Loading...
Searching...
No Matches
convergence_table.cc
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 1999 - 2024 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
16
17#include <cmath>
18
20
21void
23 const std::string &data_column_key,
24 const std::string &reference_column_key,
25 const RateMode rate_mode,
26 const unsigned int dim)
27{
28 Assert(columns.count(data_column_key), ExcColumnNotExistent(data_column_key));
29 Assert(columns.count(reference_column_key),
30 ExcColumnNotExistent(reference_column_key));
31
32 if (rate_mode == none)
33 return;
34
35 // reset the auto fill mode flag since we are going to fill columns from
36 // the top that don't yet exist
37 set_auto_fill_mode(false);
38
39 const std::vector<internal::TableEntry> &entries =
40 columns[data_column_key].entries;
41 const std::vector<internal::TableEntry> &ref_entries =
42 columns[reference_column_key].entries;
43 std::string rate_key = data_column_key + "...";
44
45 const unsigned int n = entries.size();
46 const unsigned int n_ref = ref_entries.size();
47 Assert(n == n_ref, ExcDimensionMismatch(n, n_ref));
48
49 std::vector<double> values(n);
50 std::vector<double> ref_values(n_ref);
51
52 for (unsigned int i = 0; i < n; ++i)
53 {
54 values[i] = entries[i].get_numeric_value();
55 ref_values[i] = ref_entries[i].get_numeric_value();
56 }
57
58 unsigned int no_rate_entries = 0;
59
60 switch (rate_mode)
61 {
62 // case none: already considered above
63 case reduction_rate:
64 rate_key += "red.rate";
65 no_rate_entries = columns[rate_key].entries.size();
66 // Calculate all missing rate values:
67 for (unsigned int i = no_rate_entries; i < n; ++i)
68 {
69 if (i == 0)
70 {
71 // no value available for the first row
72 add_value(rate_key, std::string("-"));
73 }
74 else
75 {
76 add_value(rate_key,
77 values[i - 1] / values[i] * ref_values[i] /
78 ref_values[i - 1]);
79 }
80 }
81 break;
83 rate_key += "red.rate.log2";
84 no_rate_entries = columns[rate_key].entries.size();
85 // Calculate all missing rate values:
86 for (unsigned int i = no_rate_entries; i < n; ++i)
87 {
88 if (i == 0)
89 {
90 // no value available for the first row
91 add_value(rate_key, std::string("-"));
92 }
93 else
94 {
95 add_value(rate_key,
96 dim * std::log(std::fabs(values[i - 1] / values[i])) /
98 std::fabs(ref_values[i] / ref_values[i - 1])));
99 }
100 }
101 break;
102 default:
104 }
105
106 Assert(columns.count(rate_key), ExcInternalError());
107 columns[rate_key].flag = 1;
108 set_precision(rate_key, 2);
109
110 const std::string &superkey = data_column_key;
111 if (supercolumns.count(superkey) == 0u)
112 {
113 add_column_to_supercolumn(data_column_key, superkey);
114 set_tex_supercaption(superkey, columns[data_column_key].tex_caption);
115 }
116
117 // only add rate_key to the supercolumn once
118 if (no_rate_entries == 0)
119 {
120 add_column_to_supercolumn(rate_key, superkey);
121 }
122}
123
124
125
126void
127ConvergenceTable::evaluate_convergence_rates(const std::string &data_column_key,
128 const RateMode rate_mode)
129{
130 Assert(columns.count(data_column_key), ExcColumnNotExistent(data_column_key));
131
132 // reset the auto fill mode flag since we are going to fill columns from
133 // the top that don't yet exist
134 set_auto_fill_mode(false);
135
136 const std::vector<internal::TableEntry> &entries =
137 columns[data_column_key].entries;
138 std::string rate_key = data_column_key + "...";
139
140 const unsigned int n = entries.size();
141
142 std::vector<double> values(n);
143 for (unsigned int i = 0; i < n; ++i)
144 values[i] = entries[i].get_numeric_value();
145
146 unsigned int no_rate_entries = 0;
147
148 switch (rate_mode)
149 {
150 case none:
151 break;
152
153 case reduction_rate:
154 rate_key += "red.rate";
155 no_rate_entries = columns[rate_key].entries.size();
156 // Calculate all missing rate values:
157 for (unsigned int i = no_rate_entries; i < n; ++i)
158 {
159 if (i == 0)
160 {
161 // no value available for the first row
162 add_value(rate_key, std::string("-"));
163 }
164 else
165 {
166 add_value(rate_key, values[i - 1] / values[i]);
167 }
168 }
169 break;
170
172 rate_key += "red.rate.log2";
173 no_rate_entries = columns[rate_key].entries.size();
174 // Calculate all missing rate values:
175 for (unsigned int i = no_rate_entries; i < n; ++i)
176 {
177 if (i == 0)
178 {
179 // no value available for the first row
180 add_value(rate_key, std::string("-"));
181 }
182 else
183 {
184 add_value(rate_key,
185 std::log(std::fabs(values[i - 1] / values[i])) /
186 std::log(2.0));
187 }
188 }
189 break;
190
191 default:
193 }
194
195 Assert(columns.count(rate_key), ExcInternalError());
196 columns[rate_key].flag = 1;
197 set_precision(rate_key, 2);
198
199 // set the superkey equal to the key
200 const std::string &superkey = data_column_key;
201 // and set the tex caption of the supercolumn to the tex caption of the
202 // data_column.
203 if (supercolumns.count(superkey) == 0u)
204 {
205 add_column_to_supercolumn(data_column_key, superkey);
206 set_tex_supercaption(superkey, columns[data_column_key].tex_caption);
207 }
208
209 // only add rate_key to the supercolumn once
210 if (no_rate_entries == 0)
211 {
212 add_column_to_supercolumn(rate_key, superkey);
213 }
214}
215
216
217
218void
220 const std::string &key)
221{
222 Assert(columns.count(key), ExcColumnNotExistent(key));
223
224 const std::map<std::string, Column>::iterator col_iter = columns.find(key);
225 col_iter->second.flag = 1;
226}
227
228
229
230void
232 const std::string &reference_column_key,
233 const RateMode rate_mode)
234{
235 for (std::map<std::string, Column>::const_iterator col_iter = columns.begin();
236 col_iter != columns.end();
237 ++col_iter)
238 if (col_iter->second.flag == 0u)
239 evaluate_convergence_rates(col_iter->first,
240 reference_column_key,
241 rate_mode);
242}
243
244
245
246void
248{
249 for (std::map<std::string, Column>::const_iterator col_iter = columns.begin();
250 col_iter != columns.end();
251 ++col_iter)
252 if (col_iter->second.flag == 0u)
253 evaluate_convergence_rates(col_iter->first, rate_mode);
254}
255
void evaluate_all_convergence_rates(const std::string &reference_column_key, const RateMode rate_mode)
void evaluate_convergence_rates(const std::string &data_column_key, const std::string &reference_column_key, const RateMode rate_mode, const unsigned int dim=2)
void omit_column_from_convergence_rate_evaluation(const std::string &key)
void set_tex_supercaption(const std::string &superkey, const std::string &tex_supercaption)
void set_auto_fill_mode(const bool state)
void add_value(const std::string &key, const T value)
std::map< std::string, std::vector< std::string > > supercolumns
void add_column_to_supercolumn(const std::string &key, const std::string &superkey)
void set_precision(const std::string &key, const unsigned int precision)
std::map< std::string, Column > columns
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:498
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:499
static ::ExceptionBase & ExcColumnNotExistent(std::string arg1)
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define AssertThrow(cond, exc)
::VectorizedArray< Number, width > log(const ::VectorizedArray< Number, width > &)