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