Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
convergence_table.cc
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 #include <deal.II/base/convergence_table.h>
17 
18 #include <cmath>
19 
20 DEAL_II_NAMESPACE_OPEN
21 
22 void
24  const std::string &data_column_key,
25  const std::string &reference_column_key,
26  const RateMode rate_mode,
27  const unsigned int dim)
28 {
29  Assert(columns.count(data_column_key), ExcColumnNotExistent(data_column_key));
30  Assert(columns.count(reference_column_key),
31  ExcColumnNotExistent(reference_column_key));
32 
33  if (rate_mode == none)
34  return;
35 
36  // reset the auto fill mode flag since we are going to fill columns from
37  // the top that don't yet exist
38  set_auto_fill_mode(false);
39 
40  std::vector<internal::TableEntry> &entries = columns[data_column_key].entries;
41  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])) /
97  std::log(
98  std::fabs(ref_values[i] / ref_values[i - 1])));
99  }
100  }
101  break;
102  default:
103  AssertThrow(false, ExcNotImplemented());
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))
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 
126 void
127 ConvergenceTable::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  std::vector<internal::TableEntry> &entries = columns[data_column_key].entries;
137  std::string rate_key = data_column_key + "...";
138 
139  const unsigned int n = entries.size();
140 
141  std::vector<double> values(n);
142  for (unsigned int i = 0; i < n; ++i)
143  values[i] = entries[i].get_numeric_value();
144 
145  unsigned int no_rate_entries = 0;
146 
147  switch (rate_mode)
148  {
149  case none:
150  break;
151 
152  case reduction_rate:
153  rate_key += "red.rate";
154  no_rate_entries = columns[rate_key].entries.size();
155  // Calculate all missing rate values:
156  for (unsigned int i = no_rate_entries; i < n; ++i)
157  {
158  if (i == 0)
159  {
160  // no value available for the first row
161  add_value(rate_key, std::string("-"));
162  }
163  else
164  {
165  add_value(rate_key, values[i - 1] / values[i]);
166  }
167  }
168  break;
169 
170  case reduction_rate_log2:
171  rate_key += "red.rate.log2";
172  no_rate_entries = columns[rate_key].entries.size();
173  // Calculate all missing rate values:
174  for (unsigned int i = no_rate_entries; i < n; ++i)
175  {
176  if (i == 0)
177  {
178  // no value available for the first row
179  add_value(rate_key, std::string("-"));
180  }
181  else
182  {
183  add_value(rate_key,
184  std::log(std::fabs(values[i - 1] / values[i])) /
185  std::log(2.0));
186  }
187  }
188  break;
189 
190  default:
191  AssertThrow(false, ExcNotImplemented());
192  }
193 
194  Assert(columns.count(rate_key), ExcInternalError());
195  columns[rate_key].flag = 1;
196  set_precision(rate_key, 2);
197 
198  // set the superkey equal to the key
199  const std::string &superkey = data_column_key;
200  // and set the tex caption of the supercolumn to the tex caption of the
201  // data_column.
202  if (!supercolumns.count(superkey))
203  {
204  add_column_to_supercolumn(data_column_key, superkey);
205  set_tex_supercaption(superkey, columns[data_column_key].tex_caption);
206  }
207 
208  // only add rate_key to the supercolumn once
209  if (no_rate_entries == 0)
210  {
211  add_column_to_supercolumn(rate_key, superkey);
212  }
213 }
214 
215 
216 
217 void
219  const std::string &key)
220 {
221  Assert(columns.count(key), ExcColumnNotExistent(key));
222 
223  const std::map<std::string, Column>::iterator col_iter = columns.find(key);
224  col_iter->second.flag = 1;
225 }
226 
227 
228 
229 void
231  const std::string &reference_column_key,
232  const RateMode rate_mode)
233 {
234  for (std::map<std::string, Column>::const_iterator col_iter = columns.begin();
235  col_iter != columns.end();
236  ++col_iter)
237  if (!col_iter->second.flag)
238  evaluate_convergence_rates(col_iter->first,
239  reference_column_key,
240  rate_mode);
241 }
242 
243 
244 
245 void
247 {
248  for (std::map<std::string, Column>::const_iterator col_iter = columns.begin();
249  col_iter != columns.end();
250  ++col_iter)
251  if (!col_iter->second.flag)
252  evaluate_convergence_rates(col_iter->first, rate_mode);
253 }
254 
255 DEAL_II_NAMESPACE_CLOSE
void set_precision(const std::string &key, const unsigned int precision)
void set_tex_supercaption(const std::string &superkey, const std::string &tex_supercaption)
std::map< std::string, Column > columns
void add_value(const std::string &key, const T value)
#define AssertThrow(cond, exc)
Definition: exceptions.h:1519
void add_column_to_supercolumn(const std::string &key, const std::string &superkey)
void omit_column_from_convergence_rate_evaluation(const std::string &key)
void evaluate_all_convergence_rates(const std::string &reference_column_key, const RateMode rate_mode)
#define Assert(cond, exc)
Definition: exceptions.h:1407
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
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)
static ::ExceptionBase & ExcColumnNotExistent(std::string arg1)
static ::ExceptionBase & ExcNotImplemented()
std::map< std::string, std::vector< std::string > > supercolumns
static ::ExceptionBase & ExcInternalError()
void set_auto_fill_mode(const bool state)