https://mooseframework.inl.gov
SteffensenSolve.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://mooseframework.inl.gov
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
10 #include "SteffensenSolve.h"
11 
12 #include "Executioner.h"
13 #include "FEProblemBase.h"
14 #include "NonlinearSystem.h"
16 #include "Console.h"
18 
21 {
23 
24  return params;
25 }
26 
28 
29 void
31 {
32  TagID fxn_m1_tagid;
33  TagID xn_m1_tagid;
34  const std::vector<PostprocessorName> * transformed_pps;
35  std::vector<std::vector<PostprocessorValue>> * transformed_pps_values;
36  if (primary)
37  {
38  xn_m1_tagid = _problem.addVectorTag("xn_m1", Moose::VECTOR_TAG_SOLUTION);
39  fxn_m1_tagid = _problem.addVectorTag("fxn_m1", Moose::VECTOR_TAG_SOLUTION);
40  transformed_pps = &_transformed_pps;
41  transformed_pps_values = &_transformed_pps_values;
42  _xn_m1_tagid = xn_m1_tagid;
43  _fxn_m1_tagid = fxn_m1_tagid;
44  }
45  else
46  {
47  xn_m1_tagid = _problem.addVectorTag("secondary_xn_m1", Moose::VECTOR_TAG_SOLUTION);
48  fxn_m1_tagid = _problem.addVectorTag("secondary_fxn_m1", Moose::VECTOR_TAG_SOLUTION);
49  transformed_pps = &_secondary_transformed_pps;
50  transformed_pps_values = &_secondary_transformed_pps_values;
51  _secondary_xn_m1_tagid = xn_m1_tagid;
52  _secondary_fxn_m1_tagid = fxn_m1_tagid;
53  }
54 
55  // Store a copy of the previous solution here
56  _solver_sys.addVector(xn_m1_tagid, false, PARALLEL);
57  _solver_sys.addVector(fxn_m1_tagid, false, PARALLEL);
58 
59  // Allocate storage for the previous postprocessor values
60  (*transformed_pps_values).resize((*transformed_pps).size());
61  for (size_t i = 0; i < (*transformed_pps).size(); i++)
62  (*transformed_pps_values)[i].resize(2);
63 }
64 
65 void
67 {
69 
71  if (!dynamic_cast<DefaultMultiAppFixedPointConvergence *>(&convergence))
72  mooseError(
73  "Only DefaultMultiAppFixedPointConvergence objects may be used for "
74  "'multiapp_fixed_point_convergence' when using the Steffensen fixed point algorithm.");
75 }
76 
77 void
79 {
80  unsigned int iteration;
81  TagID fxn_m1_tagid;
82  TagID xn_m1_tagid;
83  if (primary)
84  {
85  iteration = _fixed_point_it;
86  fxn_m1_tagid = _fxn_m1_tagid;
87  xn_m1_tagid = _xn_m1_tagid;
88  }
89  else
90  {
91  iteration = _main_fixed_point_it;
92  fxn_m1_tagid = _secondary_fxn_m1_tagid;
93  xn_m1_tagid = _secondary_xn_m1_tagid;
94  }
95 
96  // Check to make sure allocateStorage has been called
97  mooseAssert(fxn_m1_tagid != Moose::INVALID_TAG_ID,
98  "allocateStorage has not been called with primary = " + Moose::stringify(primary));
99  mooseAssert(xn_m1_tagid != Moose::INVALID_TAG_ID,
100  "allocateStorage has not been called with primary = " + Moose::stringify(primary));
101 
102  // Save previous variable values
104  NumericVector<Number> & fxn_m1 = _solver_sys.getVector(fxn_m1_tagid);
105  NumericVector<Number> & xn_m1 = _solver_sys.getVector(xn_m1_tagid);
106 
107  // What 'solution' is with regards to the Steffensen solve depends on the step
108  if (iteration % 2 == 1)
109  xn_m1 = solution;
110  else
111  fxn_m1 = solution;
112 }
113 
114 void
116 {
117  unsigned int iteration;
118  const std::vector<PostprocessorName> * transformed_pps;
119  std::vector<std::vector<PostprocessorValue>> * transformed_pps_values;
120  if (primary)
121  {
122  iteration = _fixed_point_it;
123  transformed_pps = &_transformed_pps;
124  transformed_pps_values = &_transformed_pps_values;
125  }
126  else
127  {
128  iteration = _main_fixed_point_it;
129  transformed_pps = &_secondary_transformed_pps;
130  transformed_pps_values = &_secondary_transformed_pps_values;
131  }
132 
133  // Save previous postprocessor values
134  for (size_t i = 0; i < (*transformed_pps).size(); i++)
135  {
136  if (iteration % 2 == 0)
137  (*transformed_pps_values)[i][1] = getPostprocessorValueByName((*transformed_pps)[i]);
138  else
139  (*transformed_pps_values)[i][0] = getPostprocessorValueByName((*transformed_pps)[i]);
140  }
141 }
142 
143 bool
145 {
146  // Need at least two values to compute the Steffensen update, and the update is only performed
147  // every other iteration as two evaluations of the coupled problem are necessary
148  if (primary)
149  return _fixed_point_it > 1 && (_fixed_point_it % 2 == 0);
150  else
151  return _main_fixed_point_it > 1 && (_main_fixed_point_it % 2 == 0);
152 }
153 
154 void
156 {
157  Real relaxation_factor;
158  const std::vector<PostprocessorName> * transformed_pps;
159  std::vector<std::vector<PostprocessorValue>> * transformed_pps_values;
160  if (primary)
161  {
162  relaxation_factor = _relax_factor;
163  transformed_pps = &_transformed_pps;
164  transformed_pps_values = &_transformed_pps_values;
165  }
166  else
167  {
168  relaxation_factor = _secondary_relaxation_factor;
169  transformed_pps = &_secondary_transformed_pps;
170  transformed_pps_values = &_secondary_transformed_pps_values;
171  }
172 
173  // Relax postprocessors for the main application
174  for (size_t i = 0; i < (*transformed_pps).size(); i++)
175  {
176  // Get new postprocessor value
177  const Real current_value = getPostprocessorValueByName((*transformed_pps)[i]);
178  const Real fxn_m1 = (*transformed_pps_values)[i][0];
179  const Real xn_m1 = (*transformed_pps_values)[i][1];
180 
181  // Compute and set relaxed value
182  Real new_value = current_value;
183  if (!MooseUtils::absoluteFuzzyEqual(current_value + xn_m1 - 2 * fxn_m1, 0))
184  new_value =
185  xn_m1 - (fxn_m1 - xn_m1) * (fxn_m1 - xn_m1) / (current_value + xn_m1 - 2 * fxn_m1);
186 
187  // Relax update
188  new_value = relaxation_factor * new_value + (1 - relaxation_factor) * fxn_m1;
189 
190  _problem.setPostprocessorValueByName((*transformed_pps)[i], new_value);
191  }
192 }
193 
194 void
195 SteffensenSolve::transformVariables(const std::set<dof_id_type> & transformed_dofs,
196  const bool primary)
197 {
198  Real relaxation_factor;
199  TagID fxn_m1_tagid;
200  TagID xn_m1_tagid;
201  if (primary)
202  {
203  relaxation_factor = _relax_factor;
204  fxn_m1_tagid = _fxn_m1_tagid;
205  xn_m1_tagid = _xn_m1_tagid;
206  }
207  else
208  {
209  relaxation_factor = _secondary_relaxation_factor;
210  fxn_m1_tagid = _secondary_fxn_m1_tagid;
211  xn_m1_tagid = _secondary_xn_m1_tagid;
212  }
213 
215  NumericVector<Number> & fxn_m1 = _solver_sys.getVector(fxn_m1_tagid);
216  NumericVector<Number> & xn_m1 = _solver_sys.getVector(xn_m1_tagid);
217 
218  for (const auto & dof : transformed_dofs)
219  {
220  // Avoid 0 denominator issue
221  Real new_value = solution(dof);
222  if (!MooseUtils::absoluteFuzzyEqual(solution(dof) + xn_m1(dof) - 2 * fxn_m1(dof), 0))
223  new_value = xn_m1(dof) - (fxn_m1(dof) - xn_m1(dof)) * (fxn_m1(dof) - xn_m1(dof)) /
224  (solution(dof) + xn_m1(dof) - 2 * fxn_m1(dof));
225 
226  // Relax update
227  new_value = relaxation_factor * new_value + (1 - relaxation_factor) * fxn_m1(dof);
228 
229  solution.set(dof, new_value);
230  }
231  solution.close();
233 }
234 
235 void
237  const Real initial_norm,
238  const std::vector<Real> & timestep_begin_norms,
239  const std::vector<Real> & timestep_end_norms) const
240 {
241  _console << "\n 0 Steffensen initialization |R| = "
242  << Console::outputNorm(std::numeric_limits<Real>::max(), initial_norm) << '\n';
243 
244  Real max_norm_old = initial_norm;
245  for (unsigned int i = 0; i <= _fixed_point_it; ++i)
246  {
247  Real max_norm = std::max(timestep_begin_norms[i], timestep_end_norms[i]);
248  std::stringstream steffensen_prefix;
249  if (i == 0)
250  steffensen_prefix << " Steffensen initialization |R| = ";
251  else if (i % 2 == 0)
252  steffensen_prefix << " Steffensen half-step |R| = ";
253  else
254  steffensen_prefix << " Steffensen step |R| = ";
255 
256  _console << std::setw(2) << i + 1 << steffensen_prefix.str()
257  << Console::outputNorm(max_norm_old, max_norm) << '\n';
258 
259  max_norm_old = max_norm;
260  }
261 }
std::vector< std::vector< PostprocessorValue > > _transformed_pps_values
Previous values of the relaxed postprocessors.
virtual void printFixedPointConvergenceHistory(Real initial_norm, const std::vector< Real > &timestep_begin_norms, const std::vector< Real > &timestep_end_norms) const override final
Print the convergence history of the coupling, at every fixed point iteration.
std::vector< std::vector< PostprocessorValue > > _secondary_transformed_pps_values
Previous values of the postprocessors relaxed outside of the fixed point iteration (used as a subapp)...
FEProblemBase & _problem
Reference to FEProblem.
Definition: SolveObject.h:47
SteffensenSolve(Executioner &ex)
bool absoluteFuzzyEqual(const T &var1, const T2 &var2, const T3 &tol=libMesh::TOLERANCE *libMesh::TOLERANCE)
Function to check whether two variables are equal within an absolute tolerance.
Definition: MooseUtils.h:380
unsigned int TagID
Definition: MooseTypes.h:210
NumericVector< Number > & solution()
Definition: SystemBase.h:195
PARALLEL
void setPostprocessorValueByName(const PostprocessorName &name, const PostprocessorValue &value, std::size_t t_index=0)
Set the value of a PostprocessorValue.
virtual TagID addVectorTag(const TagName &tag_name, const Moose::VectorTagType type=Moose::VECTOR_TAG_RESIDUAL)
Create a Tag.
Definition: SubProblem.C:92
virtual void transformPostprocessors(const bool primary) override final
Use the fixed point algorithm to transform the postprocessors.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
virtual void savePostprocessorValues(const bool primary) override final
Saves the current values of the postprocessors, and update the old(er) vectors.
virtual void transformVariables(const std::set< dof_id_type > &transformed_dofs, const bool primary) override final
Use the fixed point algorithm to transform the variables.
std::vector< PostprocessorName > _secondary_transformed_pps
Postprocessors to be relaxed outside of fixed point iteration (used as a subapp)
auto max(const L &left, const R &right)
const std::vector< PostprocessorName > _transformed_pps
The postprocessors (transferred or not) that are going to be relaxed.
NumericVector< Number > & addVector(const std::string &vector_name, const bool project, const libMesh::ParallelType type)
Adds a solution length vector to the system.
TagID _fxn_m1_tagid
Vector tag id for the most recent solution variable, pre-Steffensen transform, as a main app...
TagID _xn_m1_tagid
Vector tag id for the solution variable before the latest solve, as a main app.
void update()
Update the system (doing libMesh magic)
Definition: SystemBase.C:1245
const Real _relax_factor
Relaxation factor for fixed point Iteration.
virtual Convergence & getConvergence(const std::string &name, const THREAD_ID tid=0) const
Gets a Convergence object.
virtual bool useFixedPointAlgorithmUpdateInsteadOfPicard(const bool primary) override final
Use the fixed point algorithm transform instead of simply using the Picard update.
const ConvergenceName & getMultiAppFixedPointConvergenceName() const
Gets the MultiApp fixed point convergence object name.
TagID _secondary_xn_m1_tagid
Vector tag id for the solution variable before the latest solve, as a sub app.
Executioners are objects that do the actual work of solving your problem.
Definition: Executioner.h:30
virtual void saveVariableValues(const bool primary) override final
Saves the current values of the variables, and update the old(er) vectors.
Real _secondary_relaxation_factor
Relaxation factor outside of fixed point iteration (used as a subapp)
virtual void allocateStorage(const bool primary) override final
Allocate storage for the fixed point algorithm.
std::string stringify(const T &t)
conversion to string
Definition: Conversion.h:64
virtual void close()=0
virtual const PostprocessorValue & getPostprocessorValueByName(const PostprocessorName &name) const
Retrieve the value of the Postprocessor.
unsigned int _fixed_point_it
SystemBase & _solver_sys
Reference to a system for creating vectors as needed for the solve, etc.
Definition: SolveObject.h:55
static InputParameters validParams()
TagID _secondary_fxn_m1_tagid
Vector tag id for the most recent solution variable, pre-Steffensen transform, as a sub app...
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual void initialSetup() override
Method that should be executed once, before any solve calls.
static std::string outputNorm(const Real &old_norm, const Real &norm, const unsigned int precision=6)
A helper function for outputting norms in color.
Definition: Console.C:623
const TagID INVALID_TAG_ID
Definition: MooseTypes.C:23
void mooseError(Args &&... args) const
Emits an error prefixed with object name and type.
virtual void set(const numeric_index_type i, const Number value)=0
const ConsoleStream _console
An instance of helper class to write streams to the Console objects.
virtual void initialSetup() override
Method that should be executed once, before any solve calls.
virtual NumericVector< Number > & getVector(const std::string &name)
Get a raw NumericVector by name.
Definition: SystemBase.C:916
static InputParameters validParams()
unsigned int _main_fixed_point_it
fixed point iteration counter for the main app