ug4
newton.h
Go to the documentation of this file.
1 /*
2  * Copyright (c) 2010-2015: G-CSC, Goethe University Frankfurt
3  * Author: Andreas Vogel
4  *
5  * This file is part of UG4.
6  *
7  * UG4 is free software: you can redistribute it and/or modify it under the
8  * terms of the GNU Lesser General Public License version 3 (as published by the
9  * Free Software Foundation) with the following additional attribution
10  * requirements (according to LGPL/GPL v3 §7):
11  *
12  * (1) The following notice must be displayed in the Appropriate Legal Notices
13  * of covered and combined works: "Based on UG4 (www.ug4.org/license)".
14  *
15  * (2) The following notice must be displayed at a prominent place in the
16  * terminal output of covered works: "Based on UG4 (www.ug4.org/license)".
17  *
18  * (3) The following bibliography is recommended for citation and must be
19  * preserved in all covered files:
20  * "Reiter, S., Vogel, A., Heppner, I., Rupp, M., and Wittum, G. A massively
21  * parallel geometric multigrid solver on hierarchically distributed grids.
22  * Computing and visualization in science 16, 4 (2013), 151-164"
23  * "Vogel, A., Reiter, S., Rupp, M., Nägel, A., and Wittum, G. UG4 -- a novel
24  * flexible software system for simulating pde based models on high performance
25  * computers. Computing and visualization in science 16, 4 (2013), 165-179"
26  *
27  * This program is distributed in the hope that it will be useful,
28  * but WITHOUT ANY WARRANTY; without even the implied warranty of
29  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
30  * GNU Lesser General Public License for more details.
31  */
32 
33 #ifndef __H__UG__LIB_DISC__OPERATOR__NON_LINEAR_OPERATOR__NEWTON_SOLVER__NEWTON__
34 #define __H__UG__LIB_DISC__OPERATOR__NON_LINEAR_OPERATOR__NEWTON_SOLVER__NEWTON__
35 
36 #include <cmath>
37 
41 
42 // modul intern headers
46 #include "../line_search.h"
49 
50 #include "nestedNewtonRFSwitch.h"
51 
52 #if ENABLE_NESTED_NEWTON_RESOLFUNC_UPDATE
53 
54 #include "newtonUpdaterGeneric.h"
55 
56 #endif
57 
58 namespace ug {
59 
61 template <typename TAlgebra>
63  : public IOperatorInverse<typename TAlgebra::vector_type>,
64  public DebugWritingObject<TAlgebra>
65 {
66  public:
68  typedef TAlgebra algebra_type;
69 
71  typedef typename TAlgebra::vector_type vector_type;
72 
74  typedef typename TAlgebra::matrix_type matrix_type;
75 
76  public:
79 
82 
84  NewtonSolver();
85 
89  SmartPtr<ILineSearch<vector_type> > spLineSearch);
90 
93 
96 
98  void set_line_search(SmartPtr<ILineSearch<vector_type> > spLineSearch) {m_spLineSearch = spLineSearch;}
101 
103  virtual bool init(SmartPtr<IOperator<vector_type> > N);
104 
106  virtual bool prepare(vector_type& u);
107 
109  virtual bool apply(vector_type& u);
110 
112 
119  virtual std::string config_string() const;
120 
122  void print_average_convergence() const;
123 
126  size_t num_newton_steps() const;
127  int num_linsolver_calls(size_t call) const;
128  int num_linsolver_steps(size_t call) const;
129  double average_linear_steps(size_t call) const;
130  int total_linsolver_calls() const;
131  int total_linsolver_steps() const;
132  double total_average_linear_steps() const;
133  int last_num_newton_steps() const {return m_lastNumSteps;}
135 
138 
141  {m_innerStepUpdate.push_back(NU);}
142 
145  {m_innerStepUpdate.clear();}
146 
149  {m_stepUpdate.push_back(NU);}
150 
153  {m_stepUpdate.clear();}
154 
156  void set_reassemble_J_freq(int freq)
157  {m_reassembe_J_freq = freq;};
158 
159 #if ENABLE_NESTED_NEWTON_RESOLFUNC_UPDATE
160  void setNewtonUpdater( SmartPtr<NewtonUpdaterGeneric<vector_type> > nU )
161  {
162  m_newtonUpdater = nU;
163  }
164 #endif
165 
166  private:
169  void write_debug(const vector_type& vec, std::string filename);
170  void write_debug(const matrix_type& mat, std::string filename);
172 
173  private:
176 
179 
182 
184  std::vector<SmartPtr<INewtonUpdate> > m_innerStepUpdate;
185  std::vector<SmartPtr<INewtonUpdate> > m_stepUpdate;
186 
195 
199 
202  std::vector<int> m_vTotalLinSolverSteps;
203  std::vector<int> m_vLinSolverCalls;
204  std::vector<number> m_vNonLinSolverRates;
205  std::vector<number> m_vLinSolverRates;
207 
208 #if ENABLE_NESTED_NEWTON_RESOLFUNC_UPDATE
209 
211 
212 #endif
213 
214 };
215 
216 }
217 
218 #include "newton_impl.h"
219 
220 #endif /* __H__UG__LIB_DISC__OPERATOR__NON_LINEAR_OPERATOR__NEWTON_SOLVER__NEWTON__ */
Definition: smart_pointer.h:108
Definition: debug_writer.h:354
Interface providing Jacobian and Defect of a discretization.
Definition: assemble_interface.h:110
Definition: convergence_check.h:72
Definition: line_search.h:71
describes an inverse linear mapping X->Y
Definition: linear_operator_inverse.h:80
describes a mapping X->Y
Definition: operator.h:86
describes an inverse mapping X->Y
Definition: operator_inverse.h:79
linear solver using abstract preconditioner interface
Definition: linear_solver.h:57
Newton solver for assembling based discretizations.
Definition: newton.h:65
std::vector< SmartPtr< INewtonUpdate > > m_innerStepUpdate
Update.
Definition: newton.h:184
int total_linsolver_calls() const
Definition: newton_impl.h:433
void add_step_update(SmartPtr< INewtonUpdate > NU)
add outer step update (applied before every Newton step)
Definition: newton.h:148
std::vector< SmartPtr< INewtonUpdate > > m_stepUpdate
Definition: newton.h:185
void clear_average_convergence()
resets average linear solver convergence
Definition: newton_impl.h:458
void clear_inner_step_update(SmartPtr< INewtonUpdate > NU)
clears inner step update
Definition: newton.h:144
virtual bool prepare(vector_type &u)
prepare Operator
Definition: newton_impl.h:161
SmartPtr< IConvergenceCheck< vector_type > > m_spConvCheck
Convergence Check.
Definition: newton.h:178
SmartPtr< ILineSearch< vector_type > > line_search()
Definition: newton.h:100
void disable_line_search()
Definition: newton.h:99
void set_reassemble_J_freq(int freq)
sets the frequency of reassembling of the Jacobian (0 == 1 == in every step, i.e. classically)
Definition: newton.h:156
virtual bool apply(vector_type &u)
apply Operator, i.e. N^{-1}(0) = u
Definition: newton_impl.h:168
NewtonSolver()
default constructor
Definition: newton_impl.h:81
virtual bool init(SmartPtr< IOperator< vector_type > > N)
This operator inverts the Operator N: Y -> X.
Definition: newton_impl.h:149
int num_linsolver_steps(size_t call) const
Definition: newton_impl.h:421
void write_debug(const vector_type &vec, std::string filename)
Definition: newton_impl.h:467
int m_dgbCall
call counter
Definition: newton.h:197
TAlgebra::matrix_type matrix_type
Matrix type.
Definition: newton.h:74
int m_lastNumSteps
Definition: newton.h:198
int last_num_newton_steps() const
Definition: newton.h:133
void clear_step_update(SmartPtr< INewtonUpdate > NU)
clears outer step update
Definition: newton.h:152
void set_line_search(SmartPtr< ILineSearch< vector_type > > spLineSearch)
sets the line search
Definition: newton.h:98
double total_average_linear_steps() const
Definition: newton_impl.h:451
std::vector< number > m_vLinSolverRates
Definition: newton.h:205
std::vector< number > m_vNonLinSolverRates
Definition: newton.h:204
virtual std::string config_string() const
returns information about configuration parameters
Definition: newton_impl.h:489
std::vector< int > m_vTotalLinSolverSteps
Definition: newton.h:202
SmartPtr< ILinearOperatorInverse< vector_type > > m_spLinearSolver
linear solver
Definition: newton.h:175
SmartPtr< ILineSearch< vector_type > > m_spLineSearch
LineSearch.
Definition: newton.h:181
double average_linear_steps(size_t call) const
Definition: newton_impl.h:427
int m_reassembe_J_freq
how often to reassemble the Jacobian (0 == 1 == in every step, i.e. classically)
Definition: newton.h:194
TAlgebra::vector_type vector_type
Vector type.
Definition: newton.h:71
SmartPtr< AssembledLinearOperator< algebra_type > > m_J
jacobi operator
Definition: newton.h:190
SmartPtr< AssembledOperator< algebra_type > > m_N
assembling routine
Definition: newton.h:188
TAlgebra algebra_type
Algebra type.
Definition: newton.h:68
void print_average_convergence() const
prints average linear solver convergence
Definition: newton_impl.h:379
size_t num_newton_steps() const
Definition: newton_impl.h:409
void set_convergence_check(SmartPtr< IConvergenceCheck< vector_type > > spConvCheck)
sets the convergence check
Definition: newton_impl.h:140
SmartPtr< IAssemble< TAlgebra > > m_spAss
assembling
Definition: newton.h:192
std::vector< int > m_vLinSolverCalls
Definition: newton.h:203
void set_linear_solver(SmartPtr< ILinearOperatorInverse< vector_type > > LinearSolver)
sets the linear solver
Definition: newton.h:92
int num_linsolver_calls(size_t call) const
Definition: newton_impl.h:415
void add_inner_step_update(SmartPtr< INewtonUpdate > NU)
add inner step update (applied before every linear solver step)
Definition: newton.h:140
int total_linsolver_steps() const
Definition: newton_impl.h:442
Definition: newtonUpdaterGeneric.h:25
Wrapper for sequential matrices to handle them in parallel.
Definition: parallel_matrix.h:65
const NullSmartPtr SPNULL
The equivalent to NULL for smart pointers.
Definition: smart_pointer.h:90
CPUAlgebra::vector_type vector_type
the ug namespace