OGS
ConvergenceCriterionPerComponentResidual.cpp
Go to the documentation of this file.
1
12
13#include "BaseLib/ConfigTree.h"
14#include "BaseLib/Logging.h"
17
18namespace NumLib
19{
22 std::vector<double>&& absolute_tolerances,
23 std::vector<double>&& relative_tolerances,
24 const MathLib::VecNormType norm_type)
26 _abstols(std::move(absolute_tolerances)),
27 _reltols(std::move(relative_tolerances)),
28 _residual_norms_0(_abstols.size())
29{
30 if (_abstols.size() != _reltols.size())
31 {
33 "The number of absolute and relative tolerances given must be the "
34 "same.");
35 }
36
37 if (_abstols.empty())
38 {
39 OGS_FATAL("The given tolerances vector is empty.");
40 }
41}
42
44 const GlobalVector& minus_delta_x, GlobalVector const& x)
45{
46 if ((!_dof_table) || (!_mesh))
47 {
48 OGS_FATAL("D.o.f. table or mesh have not been set.");
49 }
50
51 for (unsigned global_component = 0; global_component < _abstols.size();
52 ++global_component)
53 {
54 // TODO short cut if tol <= 0.0
55 auto error_dx = norm(minus_delta_x, global_component, _norm_type,
56 *_dof_table, *_mesh);
57 auto norm_x =
58 norm(x, global_component, _norm_type, *_dof_table, *_mesh);
59
60 INFO(
61 "Convergence criterion, component {:d}: |dx|={:.4e}, |x|={:.4e}, "
62 "|dx|/|x|={:.4e}",
63 global_component, error_dx, norm_x,
64 (norm_x == 0. ? std::numeric_limits<double>::quiet_NaN()
65 : (error_dx / norm_x)));
66 }
67}
68
70 const GlobalVector& residual)
71{
72 if ((!_dof_table) || (!_mesh))
73 {
74 OGS_FATAL("D.o.f. table or mesh have not been set.");
75 }
76
77 for (unsigned global_component = 0; global_component < _abstols.size();
78 ++global_component)
79 {
80 // TODO short cut if tol <= 0.0
81 auto norm_res =
82 norm(residual, global_component, _norm_type, *_dof_table, *_mesh);
83
85 {
86 INFO("Convergence criterion, component {:d}: |r0|={:.4e}",
87 global_component, norm_res);
88 _residual_norms_0[global_component] = norm_res;
89 }
90 else
91 {
92 auto const norm_res0 = _residual_norms_0[global_component];
93 INFO(
94 "Convergence criterion, component {:d}: |r|={:.4e}, "
95 "|r0|={:.4e}, |r|/|r0|={:.4e}",
96 global_component, norm_res, norm_res0,
97 (norm_res0 == 0. ? std::numeric_limits<double>::quiet_NaN()
98 : (norm_res / norm_res0)));
99 }
100
101 bool const satisfied_abs = norm_res < _abstols[global_component];
102 // Make sure that in the first iteration the relative residual tolerance
103 // is not satisfied.
104 bool const satisfied_rel =
106 checkRelativeTolerance(_reltols[global_component], norm_res,
107 _residual_norms_0[global_component]);
108 _satisfied = _satisfied && (satisfied_abs || satisfied_rel);
109 }
110}
111
113 const LocalToGlobalIndexMap& dof_table, MeshLib::Mesh const& mesh)
114{
115 _dof_table = &dof_table;
116 _mesh = &mesh;
117
119 static_cast<int>(_abstols.size()))
120 {
121 OGS_FATAL(
122 "The number of components in the DOF table and the number of "
123 "tolerances given do not match.");
124 }
125}
126
127std::unique_ptr<ConvergenceCriterionPerComponentResidual>
129 const BaseLib::ConfigTree& config)
130{
132 config.checkConfigParameter("type", "PerComponentResidual");
133
134 auto abstols =
136 config.getConfigParameterOptional<std::vector<double>>("abstols");
137 auto reltols =
139 config.getConfigParameterOptional<std::vector<double>>("reltols");
140 auto const norm_type_str =
142 config.getConfigParameter<std::string>("norm_type");
143
144 if ((!abstols) && (!reltols))
145 {
146 OGS_FATAL(
147 "At least one of absolute or relative tolerance has to be "
148 "specified.");
149 }
150 if (!abstols)
151 {
152 abstols = std::vector<double>(reltols->size());
153 }
154 else if (!reltols)
155 {
156 reltols = std::vector<double>(abstols->size());
157 }
158
159 auto const norm_type = MathLib::convertStringToVecNormType(norm_type_str);
160
161 if (norm_type == MathLib::VecNormType::INVALID)
162 {
163 OGS_FATAL("Unknown vector norm type `{:s}'.", norm_type_str);
164 }
165
166 return std::make_unique<ConvergenceCriterionPerComponentResidual>(
167 std::move(*abstols), std::move(*reltols), norm_type);
168}
169
170} // namespace NumLib
#define OGS_FATAL(...)
Definition Error.h:26
void INFO(fmt::format_string< Args... > fmt, Args &&... args)
Definition Logging.h:35
std::optional< T > getConfigParameterOptional(std::string const &param) const
T getConfigParameter(std::string const &param) const
void checkConfigParameter(std::string const &param, std::string_view const value) const
Global vector based on Eigen vector.
Definition EigenVector.h:25
void setDOFTable(const LocalToGlobalIndexMap &dof_table, MeshLib::Mesh const &mesh) override
Sets the d.o.f. table used to extract data for a specific component.
void checkResidual(const GlobalVector &residual) override
Check if the residual satisfies the convergence criterion.
void checkDeltaX(const GlobalVector &minus_delta_x, GlobalVector const &x) override
ConvergenceCriterionPerComponentResidual(std::vector< double > &&absolute_tolerances, std::vector< double > &&relative_tolerances, const MathLib::VecNormType norm_type)
const MathLib::VecNormType _norm_type
VecNormType
Norm type. Not declared as class type in order to use the members as integers.
Definition LinAlgEnums.h:22
VecNormType convertStringToVecNormType(const std::string &str)
convert string to VecNormType
double norm(GlobalVector const &x, unsigned const global_component, MathLib::VecNormType norm_type, LocalToGlobalIndexMap const &dof_table, MeshLib::Mesh const &mesh)
bool checkRelativeTolerance(const double reltol, const double numerator, const double denominator)
std::unique_ptr< ConvergenceCriterionPerComponentResidual > createConvergenceCriterionPerComponentResidual(const BaseLib::ConfigTree &config)