1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2020 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
17 
18 
19 
20 #ifndef LIBMESH_DIFF_SOLVER_H
21 #define LIBMESH_DIFF_SOLVER_H
22 
23 // Local includes
24 #include "libmesh/libmesh_common.h"
25 #include "libmesh/reference_counted_object.h"
26 #include "libmesh/parallel_object.h"
27 
28 // C++ includes
29 #include <vector>
30 #include <memory>
31 
32 namespace libMesh
33 {
34 
35 // Forward Declarations
36 class ImplicitSystem;
37 template <typename T> class NumericVector;
38 
39 /**
40  * Functor for use as callback in solve of nonlinear solver.
41  */
42 class LinearSolutionMonitor
43 {
44 public:
45   virtual void operator() (const NumericVector<Number> & delta_u,
46                            const Real & norm_delta_u,
47                            const NumericVector<Number> & u,
48                            const Real & norm_u,
49                            const NumericVector<Number> & res,
50                            const Real & norm_res,
51                            const unsigned int iteration) = 0;
52   virtual ~LinearSolutionMonitor();
53 };
54 
~LinearSolutionMonitor()55 inline LinearSolutionMonitor::~LinearSolutionMonitor() {}
56 
57 /**
58  * This is a generic class that defines a solver to handle
59  * ImplicitSystem classes, including NonlinearImplicitSystem and
60  * DifferentiableSystem   A user can define a solver by
61  * deriving from this class and implementing certain functions.
62  *
63  * This class is part of the new DifferentiableSystem framework,
64  * which is still experimental.  Users of this framework should
65  * beware of bugs and future API changes.
66  *
67  * \author Roy H. Stogner
68  * \date 2006-2010
69  */
70 class DiffSolver : public ReferenceCountedObject<DiffSolver>,
71                    public ParallelObject
72 {
73 public:
74   /**
75    * The type of system
76    */
77   typedef ImplicitSystem sys_type;
78 
79   /**
80    * Constructor. Requires a reference to the system
81    * to be solved.
82    */
83   DiffSolver (sys_type & s);
84 
85   /**
86    * Factory method.  Requires a reference to the system
87    * to be solved.
88    *
89    * \returns A NewtonSolver by default.
90    */
91   static std::unique_ptr<DiffSolver> build(sys_type & s);
92 
93   /**
94    * Destructor.
95    */
~DiffSolver()96   virtual ~DiffSolver () {}
97 
98   /**
99    * The initialization function.  This method is used to
100    * initialize internal data structures before a simulation begins.
101    */
102   virtual void init ();
103 
104   /**
105    * The reinitialization function.  This method is used after
106    * changes in the mesh.
107    */
108   virtual void reinit ();
109 
110   /**
111    * This method performs a solve.  What occurs in
112    * this method will depend on the type of solver.  See
113    * the subclasses for more details.
114    */
115   virtual unsigned int solve () = 0;
116 
117   /**
118    * \returns The number of "outer" (e.g. quasi-Newton) iterations
119    * required by the last solve.
120    */
total_outer_iterations()121   unsigned int total_outer_iterations() { return _outer_iterations; }
122 
123   /**
124    * \returns The number of "inner" (e.g. Krylov) iterations
125    * required by the last solve.
126    */
total_inner_iterations()127   unsigned int total_inner_iterations() { return _inner_iterations; }
128 
129   /**
130    * \returns The value of the SolveResult from the last solve.
131    */
solve_result()132   unsigned int solve_result() { return _solve_result; }
133 
134   /**
135    * \returns A constant reference to the system we are solving.
136    */
system()137   const sys_type & system () const { return _system; }
138 
139   /**
140    * \returns A writable reference to the system we are solving.
141    */
system()142   sys_type & system () { return _system; }
143 
144   /**
145    * Each linear solver step should exit after \p max_linear_iterations
146    * is exceeded.
147    */
148   unsigned int max_linear_iterations;
149 
150   /**
151    * The DiffSolver should exit in failure if \p max_nonlinear_iterations
152    * is exceeded and \p continue_after_max_iterations is false, or should
153    * end the nonlinear solve if \p max_nonlinear_iterations is exceeded and \p
154    * continue_after_max_iterations is true.
155    */
156   unsigned int max_nonlinear_iterations;
157 
158   /**
159    * The DiffSolver should not print anything to libMesh::out
160    * unless quiet is set to false; default is true.
161    */
162   bool quiet;
163 
164   /**
165    * The DiffSolver may print a lot more to libMesh::out
166    * if verbose is set to true; default is false.
167    */
168   bool verbose;
169 
170   /**
171    * Defaults to true, telling the DiffSolver to continue rather than exit when
172    * a solve has reached its maximum number of nonlinear iterations.
173    */
174   bool continue_after_max_iterations;
175 
176   /**
177    * Defaults to false, telling the DiffSolver to throw an error when
178    * the backtracking scheme fails to find a descent direction.
179    */
180   bool continue_after_backtrack_failure;
181 
182   /**
183    * The DiffSolver should exit after the residual is
184    * reduced to either less than absolute_residual_tolerance
185    * or less than relative_residual_tolerance times the
186    * initial residual.
187    *
188    * Users should increase any of these tolerances that they want to use for a
189    * stopping condition.
190    */
191   Real absolute_residual_tolerance;
192   Real relative_residual_tolerance;
193 
194   /**
195    * The DiffSolver should exit after the full nonlinear step norm is
196    * reduced to either less than absolute_step_tolerance
197    * or less than relative_step_tolerance times the largest
198    * nonlinear solution which has been seen so far.
199    *
200    * Users should increase any of these tolerances that they want to use for a
201    * stopping condition.
202    */
203   Real absolute_step_tolerance;
204   Real relative_step_tolerance;
205 
206   /**
207    * Any required linear solves will at first be done with this tolerance;
208    * the DiffSolver may tighten the tolerance for later solves.
209    */
210   double initial_linear_tolerance;
211 
212   /**
213    * The tolerance for linear solves is kept above this minimum
214    */
215   double minimum_linear_tolerance;
216 
217   /**
218    * Enumeration return type for the solve() function.  Multiple SolveResults
219    * may be combined (OR'd) in the single return.  To test which ones are present,
220    * just AND the return value with any of the SolveResult flags defined below.
221    */
222   enum SolveResult {
223     /**
224      * A default or invalid solve result.  This usually means
225      * no solve has occurred yet.
226      */
227     INVALID_SOLVE_RESULT = 0,
228 
229     /**
230      * The solver converged but no
231      * particular reason is specified.
232      */
233     CONVERGED_NO_REASON = 1,
234 
235     /**
236      * The DiffSolver achieved the desired
237      * absolute residual tolerance.
238      */
239     CONVERGED_ABSOLUTE_RESIDUAL = 2,
240 
241     /**
242      * The DiffSolver achieved the desired
243      * relative residual tolerance.
244      */
245     CONVERGED_RELATIVE_RESIDUAL = 4,
246 
247     /**
248      * The DiffSolver achieved the desired
249      * absolute step size tolerance.
250      */
251     CONVERGED_ABSOLUTE_STEP = 8,
252 
253     /**
254      * The DiffSolver achieved the desired
255      * relative step size tolerance.
256      */
257     CONVERGED_RELATIVE_STEP = 16,
258 
259     /**
260      * The DiffSolver diverged but no
261      * particular reason is specified.
262      */
263     DIVERGED_NO_REASON = 32,
264 
265     /**
266      * The DiffSolver reached the maximum allowed
267      * number of nonlinear iterations before satisfying
268      * any convergence tests.
269      */
270     DIVERGED_MAX_NONLINEAR_ITERATIONS = 64,
271 
272     /**
273      * The DiffSolver failed to find a descent direction
274      * by backtracking (See newton_solver.C)
275      */
276     DIVERGED_BACKTRACKING_FAILURE = 128,
277 
278     /**
279      * The linear solver used by the DiffSolver failed to
280      * find a solution.
281      */
282     DIVERGED_LINEAR_SOLVER_FAILURE = 256
283   };
284 
285   /**
286    * Pointer to functor which is called right after each linear solve
287    */
288   std::unique_ptr<LinearSolutionMonitor> linear_solution_monitor;
289 
290 protected:
291 
292   /**
293    * The largest solution norm which the DiffSolver has yet seen will be stored
294    * here, to be used for stopping criteria based on relative_step_tolerance
295    */
296   Real max_solution_norm;
297 
298   /**
299    * The largest nonlinear residual which the DiffSolver has yet seen will be
300    * stored here, to be used for stopping criteria based on
301    * relative_residual_tolerance
302    */
303   Real max_residual_norm;
304 
305   /**
306    * The number of outer iterations used by the last solve
307    */
308   unsigned int _outer_iterations;
309 
310   /**
311    * The number of inner iterations used by the last solve
312    */
313   unsigned int _inner_iterations;
314 
315   /**
316    * A reference to the system we are solving.
317    */
318   sys_type & _system;
319 
320   /**
321    * Initialized to zero.  solve_result is typically set internally in
322    * the solve() function before it returns.  When non-zero,
323    * solve_result tells the result of the latest solve.  See enum
324    * definition for description.
325    */
326   unsigned int _solve_result;
327 };
328 
329 
330 } // namespace libMesh
331 
332 
333 #endif // LIBMESH_DIFF_SOLVER_H
334