1 | // Copyright (C) 2004, 2005 International Business Machines and others. |
---|
2 | // All Rights Reserved. |
---|
3 | // This code is published under the Common Public License. |
---|
4 | // |
---|
5 | // $Id: IpPDSystemSolver.hpp 501 2005-08-26 15:43:07Z claird $ |
---|
6 | // |
---|
7 | // Authors: Carl Laird, Andreas Waechter IBM 2004-08-13 |
---|
8 | |
---|
9 | #ifndef __IPPDSYSTEMSOLVER_HPP__ |
---|
10 | #define __IPPDSYSTEMSOLVER_HPP__ |
---|
11 | |
---|
12 | #include "IpUtils.hpp" |
---|
13 | #include "IpSymMatrix.hpp" |
---|
14 | #include "IpAlgStrategy.hpp" |
---|
15 | #include "IpIteratesVector.hpp" |
---|
16 | |
---|
17 | namespace Ipopt |
---|
18 | { |
---|
19 | |
---|
20 | /** Pure Primal Dual System Solver Base Class. |
---|
21 | * This is the base class for all derived Primal-Dual System Solver Types. |
---|
22 | * |
---|
23 | * Here, we understand the primal-dual system as the following linear |
---|
24 | * system: |
---|
25 | * |
---|
26 | * \f$ |
---|
27 | * \left[\begin{array}{cccccccc} |
---|
28 | * W & 0 & J_c^T & J_d^T & -P^x_L & P^x_U & 0 & 0 \\ |
---|
29 | * 0 & 0 & 0 & -I & 0 & 0 & -P_L^d & P_U^d \\ |
---|
30 | * J_c & 0 & 0 & 0 & 0 & 0 & 0 & 0\\ |
---|
31 | * J_d & -I & 0 & 0 & 0 & 0 & 0 & 0\\ |
---|
32 | * Z_L(P_L^x)^T & 0 & 0 & 0 & Sl^x_L & 0 & 0 & 0\\ |
---|
33 | * -Z_U(P_U^x)^T & 0 & 0 & 0 & 0 & Sl^x_U & 0 & 0\\ |
---|
34 | * 0 & V_L(P_L^d)^T & 0 & 0 & 0 & 0 & Sl^s_L & 0 \\ |
---|
35 | * 0 & -V_U(P_U^d)^T & 0 & 0 & 0 & 0 & 0 & Sl^s_U \\ |
---|
36 | * \end{array}\right] |
---|
37 | * \left(\begin{array}{c} |
---|
38 | * sol_x\\ sol_s\\ sol_c\\ sol_d\\ sol^z_L\\ sol^z_U\\ sol^v_L\\ |
---|
39 | * sol^v_U |
---|
40 | * \end{array}\right) = |
---|
41 | * \left(\begin{array}{c} |
---|
42 | * rhs_x\\ rhs_s\\ rhs_c\\ rhs_d\\ rhs^z_L\\ rhs^z_U\\ rhs^v_L\\ |
---|
43 | * rhs^v_U |
---|
44 | * \end{array}\right) |
---|
45 | * \f$ |
---|
46 | * |
---|
47 | * Here, \f$Sl^x_L = (P^x_L)^T x - x_L\f$, |
---|
48 | * \f$Sl^x_U = x_U - (P^x_U)^T x\f$, \f$Sl^d_L = (P^d_L)^T d(x) - d_L\f$, |
---|
49 | * \f$Sl^d_U = d_U - (P^d_U)^T d(x)\f$. The results returned to the |
---|
50 | * caller is \f$res = \alpha * sol + \beta * res\f$. |
---|
51 | * |
---|
52 | * The solution of this linear system (in order to compute the search |
---|
53 | * direction of the algorthim) usually requires a considerable amount of |
---|
54 | * computation time. Therefore, it is important to tailor the solution |
---|
55 | * of this system to the characteristics of the problem. The purpose of |
---|
56 | * this base class is to provide a generic interface to the algorithm |
---|
57 | * that it can use whenever it requires a solution of the above system. |
---|
58 | * Particular implementation can then be written to provide the methods |
---|
59 | * defined here. |
---|
60 | * |
---|
61 | * It is implicitly assumed here, that the upper left 2 by 2 block |
---|
62 | * is possibly modified (implicitly or explicitly) so that its |
---|
63 | * projection onto the null space of the overall constraint |
---|
64 | * Jacobian \f$\left[\begin{array}{cc}J_c & 0\\J_d & |
---|
65 | * -I\end{array}\right]\f$ is positive definite. This is necessary |
---|
66 | * to guarantee certain descent properties of the resulting search |
---|
67 | * direction. For example, in the full space implementation, a |
---|
68 | * multiple of the identity might be added to the upper left 2 by 2 |
---|
69 | * block. |
---|
70 | * |
---|
71 | * Note that the Solve method might be called several times for different |
---|
72 | * right hand sides, but with identical data. Therefore, if possible, |
---|
73 | * an implemetation of PDSystem should check whether the incoming data has |
---|
74 | * changed, and not redo factorization etc. unless necessary. |
---|
75 | */ |
---|
76 | class PDSystemSolver: public AlgorithmStrategyObject |
---|
77 | { |
---|
78 | public: |
---|
79 | /** @name /Destructor */ |
---|
80 | //@{ |
---|
81 | /** Default Constructor */ |
---|
82 | PDSystemSolver() |
---|
83 | {} |
---|
84 | |
---|
85 | /** Default destructor */ |
---|
86 | virtual ~PDSystemSolver() |
---|
87 | {} |
---|
88 | //@} |
---|
89 | |
---|
90 | /** overloaded from AlgorithmStrategyObject */ |
---|
91 | virtual bool InitializeImpl(const OptionsList& options, |
---|
92 | const std::string& prefix) = 0; |
---|
93 | |
---|
94 | /** Solve the primal dual system, given one right hand side. For |
---|
95 | * the time being, we don't pass a state object, but all items |
---|
96 | * individually, including the calculated quantity \f$\Sigma = |
---|
97 | * P_L^TS_L^{-1}Z_L + P_U^TS_U^{-1}Z_U\f$. If the flag |
---|
98 | * allow_inexact is set to true, it is not necessary to solve the |
---|
99 | * system to best accuracy; for example, we don't want iterative |
---|
100 | * refinement during the computation of the second order |
---|
101 | * correction. |
---|
102 | */ |
---|
103 | virtual void Solve(Number alpha, |
---|
104 | Number beta, |
---|
105 | const IteratesVector& rhs, |
---|
106 | IteratesVector& res, |
---|
107 | bool allow_inexact=false) =0; |
---|
108 | |
---|
109 | private: |
---|
110 | /**@name Default Compiler Generated Methods |
---|
111 | * (Hidden to avoid implicit creation/calling). |
---|
112 | * These methods are not implemented and |
---|
113 | * we do not want the compiler to implement |
---|
114 | * them for us, so we declare them private |
---|
115 | * and do not define them. This ensures that |
---|
116 | * they will not be implicitly created/called. */ |
---|
117 | //@{ |
---|
118 | /** Overloaded Equals Operator */ |
---|
119 | PDSystemSolver& operator=(const PDSystemSolver&); |
---|
120 | //@} |
---|
121 | }; |
---|
122 | |
---|
123 | |
---|
124 | } // namespace Ipopt |
---|
125 | |
---|
126 | #endif |
---|