00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011 #include "CoinHelperFunctions.hpp"
00012 #include "CouenneProblem.hpp"
00013
00014
00015 bool CouenneProblem::checkNLP (const double *solution, double &obj, bool recompute) const {
00016
00017 if (Jnlst () -> ProduceOutput (Ipopt::J_ALL, J_PROBLEM)) {
00018
00019 printf ("checking solution: [%.12e] ", obj);
00020
00021 for (int i=0; i<nOrigVars_; i++)
00022 printf ("%.12e ", solution [i]);
00023 printf ("\n");
00024 }
00025
00026
00027
00028 for (register int i=0; i < nOrigVars_; i++) {
00029
00030 CouNumber val = solution [i];
00031
00032
00033
00034 exprVar *v = variables_ [i];
00035
00036 if ((v -> Type () == VAR) &&
00037 (v -> Multiplicity () > 0) &&
00038 (v -> isInteger ()) &&
00039 (fabs (val - COUENNE_round (val)) > feas_tolerance_)) {
00040
00041 Jnlst()->Printf(Ipopt::J_MOREVECTOR, J_PROBLEM,
00042 "checkNLP: integrality %d violated: %.6f [%g,%g]\n",
00043 i, val, domain_.lb (i), domain_.ub (i));
00044
00045 return false;
00046 }
00047 }
00048
00049 const int infeasible = 1;
00050 const int wrong_obj = 2;
00051
00052 CouNumber *sol = new CouNumber [nVars ()];
00053
00054
00055
00056 CoinCopyN (solution, nOrigVars_, sol);
00057 getAuxs (sol);
00058 CoinCopyN (solution, nOrigVars_, sol);
00059
00060
00061 domain_.push (nVars (), sol, domain_.lb (), domain_.ub (), false);
00062
00063
00064 if (Jnlst () -> ProduceOutput (Ipopt::J_ALL, J_PROBLEM)) {
00065 printf ("checknlp: %d vars -------------------\n", domain_.current () -> Dimension ());
00066 for (int i=0; i<domain_.current () -> Dimension (); i++)
00067 printf ("%4d %.12e [%.12e %.12e]\n",
00068 i, domain_.x (i), domain_.lb (i), domain_.ub (i));
00069 }
00070
00071 expression *objBody = Obj (0) -> Body ();
00072
00073
00074
00075
00076
00077
00078 CouNumber realobj = obj;
00079
00080 if (objBody)
00081 realobj =
00082 (objBody -> Index () >= 0) ?
00083 sol [objBody -> Index ()] :
00084 (*(objBody -> Image () ? objBody -> Image () : objBody)) ();
00085
00086 if (Jnlst () -> ProduceOutput (Ipopt::J_ALL, J_PROBLEM)) {
00087 printf ("%.12e %.12e %.12e ------------------------------\n",
00088 realobj, sol [objBody -> Index ()],
00089 (*(objBody -> Image () ? objBody -> Image () : objBody)) ());
00090 }
00091
00092 bool retval = true;
00093
00094 try {
00095
00096
00097
00098 if (fabs (realobj - obj) / (1. + fabs (realobj)) > feas_tolerance_) {
00099
00100 Jnlst()->Printf(Ipopt::J_MOREVECTOR, J_PROBLEM,
00101 "checkNLP: false objective. %g != %g (diff. %g)\n",
00102 realobj, obj, realobj - obj);
00103
00104 if (!recompute)
00105 throw wrong_obj;
00106 }
00107
00108 if (recompute)
00109 obj = realobj;
00110
00111 if (Jnlst () -> ProduceOutput (Ipopt::J_ALL, J_PROBLEM))
00112 printf ("recomputed: %.12e\n", obj);
00113
00114 for (int i=0; i < nOrigVars_; i++) {
00115
00116 if (variables_ [i] -> Multiplicity () <= 0)
00117 continue;
00118
00119 CouNumber val = domain_.x (i);
00120
00121
00122
00123 if ((val > domain_.ub (i) + feas_tolerance_) ||
00124 (val < domain_.lb (i) - feas_tolerance_)) {
00125
00126 Jnlst()->Printf(Ipopt::J_MOREVECTOR, J_PROBLEM,
00127 "checkNLP: variable %d out of bounds: %.6f [%g,%g] (diff %g)\n",
00128 i, val, domain_.lb (i), domain_.ub (i),
00129 CoinMax (fabs (val - domain_.lb (i)),
00130 fabs (val - domain_.ub (i))));
00131 throw infeasible;
00132 }
00133
00134
00135
00136 if (variables_ [i] -> isInteger () &&
00137 (fabs (val - COUENNE_round (val)) > feas_tolerance_)) {
00138
00139 Jnlst()->Printf(Ipopt::J_MOREVECTOR, J_PROBLEM,
00140 "checkNLP: integrality %d violated: %.6f [%g,%g]\n",
00141 i, val, domain_.lb (i), domain_.ub (i));
00142
00143 throw infeasible;
00144 }
00145 }
00146
00147
00148
00149 for (int i=0; i < nVars (); i++) {
00150
00151 if (Jnlst () -> ProduceOutput (Ipopt::J_ALL, J_PROBLEM)) {
00152 if (variables_ [i] -> Type () == AUX) {
00153 printf ("checking aux ");
00154 variables_ [i] -> print (); printf (" := ");
00155 variables_ [i] -> Image () -> print ();
00156 printf (" --- %.12e = %.12e [%.12e]; {",
00157 (*(variables_ [i])) (),
00158 (*(variables_ [i] -> Image ())) (),
00159 (*(variables_ [i])) () -
00160 (*(variables_ [i] -> Image ())) ());
00161
00162
00163 printf ("}\n");
00164 }
00165 }
00166
00167 CouNumber delta;
00168
00169
00170
00171 if ((variables_ [i] -> Type () == AUX) &&
00172 (fabs (delta = (*(variables_ [i])) () -
00173 (*(variables_ [i] -> Image ())) ()) > feas_tolerance_)) {
00174
00175 Jnlst()->Printf(Ipopt::J_MOREVECTOR, J_PROBLEM,
00176 "checkNLP: auxiliarized %d violated (%g)\n", i, delta);
00177
00178 throw infeasible;
00179 }
00180 }
00181
00182
00183
00184 for (int i=0; i < nCons (); i++) {
00185
00186 CouenneConstraint *c = Con (i);
00187
00188 CouNumber
00189 body = (*(c -> Body ())) (),
00190 lhs = (*(c -> Lb ())) (),
00191 rhs = (*(c -> Ub ())) ();
00192
00193 if (((rhs < COUENNE_INFINITY) &&
00194 (body > rhs + feas_tolerance_ * (1 + CoinMax (fabs (body), fabs (rhs))))) ||
00195 ((lhs > -COUENNE_INFINITY) &&
00196 (body < lhs - feas_tolerance_ * (1 + CoinMax (fabs (body), fabs (lhs)))))) {
00197
00198 if (Jnlst()->ProduceOutput(Ipopt::J_MOREVECTOR, J_PROBLEM)) {
00199
00200 printf ("checkNLP: constraint %d violated (lhs=%+e body=%+e rhs=%+e, violation %g): ",
00201 i, lhs, body, rhs, CoinMax (lhs-body, body-rhs));
00202
00203 c -> print ();
00204 }
00205
00206 throw infeasible;
00207 }
00208 }
00209 }
00210
00211 catch (int exception) {
00212
00213 switch (exception) {
00214
00215 case wrong_obj:
00216 retval = false;
00217 break;
00218
00219 case infeasible:
00220 default:
00221 retval = false;
00222 break;
00223 }
00224 }
00225
00226 delete [] sol;
00227 domain_.pop ();
00228
00229 return retval;
00230 }