00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011 #include "BonminConfig.h"
00012
00013 #include "BonBqpdSolver.hpp"
00014 #include "BonBqpdWarmStart.hpp"
00015
00016 #include "CoinTime.hpp"
00017
00018 #define InitializeAll
00019
00020 typedef Bonmin::BqpdSolver::fint fint;
00021 typedef Bonmin::BqpdSolver::real real;
00022
00023 extern "C"
00024 {
00025 void F77_FUNC(bqpd,BQPD)(fint* n, fint* m, fint* k, fint* kmax,
00026 real* a, fint* la, real* x, real* bl, real* bu,
00027 real* f, real* fmin, real* g, real* r, real* w,
00028 real* e, fint* ls, real* alp, fint* lp, fint* mlp,
00029 fint* peq, real* ws, fint* lws, fint* m0de,
00030 fint* ifail, fint* info, fint* iprint, fint* nout);
00031
00032
00033 extern struct {
00034 fint kk,ll,kkk,lll,mxws,mxlws;
00035 }
00036 F77_FUNC(wsc,WSC);
00037
00038 extern struct {
00039 real eps,tol,emin;
00040 }
00041 F77_FUNC(epsc,EPSC);
00042
00043 extern struct {
00044 real sgnf;
00045 fint nrep,npiv,nres;
00046 }
00047 F77_FUNC(repc,REPC);
00048
00049 extern struct {
00050 fint nup,nfreq;
00051 }
00052 F77_FUNC(refactorc,REFACTORC);
00053
00054 extern struct {
00055 real vstep;
00056 }
00057 F77_FUNC(vstepc,VSTEPC);
00058
00059 extern struct {
00060 fint phl, phr, phc;
00061 }
00062 F77_FUNC(hessc,HESSC);
00063
00064 extern struct {
00065 fint scale_mode, phe;
00066 }
00067 F77_FUNC(scalec,SCALEC);
00068
00069 extern struct {
00070 fint irh1,na,na1,nb,nb1,ka1,kb1,kc1,irg1,lu1,lv,lv1,ll1;
00071 }
00072 F77_FUNC(bqpdc,BQPDC);
00073
00074 extern struct {
00075 real alpha;
00076 }
00077 F77_FUNC(alphac,ALPHAC);
00078
00079 extern struct {
00080 fint ns,ns1,nt,nt1,nu,nu1,nx,nx1,np,np1,nprof,lc;
00081 fint lc1,li,li1,lm,lm1,lp_,lp1,lq,lq1,lr,lr1,ls_,ls1,lt,lt1;
00082 }
00083 F77_FUNC(sparsec,SPARSEC);
00084
00085 extern struct {
00086 fint m1,m2,mp,mq,lastr,irow;
00087 }
00088 F77_FUNC(factorc,FACTORC);
00089
00090 extern struct {
00091 fint mxm1;
00092 }
00093 F77_FUNC(mxm1c,MXM1C);
00094
00095 extern struct {
00096 real c;
00097 }
00098 F77_FUNC(minorc,MINORS);
00099 }
00100
00101 namespace Bonmin
00102 {
00103
00104 std::string BqpdSolver::solverName_ = "Bqpd QP";
00105
00106 void BqpdSolver::
00107 registerOptions(Ipopt::SmartPtr<Bonmin::RegisteredOptions> roptions)
00108 {
00109 roptions->SetRegisteringCategory("Bqpd options",RegisteredOptions::BqpdCategory);
00110 roptions->AddLowerBoundedNumberOption("qp_fillin_factor", "Factor for estimating fill-in for factorization method in Bqpd", 0., true, 4.);
00111 }
00112
00113 BqpdSolver::BqpdSolver(bool createEmpty )
00114 :
00115 TNLPSolver(),
00116 cached_(NULL)
00117 {
00118 if (createEmpty) return;
00119 }
00120
00121 BqpdSolver::BqpdSolver(Ipopt::SmartPtr<Bonmin::RegisteredOptions> roptions,
00122 Ipopt::SmartPtr<Ipopt::OptionsList> options,
00123 Ipopt::SmartPtr<Ipopt::Journalist> journalist)
00124 :
00125 TNLPSolver(roptions, options, journalist),
00126 cached_(NULL)
00127 {
00128 options->GetNumericValue("qp_fillin_factor", fillin_factor_, "bqpd.");
00129 options->GetIntegerValue("kmax", kmax_ipt_, "bqpd.");
00130 options->GetIntegerValue("mlp", mlp_ipt_,"bqpd.");
00131 }
00132
00133 Ipopt::SmartPtr <TNLPSolver>
00134 BqpdSolver::clone()
00135 {
00136 Ipopt::SmartPtr<BqpdSolver> retval = new BqpdSolver(true);
00137 *retval->options_ = *options_;
00138 retval->roptions_ = roptions_;
00139 retval->journalist_ = journalist_;
00140 retval->fillin_factor_ = fillin_factor_;
00141 retval->kmax_ipt_ = kmax_ipt_;
00142 retval->mlp_ipt_ = mlp_ipt_;
00143 return GetRawPtr(retval);
00144 }
00145
00146 BqpdSolver::~BqpdSolver()
00147 {}
00148
00149 bool
00150 BqpdSolver::Initialize(std::string optFile)
00151 {
00152 std::ifstream is;
00153 if (optFile != "") {
00154 try {
00155 is.open(optFile.c_str());
00156 }
00157 catch (std::bad_alloc) {
00158 journalist_->Printf(Ipopt::J_SUMMARY, Ipopt::J_MAIN, "\nEXIT: Not enough memory.\n");
00159 return false;
00160 }
00161 catch (...) {
00162 Ipopt::IpoptException E("Unknown Exception caught in ipopt", "Unknown File", -1);
00163 E.ReportException(*journalist_);
00164 return false;
00165 }
00166 }
00167 bool retval = Initialize(is);
00168 if (is) {
00169 is.close();
00170 }
00171 return retval;
00172 }
00173
00174 bool
00175 BqpdSolver::Initialize(std::istream &is)
00176 {
00177 if (is.good()) {
00178 options_->ReadFromStream(*journalist_, is);
00179 }
00180 return true;
00181 }
00182
00184 TNLPSolver::ReturnStatus
00185 BqpdSolver::OptimizeTNLP(const Ipopt::SmartPtr<Ipopt::TNLP>& tnlp)
00186 {
00187 BranchingTQP* tqp = dynamic_cast<BranchingTQP*>(GetRawPtr(tnlp));
00188 if (!tqp) {
00189 Ipopt::IpoptException E("BqpdSolver called with object other than a BranchingTQP",
00190 "BonBqpdSolver.cpp",-1);
00191 throw E;
00192 }
00193 if (IsNull(cached_) || !cached_->use_warm_start_in_cache_) {
00194 cached_ = new cachedInfo(tqp, options_, kmax_ipt_, mlp_ipt_,
00195 &fillin_factor_);
00196 }
00197 return callOptimizer();
00198 }
00199
00201 TNLPSolver::ReturnStatus
00202 BqpdSolver::ReOptimizeTNLP(const Ipopt::SmartPtr<Ipopt::TNLP> & tnlp)
00203 {
00204 #ifndef NDEBUG
00205 BranchingTQP* tqp = dynamic_cast<BranchingTQP*>(GetRawPtr(tnlp));
00206 assert(tqp == GetRawPtr(cached_->tqp_));
00207 #endif
00208 int n = cached_->n;
00209 int m = cached_->m;
00210 tnlp->get_bounds_info(n, cached_->bl, cached_->bu,
00211 m, cached_->bl+n, cached_->bu+n);
00212
00213 for (int i=0; i<n+m; i++) {
00214 cached_->bl[i] = Max(cached_->bl[i], -1e50);
00215 cached_->bu[i] = Min(cached_->bu[i], 1e50);
00216 }
00217
00218 cached_->use_warm_start_in_cache_ = true;
00219 return callOptimizer();
00220 }
00221
00222 void
00223 BqpdSolver::cachedInfo::initialize(const Ipopt::SmartPtr<BranchingTQP> & tqp,
00224 Ipopt::SmartPtr<Ipopt::OptionsList>& options,
00225 int kmax_ipt,
00226 int mlp_ipt,
00227 double* fillin_factor)
00228 {
00229
00230 fillin_factor_ = fillin_factor;
00231
00232
00233
00234 F77_FUNC(epsc,EPSC).eps = 1111e-19;
00235 F77_FUNC(epsc,EPSC).tol = 1e-12;
00236 F77_FUNC(epsc,EPSC).emin = 1.;
00237 F77_FUNC(repc,REPC).sgnf = 1e-4;
00238 F77_FUNC(repc,REPC).nrep = 2;
00239 F77_FUNC(repc,REPC).npiv = 3;
00240 F77_FUNC(repc,REPC).nres = 2;
00241 F77_FUNC(refactorc,REFACTORC).nfreq = 500;
00242
00243 Ipopt::TNLP::IndexStyleEnum index_style;
00244 Index nv, nc, nnz_jac_g, nnz_hess;
00245 tqp->get_nlp_info(nv, nc, nnz_jac_g, nnz_hess, index_style);
00246 n = nv;
00247 m = nc;
00248
00249 if (kmax_ipt == -1) {
00250 kmax = n;
00251 }
00252 else {
00253 kmax = kmax_ipt;
00254 kmax = std::min(kmax,n);
00255 }
00256 mlp = mlp_ipt;
00257
00258 use_warm_start_in_cache_ = false;
00259
00260
00261 x = new real[n];
00262 bl = new real[n+m];
00263 bu = new real[n+m];
00264 g = new real[n];
00265 r = new real[n+m];
00266 w = new real[n+m];
00267 e = new real[n+m];
00268 ls = new fint[n+m];
00269 alp = new real[mlp];
00270 lp = new fint[mlp];
00271
00272
00273 tqp->get_bounds_info(n, bl, bu, m, bl+n, bu+n);
00274
00275
00276 for (int i=0; i<n+m; i++) {
00277 bl[i] = Max(bl[i], -1e50);
00278 bu[i] = Min(bu[i], 1e50);
00279 }
00280
00281
00282
00283 const Number* obj_grad = tqp->ObjGrad();
00284 amax_ = nnz_jac_g;
00285 for (int i = 0; i<n; i++) {
00286 if (obj_grad[i]!=0.) {
00287 amax_++;
00288 }
00289 }
00290 int lamax = amax_+m+2;
00291 a = new real[amax_];
00292 la = new fint [1+lamax];
00293
00294
00295 int nnz_grad = 0;
00296 for (int i = 0; i < n ; i++) {
00297 if (obj_grad[i]!=0.) {
00298 a[nnz_grad] = obj_grad[i];
00299 la[++nnz_grad] = i+1;
00300 }
00301 }
00302 la[amax_+1] = 1;
00303
00304
00305 const Number* JacVals = tqp->ConstrJacVals();
00306 const Index* RowJac = tqp->ConstrJacIRow();
00307 const Index* ColJac = tqp->ConstrJacJCol();
00308
00309 int* permutationJac = new int [nnz_jac_g];
00310 FilterSolver::TMat2RowPMat(false, n, m, nnz_jac_g, RowJac, ColJac, permutationJac,
00311 la, nnz_grad, 1, TNLP::C_STYLE);
00312 for (int i=0; i<nnz_jac_g; i++) {
00313 const int& indice = permutationJac[i];
00314 a[nnz_grad+i] = JacVals[indice];
00315 }
00316 delete [] permutationJac;
00317 #if 0
00318
00319 printf("nnz_grad = %d nnz_jac = %d\n", nnz_grad, nnz_jac_g);
00320 for (int i=0; i<1+lamax; i++) printf("la[%2d] = %d\n", i,la[i]);
00321 for (int i=0; i<amax_; i++) printf("a[%3d] = %e\n",i,a[i]);
00322 #endif
00323
00324
00325 mxws = nnz_hess + ((kmax*(kmax+9))/2 + 2*n + m) + (5*n + (int)(*fillin_factor_*(double)amax_));
00326 mxlws = (nnz_hess + n + 2) + kmax + (9*n + m);
00327
00328 ws = new real[mxws];
00329 lws = new fint[mxlws];
00330
00331 #ifdef InitializeAll
00332 for (int i=0; i<mxws; i++) {
00333 ws[i] = 42.;
00334 }
00335 for (int i=0; i<mxlws; i++) {
00336 lws[i] = 55;
00337 }
00338 #endif
00339
00340
00341 const Number* HessVals = tqp->ObjHessVals();
00342 const Index* RowHess = tqp->ObjHessIRow();
00343 const Index* ColHess = tqp->ObjHessJCol();
00344
00345 kk = nnz_hess;
00346 ll = nnz_hess + n + 2;
00347 int* permutationHess = new int[nnz_hess];
00348
00349 FilterSolver::TMat2RowPMat(true, n, n, nnz_hess, RowHess, ColHess,
00350 permutationHess, lws, 0, 0, TNLP::C_STYLE);
00351 for (int i=0; i<nnz_hess; i++) {
00352 ws[i] = HessVals[permutationHess[i]];
00353 }
00354 delete [] permutationHess;
00355 #if 0
00356
00357 printf("nnz_hess = %d\n", nnz_hess);
00358 for (int i=0; i<ll; i++) printf("hess lws[%2d] = %d\n", i,lws[i]);
00359 for (int i=0; i<kk; i++) printf("hess ws[%3d] = %e\n",i,ws[i]);
00360 #endif
00361
00362 Index bufy;
00363 options->GetIntegerValue("iprint",bufy, "bqpd.");
00364 iprint = bufy;
00365 nout = 6;
00366
00367 tqp_ = tqp;
00368 }
00369
00371 TNLPSolver::ReturnStatus
00372 BqpdSolver::callOptimizer()
00373 {
00374 cached_->optimize();
00375
00376 TNLPSolver::ReturnStatus optimizationStatus = TNLPSolver::exception;
00377 Ipopt::SolverReturn status = Ipopt::INTERNAL_ERROR;
00378 fint ifail = cached_->ifail;
00379 switch (ifail) {
00380 case 0:
00381 optimizationStatus = TNLPSolver::solvedOptimal;
00382 status = Ipopt::SUCCESS;
00383 break;
00384 case 1:
00385 optimizationStatus = TNLPSolver::unbounded;
00386 status = Ipopt::DIVERGING_ITERATES;
00387 case 2:
00388 case 3:
00389 optimizationStatus = TNLPSolver::provenInfeasible;
00390 status = Ipopt::LOCAL_INFEASIBILITY;
00391 break;
00392 }
00393
00394 Index dummy_len = Ipopt::Max(cached_->n,cached_->m);
00395 Number* dummy = new Number[dummy_len];
00396 for (int i=0; i<dummy_len; i++) {
00397 dummy[i] = 0.;
00398 }
00399
00400 cached_->tqp_->finalize_solution(status, cached_->n,
00401 cached_->x, dummy, dummy,
00402 cached_->m, dummy, dummy,
00403 cached_->f, NULL, NULL);
00404 delete [] dummy;
00405 return optimizationStatus;
00406 }
00407
00408 BqpdSolver::cachedInfo::~cachedInfo()
00409 {
00410 if (haveHotStart_) {
00411 unmarkHotStart();
00412 }
00413 delete [] a;
00414 delete [] la;
00415 delete [] x;
00416 delete [] bl;
00417 delete [] bu;
00418 delete [] g;
00419 delete [] r;
00420 delete [] w;
00421 delete [] e;
00422 delete [] ls;
00423 delete [] alp;
00424 delete [] lp;
00425 delete [] ws;
00426 delete [] lws;
00427 }
00428
00429 bool
00430 BqpdSolver::cachedInfo::markHotStart()
00431 {
00432 haveHotStart_ = true;
00433 irh1 = F77_FUNC(bqpdc,BQPDC).irh1;
00434 na = F77_FUNC(bqpdc,BQPDC).na;
00435 na1 = F77_FUNC(bqpdc,BQPDC).na1;
00436 nb = F77_FUNC(bqpdc,BQPDC).nb;
00437 nb1 = F77_FUNC(bqpdc,BQPDC).nb1;
00438 ka1 = F77_FUNC(bqpdc,BQPDC).ka1;
00439 kb1 = F77_FUNC(bqpdc,BQPDC).kb1;
00440 kc1 = F77_FUNC(bqpdc,BQPDC).kc1;
00441 irg1 = F77_FUNC(bqpdc,BQPDC).irg1;
00442 lu1 = F77_FUNC(bqpdc,BQPDC).lu1;
00443 lv = F77_FUNC(bqpdc,BQPDC).lv;
00444 lv1 = F77_FUNC(bqpdc,BQPDC).lv1;
00445 ll1 = F77_FUNC(bqpdc,BQPDC).ll1;
00446 eps = F77_FUNC(epsc,EPSC).eps;
00447 tol = F77_FUNC(epsc,EPSC).tol;
00448 emin = F77_FUNC(epsc,EPSC).emin;
00449 vstep = F77_FUNC(vstepc,VSTEPC).vstep;
00450 sgnf = F77_FUNC(repc,REPC).sgnf;
00451 nrep = F77_FUNC(repc,REPC).nrep;
00452 npiv = F77_FUNC(repc,REPC).npiv;
00453 nres = F77_FUNC(repc,REPC).nres;
00454 nup = F77_FUNC(refactorc,REFACTORC).nup;
00455 nfreq = F77_FUNC(refactorc,REFACTORC).nfreq;
00456 alpha = F77_FUNC(alphac,ALPHAC).alpha;
00457
00458 ns = F77_FUNC(sparsec,SPARSEC).ns;
00459 ns1 = F77_FUNC(sparsec,SPARSEC).ns1;
00460 nt = F77_FUNC(sparsec,SPARSEC).nt;
00461 nt1 = F77_FUNC(sparsec,SPARSEC).nt1;
00462 nu = F77_FUNC(sparsec,SPARSEC).nu;
00463 nu1 = F77_FUNC(sparsec,SPARSEC).nu1;
00464 nx = F77_FUNC(sparsec,SPARSEC).nx;
00465 nx1 = F77_FUNC(sparsec,SPARSEC).nx1;
00466 np = F77_FUNC(sparsec,SPARSEC).np;
00467 np1 = F77_FUNC(sparsec,SPARSEC).np1;
00468 nprof = F77_FUNC(sparsec,SPARSEC).nprof;
00469 lc = F77_FUNC(sparsec,SPARSEC).lc;
00470 lc1 = F77_FUNC(sparsec,SPARSEC).lc1;
00471 li = F77_FUNC(sparsec,SPARSEC).li;
00472 li1 = F77_FUNC(sparsec,SPARSEC).li1;
00473 lm = F77_FUNC(sparsec,SPARSEC).lm;
00474 lm1 = F77_FUNC(sparsec,SPARSEC).lm1;
00475 lp_ = F77_FUNC(sparsec,SPARSEC).lp_;
00476 lp1 = F77_FUNC(sparsec,SPARSEC).lp1;
00477 lq = F77_FUNC(sparsec,SPARSEC).lq;
00478 lq1 = F77_FUNC(sparsec,SPARSEC).lq1;
00479 lr = F77_FUNC(sparsec,SPARSEC).lr;
00480 lr1 = F77_FUNC(sparsec,SPARSEC).lr1;
00481 ls_ = F77_FUNC(sparsec,SPARSEC).ls_;
00482 ls1 = F77_FUNC(sparsec,SPARSEC).ls1;
00483 lt = F77_FUNC(sparsec,SPARSEC).lt;
00484 lt1 = F77_FUNC(sparsec,SPARSEC).lt1;
00485
00486 m1 = F77_FUNC(factorc,FACTORC).m1;
00487 m2 = F77_FUNC(factorc,FACTORC).m2;
00488 mp = F77_FUNC(factorc,FACTORC).mp;
00489 mq = F77_FUNC(factorc,FACTORC).mq;
00490 lastr = F77_FUNC(factorc,FACTORC).lastr;
00491 irow = F77_FUNC(factorc,FACTORC).irow;
00492
00493 mxm1 = F77_FUNC(mxm1c,MXM1c).mxm1;
00494 c = F77_FUNC(minorc,MINORC).c;
00495 kkkHot = F77_FUNC(wsc,WSC).kkk;
00496 lllHot = F77_FUNC(wsc,WSC).lll;
00497
00498 kHot = k;
00499 xHot = CoinCopyOfArray(x,n);
00500 fHot = f;
00501 gHot = CoinCopyOfArray(g,n);
00502 rHot = CoinCopyOfArray(r,n+m);
00503 wHot = CoinCopyOfArray(w,n+m);
00504 eHot = CoinCopyOfArray(e,n+m);
00505 lsHot = CoinCopyOfArray(ls,n+m);
00506 alpHot = CoinCopyOfArray(alp,mlp);
00507 lpHot = CoinCopyOfArray(lp,mlp);
00508 peqHot = peq;
00509 wsHot = CoinCopyOfArray(ws,mxws);
00510 lwsHot = CoinCopyOfArray(lws,mxlws);
00511 infoHot[0] = info[0];
00512
00513 return true;
00514 }
00515
00516 void
00517 BqpdSolver::cachedInfo::copyFromHotStart()
00518 {
00519 F77_FUNC(bqpdc,BQPDC).irh1 = irh1;
00520 F77_FUNC(bqpdc,BQPDC).na = na;
00521 F77_FUNC(bqpdc,BQPDC).na1 = na1;
00522 F77_FUNC(bqpdc,BQPDC).nb = nb;
00523 F77_FUNC(bqpdc,BQPDC).nb1 = nb1;
00524 F77_FUNC(bqpdc,BQPDC).ka1 = ka1;
00525 F77_FUNC(bqpdc,BQPDC).kb1 = kb1;
00526 F77_FUNC(bqpdc,BQPDC).kc1 = kc1;
00527 F77_FUNC(bqpdc,BQPDC).irg1 = irg1;
00528 F77_FUNC(bqpdc,BQPDC).lu1 = lu1;
00529 F77_FUNC(bqpdc,BQPDC).lv = lv;
00530 F77_FUNC(bqpdc,BQPDC).lv1 = lv1;
00531 F77_FUNC(bqpdc,BQPDC).ll1 = ll1;
00532 F77_FUNC(epsc,EPSC).eps = eps;
00533 F77_FUNC(epsc,EPSC).tol = tol;
00534 F77_FUNC(epsc,EPSC).emin = emin;
00535 F77_FUNC(vstepc,VSTEPC).vstep = vstep;
00536 F77_FUNC(repc,REPC).sgnf = sgnf;
00537 F77_FUNC(repc,REPC).nrep = nrep;
00538 F77_FUNC(repc,REPC).npiv = npiv;
00539 F77_FUNC(repc,REPC).nres = nres;
00540 F77_FUNC(refactorc,REFACTORC).nup = nup;
00541 F77_FUNC(refactorc,REFACTORC).nfreq = nfreq;
00542 F77_FUNC(alphac,ALPHAC).alpha = alpha;
00543
00544 F77_FUNC(sparsec,SPARSEC).ns = ns;
00545 F77_FUNC(sparsec,SPARSEC).ns1 = ns1;
00546 F77_FUNC(sparsec,SPARSEC).nt = nt;
00547 F77_FUNC(sparsec,SPARSEC).nt1 = nt1;
00548 F77_FUNC(sparsec,SPARSEC).nu = nu;
00549 F77_FUNC(sparsec,SPARSEC).nu1 = nu1;
00550 F77_FUNC(sparsec,SPARSEC).nx = nx;
00551 F77_FUNC(sparsec,SPARSEC).nx1 = nx1;
00552 F77_FUNC(sparsec,SPARSEC).np = np;
00553 F77_FUNC(sparsec,SPARSEC).np1 = np1;
00554 F77_FUNC(sparsec,SPARSEC).nprof = nprof;
00555 F77_FUNC(sparsec,SPARSEC).lc = lc;
00556 F77_FUNC(sparsec,SPARSEC).lc1 = lc1;
00557 F77_FUNC(sparsec,SPARSEC).li = li;
00558 F77_FUNC(sparsec,SPARSEC).li1 = li1;
00559 F77_FUNC(sparsec,SPARSEC).lm = lm;
00560 F77_FUNC(sparsec,SPARSEC).lm1 = lm1;
00561 F77_FUNC(sparsec,SPARSEC).lp_ = lp_;
00562 F77_FUNC(sparsec,SPARSEC).lp1 = lp1;
00563 F77_FUNC(sparsec,SPARSEC).lq = lq;
00564 F77_FUNC(sparsec,SPARSEC).lq1 = lq1;
00565 F77_FUNC(sparsec,SPARSEC).lr = lr;
00566 F77_FUNC(sparsec,SPARSEC).lr1 = lr1;
00567 F77_FUNC(sparsec,SPARSEC).ls_ = ls_;
00568 F77_FUNC(sparsec,SPARSEC).ls1 = ls1;
00569 F77_FUNC(sparsec,SPARSEC).lt = lt;
00570 F77_FUNC(sparsec,SPARSEC).lt1 = lt1;
00571
00572 F77_FUNC(factorc,FACTORC).m1 = m1;
00573 F77_FUNC(factorc,FACTORC).m2 = m2;
00574 F77_FUNC(factorc,FACTORC).mp = mp;
00575 F77_FUNC(factorc,FACTORC).mq = mq;
00576 F77_FUNC(factorc,FACTORC).lastr = lastr;
00577 F77_FUNC(factorc,FACTORC).irow = irow;
00578
00579 F77_FUNC(mxm1c,MXM1c).mxm1 = mxm1;
00580 F77_FUNC(minorc,MINORC).c = c;
00581
00582 F77_FUNC(wsc,WSC).kkk = kkkHot;
00583 F77_FUNC(wsc,WSC).lll = lllHot;
00584
00585 k = kHot;
00586 CoinCopyN(xHot,n,x);
00587 f = fHot;
00588 CoinCopyN(gHot,n,g);
00589 CoinCopyN(rHot,n+m,r);
00590 CoinCopyN(wHot,n+m,w);
00591 CoinCopyN(eHot,n+m,e);
00592 CoinCopyN(lsHot,n+m,ls);
00593 CoinCopyN(alpHot,mlp,alp);
00594 CoinCopyN(lpHot,mlp,lp);
00595 peq = peqHot;
00596 CoinCopyN(wsHot,mxws,ws);
00597 CoinCopyN(lwsHot,mxlws,lws);
00598 info[0] = infoHot[0];
00599 }
00600
00601 void
00602 BqpdSolver::cachedInfo::unmarkHotStart()
00603 {
00604 delete [] xHot;
00605 delete [] gHot;
00606 delete [] rHot;
00607 delete [] wHot;
00608 delete [] eHot;
00609 delete [] lsHot;
00610 delete [] alpHot;
00611 delete [] lpHot;
00612 delete [] wsHot;
00613 delete [] lwsHot;
00614 haveHotStart_ = false;
00615 }
00616
00618 void
00619 BqpdSolver::cachedInfo::optimize()
00620 {
00621
00622 F77_FUNC(scalec,SCALEC).scale_mode = 0;
00623 F77_FUNC(scalec,SCALEC).phe = 0;
00624 F77_FUNC(hessc,HESSC).phl = 1;
00625 F77_FUNC(wsc,WSC).kk = kk;
00626 F77_FUNC(wsc,WSC).ll = ll;
00627 F77_FUNC(wsc,WSC).mxws = mxws;
00628 F77_FUNC(wsc,WSC).mxlws = mxlws;
00629
00630 if (use_warm_start_in_cache_ && !bad_warm_start_info_) {
00631 m0de = 6;
00632 ifail = 0;
00633 use_warm_start_in_cache_ = false;
00634 if (haveHotStart_) copyFromHotStart();
00635 }
00636 else {
00637 m0de = 0;
00638 tqp_->get_starting_point(n, 1, x, 0, NULL, NULL, m, 0, NULL);
00639 ifail = 0;
00640 bad_warm_start_info_ = false;
00641 }
00642
00643 #if 0
00644 printf("========= 222222222222 =============\n");
00645 printf("kk = %d ll = %d mxws = %d mxlws = %d\n", kk, ll, mxws, mxlws);
00646 for (int i=0; i<n; i++) {
00647 printf("xL[%3d] = %15.8e xU[%3d] = %15.8e\n", i, bl[i], i, bu[i]);
00648 }
00649 for (int i=0; i<m; i++) {
00650 printf("gL[%3d] = %15.8e gU[%3d] = %15.8e\n", i, bl[n+i], i, bu[n+i]);
00651 }
00652 #endif
00653 cpuTime_ = - CoinCpuTime();
00654 real fmin = -1e100;
00655 #if 0
00656 for (int i=0; i<n; i++) {
00657 printf("qxstart[%2d] = %23.16e\n", i, x[i]);
00658 }
00659 #endif
00660 F77_FUNC(bqpd,BQPD)(&n, &m, &k, &kmax, a, la, x, bl, bu, &f, &fmin,
00661 g, r, w, e, ls, alp, lp, &mlp, &peq, ws, lws,
00662 &m0de, &ifail, info, &iprint, &nout);
00663 if (ifail == 8 && m0de == 6) {
00664 printf("Restarting Bqpd...");
00665 m0de = 0;
00666 tqp_->get_starting_point(n, 1, x, 0, NULL, NULL, m, 0, NULL);
00667 ifail = 0;
00668 F77_FUNC(bqpd,BQPD)(&n, &m, &k, &kmax, a, la, x, bl, bu, &f, &fmin,
00669 g, r, w, e, ls, alp, lp, &mlp, &peq, ws, lws,
00670 &m0de, &ifail, info, &iprint, &nout);
00671 printf("new ifail = %d\n", ifail);
00672 }
00673 while (ifail == 7) {
00674
00675 if (haveHotStart_) unmarkHotStart();
00676
00677 printf("Increasing fillin_factor from %e ", *fillin_factor_);
00678 *fillin_factor_ *= 2.;
00679 printf("to %e\n", *fillin_factor_);
00680 int mxws_new = kk + F77_FUNC(wsc,WSC).kkk + (5*n + (int)(*fillin_factor_*(double)amax_));
00681 real* ws_new = new real[mxws_new];
00682 #ifdef InitializeAll
00683 for (int i=0; i<mxws_new; i++) {
00684 ws_new[i] = 42.;
00685 }
00686 #endif
00687 CoinCopyN(ws, kk, ws_new);
00688 delete [] ws;
00689 ws = ws_new;
00690 mxws = mxws_new;
00691 F77_FUNC(wsc,WSC).mxws = mxws;
00692
00693 m0de = 0;
00694 tqp_->get_starting_point(n, 1, x, 0, NULL, NULL, m, 0, NULL);
00695 ifail = 0;
00696 F77_FUNC(bqpd,BQPD)(&n, &m, &k, &kmax, a, la, x, bl, bu, &f, &fmin,
00697 g, r, w, e, ls, alp, lp, &mlp, &peq, ws, lws,
00698 &m0de, &ifail, info, &iprint, &nout);
00699 }
00700 if (ifail == 8) bad_warm_start_info_ = true;
00701 #if 0
00702 for (int i=0; i<n; i++) {
00703 printf("qxsol[%2d] = %23.16e\n", i, x[i]);
00704 }
00705 #endif
00706 #if 0
00707 printf("ifail = %d\n", ifail);
00708 printf("final f = %e\n", f);
00709 printf("final f + obj_val = %e\n", f+tqp_->ObjVal());
00710 #endif
00711 #if 0
00712 int kkk = F77_FUNC(wsc,WSC).kkk;
00713 int lll = F77_FUNC(wsc,WSC).lll;
00714 printf("========= 3333333333 =============\n");
00715 printf("kk = %d ll = %d kkk = %d lll = %d mxws = %d mxlws = %d\n", kk, ll, kkk, lll, mxws, mxlws);
00716 for (int i=0; i<kk+kkk; i++) {
00717 printf("ws[%3d] = %15.8e\n ", i, ws[i]);
00718 }
00719 printf("--------\n");
00720 for (int i=kk+kkk; i<mxws; i++) {
00721 printf("ws[%3d] = %15.8e\n ", i, ws[i]);
00722 }
00723 for (int i=0; i<ll+lll; i++) {
00724 printf("lws[%5d] = %8d\n", i, lws[i]);
00725 }
00726 printf("------\n");
00727 for (int i=ll+lll; i<mxlws; i++) {
00728 printf("lws[%5d] = %8d\n", i, lws[i]);
00729 }
00730 #endif
00731 cpuTime_ += CoinCpuTime();
00732 }
00733
00734 std::string
00735 BqpdSolver::UnsolvedBqpdError::errorNames_[1] =
00736 {"Internal error in Filter SQP."};
00737
00738 std::string
00739 BqpdSolver::UnsolvedBqpdError::solverName_ =
00740 "filterSqp";
00741
00742 const std::string&
00743 BqpdSolver::UnsolvedBqpdError::errorName() const
00744 {
00745 return errorNames_[0];
00746 }
00747
00748 const std::string&
00749 BqpdSolver::UnsolvedBqpdError::solverName() const
00750 {
00751 return solverName_;
00752 }
00753
00754 bool
00755 BqpdSolver::setWarmStart(const CoinWarmStart * warm,
00756 Ipopt::SmartPtr<TMINLP2TNLP> tnlp)
00757 {
00758 #if 0
00759 if (IsNull(cached_)) {
00760 cached_ = new cachedInfo(GetRawPtr(tnlp), options_);
00761 }
00762
00763 const FilterWarmStart * warmF = dynamic_cast<const FilterWarmStart *> (warm);
00764
00765 const fint xsize = warmF->xSize();
00766 real* x = cached_->x;
00767 const real* xarray = warmF->xArray();
00768 for (int i = 0; i<xsize; i++) {
00769 x[i] = xarray[i];
00770 }
00771 CoinCopyN(warmF->lamArray(), warmF->lamSize(), cached_->lam);
00772 CoinCopyN(warmF->lwsArray(), warmF->lwsSize(), cached_->lws);
00773 for (int i = 0 ; i < 14 ; i ++) {
00774 cached_->istat[i] = warmF->istat()[i];
00775 }
00776 cached_->use_warm_start_in_cache_ = true;
00777 #endif
00778 printf("BqpdSolver::setWarmStart called!\n");
00779 return true;
00780 }
00781
00782 CoinWarmStart *
00783 BqpdSolver::getWarmStart(Ipopt::SmartPtr<TMINLP2TNLP> tnlp) const
00784 {
00785 #if 0
00786 return new FilterWarmStart(cached_->n, cached_->x,
00787 cached_->n+cached_->m, cached_->lam,
00788 cached_->maxiWk, cached_->lws, cached_->istat);
00789 #endif
00790 printf("BqpdSolver::getWarmStart called!\n");
00791 return NULL;
00792 }
00793
00794 CoinWarmStart *
00795 BqpdSolver::getEmptyWarmStart() const
00796 {
00797 #if 0
00798 return new FilterWarmStart;
00799 #endif
00800 printf("BqpdSolver::getEmptyWarmStart called \n");
00801 return NULL;
00802 }
00803
00805 bool
00806 BqpdSolver::warmStartIsValid(const CoinWarmStart * ws) const{
00807 const BqpdWarmStart* bws = dynamic_cast<const BqpdWarmStart*>(ws);
00808 if (bws && ! bws->empty()) {
00809 return true;
00810 }
00811 return false;
00812 }
00813
00814 }