1 // Copyright (C) 2008 International Business Machines and others. 2 // All Rights Reserved. 3 // This code is published under the Eclipse Public License. 4 // 5 // $Id$ 6 // 7 // Authors: Andreas Waechter IBM 2008-09-19 8 9 #include "IpInexactNormalTerminationTester.hpp" 10 #include "IpBlas.hpp" 11 12 #ifdef HAVE_CMATH 13 # include <cmath> 14 #else 15 # ifdef HAVE_MATH_H 16 # include <math.h> 17 # else 18 # error "don't have header file for math" 19 # endif 20 #endif 21 22 namespace Ipopt 23 { 24 #if COIN_IPOPT_VERBOSITY > 0 25 static const Index dbg_verbosity = 0; 26 #endif 27 InexactNormalTerminationTester()28 InexactNormalTerminationTester::InexactNormalTerminationTester() 29 { 30 DBG_START_METH("InexactNormalTerminationTester::InexactNormalTerminationTester",dbg_verbosity); 31 } 32 ~InexactNormalTerminationTester()33 InexactNormalTerminationTester::~InexactNormalTerminationTester() 34 { 35 DBG_START_METH("InexactNormalTerminationTester::~InexactNormalTerminationTester()",dbg_verbosity); 36 } 37 RegisterOptions(SmartPtr<RegisteredOptions> roptions)38 void InexactNormalTerminationTester::RegisterOptions(SmartPtr<RegisteredOptions> roptions) 39 { 40 roptions->AddLowerBoundedNumberOption( 41 "inexact_normal_tol", 42 "Desired relative residual tolerance for iterative solver during normal step computation.", 43 0.0, true, 1e-3, 44 ""); 45 roptions->AddLowerBoundedIntegerOption( 46 "inexact_normal_max_iter", 47 "Maximal number of iterative solver iterations during normal step computation", 48 0, 200, 49 ""); 50 } 51 52 InitializeImpl(const OptionsList & options,const std::string & prefix)53 bool InexactNormalTerminationTester::InitializeImpl(const OptionsList& options, 54 const std::string& prefix) 55 { 56 options.GetNumericValue("inexact_normal_tol", inexact_normal_tol_, prefix); 57 options.GetIntegerValue("inexact_normal_max_iter", 58 inexact_normal_max_iter_, prefix); 59 60 std::string inexact_linear_system_scaling; 61 options.GetStringValue("inexact_linear_system_scaling", 62 inexact_linear_system_scaling, prefix); 63 if (inexact_linear_system_scaling=="slack-based") { 64 requires_scaling_ = true; 65 } 66 else { 67 requires_scaling_ = false; 68 } 69 70 c_Avc_norm_cauchy_ = -1; 71 return true; 72 } 73 InitializeSolve()74 bool InexactNormalTerminationTester::InitializeSolve() 75 { 76 DBG_START_METH("InexactNormalTerminationTester::InitializeSolve", 77 dbg_verbosity); 78 79 80 return true; 81 } 82 83 InexactNormalTerminationTester::ETerminationTest 84 InexactNormalTerminationTester:: TestTermination(Index ndim,const Number * sol,const Number * resid,Index iter,Number norm2_rhs)85 TestTermination(Index ndim, const Number* sol, const Number* resid, 86 Index iter, Number norm2_rhs) 87 { 88 DBG_START_METH("InexactNormalTerminationTester::TestTermination", 89 dbg_verbosity); 90 91 last_iter_ = iter; 92 93 DBG_ASSERT(c_Avc_norm_cauchy_ >= 0.); 94 95 ETerminationTest retval = CONTINUE; 96 97 double norm2_resid = IpBlasDnrm2(ndim, resid, 1); 98 Jnlst().Printf(J_DETAILED, J_LINEAR_ALGEBRA, 99 "TTNormal: iter = %d ||resid|| = %23.16e ||rhs|| = %23.16e\n", 100 iter, norm2_resid, norm2_rhs); 101 102 if (iter > inexact_normal_max_iter_) { 103 IpData().Append_info_string("N+"); 104 retval = OTHER_SATISFIED; 105 return retval; 106 } 107 108 //if (Min(norm2_resid/norm2_rhs,norm2_resid) > inexact_normal_tol_) { 109 if (norm2_resid/norm2_rhs > inexact_normal_tol_ && norm2_resid > 1e-10) { 110 return retval; 111 } 112 113 // compute ||c+Av|| for current iterative solution v 114 // note that the sol_x and sol_s have the wrong sign 115 SmartPtr<const Vector> sol_x; 116 SmartPtr<const Vector> sol_s; 117 SmartPtr<const Vector> sol_c; 118 SmartPtr<const Vector> sol_d; 119 GetVectors(ndim, sol, sol_x, sol_s, sol_c, sol_d); 120 121 if (requires_scaling_) { 122 SmartPtr<const Vector> scaling_vec = InexCq().curr_scaling_slacks(); 123 SmartPtr<Vector> tmp = sol_s->MakeNewCopy(); 124 tmp->ElementWiseMultiply(*scaling_vec); 125 sol_s = ConstPtr(tmp); 126 } 127 128 SmartPtr<Vector> inf_c = IpCq().curr_c()->MakeNewCopy(); 129 IpCq().curr_jac_c()->MultVector(-1., *sol_x, 1., *inf_c); 130 SmartPtr<const Vector> curr_d_minus_s = IpCq().curr_d_minus_s(); 131 SmartPtr<Vector> inf_d = curr_d_minus_s->MakeNew(); 132 inf_d->AddTwoVectors(1., *curr_d_minus_s, 1., *sol_s, 0.); 133 IpCq().curr_jac_d()->MultVector(-1., *sol_x, 1., *inf_d); 134 135 Number trial_c_Av_norm = IpCq().CalcNormOfType(NORM_2, *inf_c, *inf_d); 136 137 Jnlst().Printf(J_MOREDETAILED, J_LINEAR_ALGEBRA, 138 "TTNormal: c_Avc_norm_cauchy = %23.16e trial_c_Av_norm = %23.16e\n", c_Avc_norm_cauchy_, trial_c_Av_norm); 139 Number BasVal = Max(1., IpCq().curr_primal_infeasibility(NORM_2)); 140 141 if (Compare_le(trial_c_Av_norm, c_Avc_norm_cauchy_, BasVal)) { 142 retval = OTHER_SATISFIED; 143 } 144 145 return retval; 146 } 147 148 void Clear()149 InexactNormalTerminationTester::Clear() 150 { 151 DBG_START_METH("InexactNormalTerminationTester::Clear", 152 dbg_verbosity); 153 } 154 155 } // namespace Ipopt 156