Line data Source code
1 : /*
2 : * Copyright (c) 2013-2015: G-CSC, Goethe University Frankfurt
3 : * Author: Raphael Prohl
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 : /*
34 : * (main parts are based on the structure of
35 : * newton.h by Andreas Vogel)
36 : */
37 :
38 : #ifndef __H__UG__LIB_DISC__OPERATOR__NON_LINEAR_OPERATOR__NL_JACOBI__NL_JACOBIL_H_
39 : #define __H__UG__LIB_DISC__OPERATOR__NON_LINEAR_OPERATOR__NL_JACOBI__NL_JACOBIL_H_
40 :
41 : #include "lib_algebra/operator/interface/operator_inverse.h"
42 :
43 : // modul intern headers
44 : #include "lib_disc/assemble_interface.h"
45 : #include "lib_disc/operator/non_linear_operator/assembled_non_linear_operator.h"
46 : #include "lib_disc/operator/linear_operator/assembled_linear_operator.h"
47 :
48 : namespace ug {
49 :
50 : /// Nonlinear Jacobi-method
51 : /**
52 : * Let L(u) denote a nonlinear functional of n components (l_1,...,l_n).
53 : * Then the basic step of the nonlinear Jacobi method is to solve the
54 : * i-th equation
55 : *
56 : * l_i(u_1^{k},...,u_{i-1}^{k},u_i,u_{i+1}^{k},...,u_{n}^{k}) = 0
57 : *
58 : * for u_i and to set u_i^{k+1} = u_i. Here k denotes the iteration-index.
59 : * Thus, in order to obtain u^{k+1} from u^k, we solve successively the n
60 : * dimensional nonlinear equations for i = 1,...,n. Here this is done
61 : * by a scalar newton step for every i. But every other scalar nonlinear
62 : * method could be applied as well.
63 : *
64 : * Using a damped version of the nonlinear jacobi method results in the
65 : * following update of the variables
66 : *
67 : * u_i^{k+1} = u_i^k + damp * (u_i -u_i^k).
68 : *
69 : * References:
70 : * <ul>
71 : * <li> J. M. Ortega and W. C. Rheinbolt. Iterative Solution of nonlinear equations in several variables.(1970)
72 : * </ul>
73 : *
74 : * \tparam TAlgebra Algebra type
75 : */
76 : template <typename TAlgebra>
77 : class NLJacobiSolver
78 : : public IOperatorInverse<typename TAlgebra::vector_type>,
79 : public DebugWritingObject<TAlgebra>
80 : {
81 : public:
82 : /// Algebra type
83 : typedef TAlgebra algebra_type;
84 :
85 : /// Vector type
86 : typedef typename TAlgebra::vector_type vector_type;
87 :
88 : /// Matrix type
89 : typedef typename TAlgebra::matrix_type matrix_type;
90 :
91 : protected:
92 : typedef DebugWritingObject<TAlgebra> base_writer_type;
93 :
94 : public:
95 : /// default constructor
96 : NLJacobiSolver();
97 :
98 : /// constructor
99 : NLJacobiSolver(SmartPtr<IConvergenceCheck<vector_type> > spConvCheck);
100 :
101 : void set_convergence_check(SmartPtr<IConvergenceCheck<vector_type> > spConvCheck);
102 :
103 0 : void set_damp(number damp) {m_damp = damp;}
104 :
105 : /// returns information about configuration parameters
106 0 : virtual std::string config_string() const
107 : {
108 0 : std::stringstream ss;
109 0 : ss << "NonlinearJacobiSolver( damp = " << m_damp << ")\n";
110 0 : ss << " ConvergenceCheck: ";
111 0 : if(m_spConvCheck.valid()) ss << ConfigShift(m_spConvCheck->config_string()) << "\n";
112 0 : else ss << " NOT SET!\n";
113 :
114 0 : return ss.str();
115 :
116 0 : }
117 :
118 : ///////////////////////////////////////////////////////////////////////////
119 : // OperatorInverse interface methods
120 : ///////////////////////////////////////////////////////////////////////////
121 :
122 : /// This operator inverts the Operator op: Y -> X
123 : virtual bool init(SmartPtr<IOperator<vector_type> > op);
124 :
125 : /// prepare Operator
126 : virtual bool prepare(vector_type& u);
127 :
128 : /// apply Operator, i.e. op^{-1}(0) = u
129 : virtual bool apply(vector_type& u);
130 :
131 : private:
132 : /// help functions for debug output
133 : /// \{
134 : void write_debug(const vector_type& vec, const char* filename);
135 : void write_debug(const matrix_type& mat, const char* filename);
136 : /// \}
137 :
138 : private:
139 : SmartPtr<IConvergenceCheck<vector_type> > m_spConvCheck;
140 :
141 : /// damping factor
142 : number m_damp;
143 :
144 : SmartPtr<AssembledOperator<algebra_type> > m_spAssOp;
145 : SmartPtr<AssembledLinearOperator<algebra_type> > m_spJ;
146 : SmartPtr<IAssemble<TAlgebra> > m_spAss;
147 :
148 : /// call counter
149 : int m_dgbCall;
150 : };
151 :
152 : }
153 :
154 : #include "nl_jacobi_impl.h"
155 :
156 : #endif /* __H__UG__LIB_DISC__OPERATOR__NON_LINEAR_OPERATOR__NL_JACOBI__NL_JACOBIL_H_ */
|
