44 #ifndef ROL_BUNDLE_TT_H 45 #define ROL_BUNDLE_TT_H 52 #include "Teuchos_RCP.hpp" 61 #include "Teuchos_SerialDenseMatrix.hpp" 62 #include "Teuchos_SerialDenseVector.hpp" 63 #include "Teuchos_LAPACK.hpp" 68 #define FIRST_VIOLATED 0 84 Bundle_TT (
const unsigned maxSize = 10,
const Real coeff = 0.0,
const unsigned remSize = 2)
85 :
Bundle<Real>(maxSize,coeff,remSize), maxSize_(maxSize) {
86 INF = std::numeric_limits<double>::max();
87 maxind = std::numeric_limits<int>::max();
88 Id.reshape(maxSize_,maxSize_);
93 Real
GiTGj(
const int i,
const int j){
111 Teuchos::SerialDenseMatrix<int, Real>
L;
112 Teuchos::SerialDenseMatrix<int, Real>
Id;
113 Teuchos::SerialDenseVector<int, Real>
tempv;
114 Teuchos::SerialDenseVector<int, Real>
tempw1;
115 Teuchos::SerialDenseVector<int, Real>
tempw2;
116 Teuchos::SerialDenseVector<int, Real>
lh;
117 Teuchos::SerialDenseVector<int, Real>
lj;
118 Teuchos::SerialDenseVector<int, Real>
z1;
119 Teuchos::SerialDenseVector<int, Real>
z2;
129 void swapRowsL(
unsigned ind1,
unsigned ind2,
bool trans=
false){
136 for (
unsigned n=ind1+1; n<=ind2; n++){
137 Teuchos::SerialDenseMatrix<int, Real> Id_n(Teuchos::Copy,Id,currSize_,currSize_);
138 Id_n(dd,dd) = 0; Id_n(dd,n) = 1.0;
139 Id_n(n,dd) = 1.0; Id_n(n,n) = 0;
140 Teuchos::SerialDenseMatrix<int, Real> prod(currSize_,currSize_);
142 prod.multiply(Teuchos::NO_TRANS,Teuchos::NO_TRANS,1.0,Id_n,L,0.0);
144 prod.multiply(Teuchos::NO_TRANS,Teuchos::NO_TRANS,1.0,L,Id_n,0.0);
151 if (currSize_ <= dependent_) {
155 Real tmpdiagMax = -
INF;
156 Real tmpdiagMin =
INF;
157 for (
unsigned j=0;j<currSize_-
dependent_;j++){
158 if(
L(j,j) > tmpdiagMax ){
162 if(
L(j,j) < tmpdiagMin ){
167 kappa = tmpdiagMax/tmpdiagMin;
175 if(dependent_ && (ind == currSize_-1)){
178 tmp = base[currSize_-2];
179 base[currSize_-2] = base[currSize_-1];
180 base[currSize_-1] = tmp;
183 std::cout <<
"Swapped last two rows of L " << std::endl;
190 unsigned zsize = ind+1;
191 z1.resize(zsize); z2.resize(zsize);
192 z1[ind] = ( 1.0 -
lhTz1 ) / delta;
193 z2[ind] = ( this->
alpha(base[ind]) -
lhTz2 ) / delta;
197 if(delta >
L(LiMax,LiMax)){
199 kappa = delta/
L(LiMin,LiMin);
201 if(delta <
L(LiMin,LiMin)){
203 kappa =
L(LiMax,LiMax)/delta;
209 if (ind >= currSize_-dependent_){
211 std::cout <<
"Eliminating dependent item " << base[ind] << std::endl;
214 if (ind < currSize_-1){
216 std::cout <<
"Eliminating Lh but keeping Lj" << std::endl;
220 base[ind] = base[currSize_-1];
228 L.reshape(currSize_,currSize_);
229 base.resize(currSize_);
232 std::cout <<
"New base = " << std::endl;
233 for (
unsigned kk=0;kk<base.size();kk++){
234 std::cout << base[kk] << std::endl;
237 std::cout <<
"New L " << std::endl;
238 std::cout << L << std::endl;
246 std::cout <<
"Eliminating INdependent item " << base[baseitem] << std::endl;
266 for (
unsigned j=ind+1;j<currSize_-
dependent_;j++){
268 if (std::abs(ai) <= tol*currSize_)
274 if (std::abs(aj) <= tol*currSize_){
275 Gc = 0.0; d = std::abs(ai); Gs = (ai < 0) - (ai > 0);
277 else if ( std::abs(ai) > std::abs(aj) ){
279 Real sgnb = (ai > 0) - (ai < 0);
280 Real u = sgnb * std::sqrt(1.0 + t*t );
287 Real sgna = (aj > 0) - (aj < 0);
288 Real u = sgna * std::sqrt(1.0 + t*t );
299 L(j,j) = d;
L(j,ind) = 0.0;
302 Real tmp1 =
L(h,ind);
304 L(h,ind) = Gc*tmp1 + Gs*tmp2;
305 L(h,j) = Gc*tmp2 - Gs*tmp1;
312 z1[ind] = Gc*tmp1 + Gs*tmp2;
313 z1[j] = Gc*tmp2 - Gs*tmp1;
314 z2[ind] = Gc*tmp3 + Gs*tmp4;
315 z2[j] = Gc*tmp4 - Gs*tmp3;
319 std::cout <<
"After Givens: L,z1,z2 " << std::endl;
320 std::cout << L << std::endl;
321 std::cout << z1 << std::endl;
322 std::cout << z2 << std::endl;
326 deltaLh =
L(currSize_-dependent_,ind);
328 deltaLj =
L(currSize_-1,ind);
334 L.reshape(currSize_-1,currSize_-1);
339 for(
unsigned m=ind; m<zsize; m++ ){
347 std::cout <<
"After elimination: L,z1,z2 " << std::endl;
348 std::cout << L << std::endl;
349 std::cout << z1 << std::endl;
350 std::cout << z2 << std::endl;
354 base.erase(base.begin()+ind);
357 std::cout <<
"New base = " << std::endl;
358 for (
unsigned kk=0;kk<base.size();kk++){
359 std::cout << base[kk] << std::endl;
371 std::cout <<
"deltaLh = " << deltaLh << std::endl;
374 Real ghNorm =
GiTGj(base[currSize_-dependent_],base[currSize_-dependent_]);
377 for (
unsigned ii=0;ii<currSize_-
dependent_;ii++){
378 lhnrm +=
L(currSize_-dependent_,ii)*
L(currSize_-dependent_,ii);
380 deltaLh = std::abs(ghNorm - lhnrm);
382 Real sgn1 = (deltaLh > 0) ? 1 : ((deltaLh < 0) ? -1 : 0);
383 Real sgn2 = (deltaLj > 0) ? 1 : ((deltaLj < 0) ? -1 : 0);
384 Real signum = sgn1 * sgn2;
385 deltaLh = std::abs( ghNorm - lhNorm + deltaLh * deltaLh);
389 std::cout <<
"ghNorm = " << ghNorm << std::endl;
390 std::cout <<
"lhNorm (exact) = " << lhnrm << std::endl;
391 std::cout <<
"lhNorm = " << lhNorm << std::endl;
392 std::cout <<
"deltaLh = " << std::sqrt(deltaLh) << std::endl;
393 std::cout <<
"kappa = " << kappa << std::endl;
396 if( std::sqrt(deltaLh) > tol*kappa*std::max(1.0,ghNorm) ){
398 std::cout <<
"Lh has become lin. INdependent" << std::endl;
406 for (
unsigned ii=0;ii<newind;ii++){
407 lh[ii] =
L(newind,ii);
408 lhTz1 += lh[ii]*z1[ii];
409 lhTz2 += lh[ii]*z2[ii];
411 deltaLh = std::sqrt(deltaLh);
416 Real gjNorm =
GiTGj(base[currSize_-1],base[currSize_-1]);
420 deltaLj = std::abs(gjNorm - ljNorm);
422 deltaLj = - std::sqrt( deltaLj );
424 deltaLj = std::sqrt( deltaLj );
427 Real gjTgh =
GiTGj(base[currSize_-1],base[currSize_-2]);
430 ljTlh +=
L(currSize_-1,ii)*
L(currSize_-2,ii);
432 deltaLj = (gjTgh - ljTlh) / deltaLh;
435 L(currSize_-1,currSize_-2) =
deltaLj;
439 std::cout <<
"Updated L, z1, z2: " << std::endl;
440 std::cout << L << std::endl;
441 std::cout << z1 << std::endl;
442 std::cout << z2 << std::endl;
443 std::cout <<
"kappa = " << kappa << std::endl;
450 Real gjNorm =
GiTGj(base[currSize_-1],base[currSize_-1]);
454 ljnrm +=
L(currSize_-1,ii)*
L(currSize_-1,ii);
456 deltaLj = std::abs(gjNorm - ljnrm);
458 deltaLj = std::abs( gjNorm - ljNorm + deltaLj * deltaLj);
462 std::cout <<
"gjNorm = " << gjNorm << std::endl;
463 std::cout <<
"ljNorm (exact) = " << ljnrm << std::endl;
464 std::cout <<
"ljNorm = " << ljNorm << std::endl;
465 std::cout <<
"deltaLj = " << std::sqrt(deltaLj) << std::endl;
466 std::cout <<
"kappa = " << kappa << std::endl;
469 if( std::sqrt(deltaLj) > tol*kappa*std::max(1.0,gjNorm) ){
471 std::cout <<
"Lj has become lin. INdependent" << std::endl;
473 unsigned newind = currSize_-1;
479 for (
unsigned ii=0;ii<newind-1;ii++){
480 lj[ii] =
L(newind,ii);
481 ljTz1 += lj[ii]*z1[ii];
482 ljTz2 += lj[ii]*z2[ii];
484 deltaLj = std::sqrt(deltaLj);
487 deltaLh =
GiTGj(base[currSize_-2],base[currSize_-1]);
488 for (
unsigned ii=0;ii<currSize_-1;ii++){
489 deltaLh -=
L(currSize_-2,ii)*
L(currSize_-1,ii);
494 deltaLh = - std::sqrt( deltaLh );
496 deltaLh = std::sqrt ( deltaLh );
499 L(currSize_-1,currSize_-2) =
deltaLh;
509 this->
gx_->zero(); this->
eG_->zero();
510 for (
unsigned i = 0; i < this->
size(); i++) {
512 this->
tG_->set(*this->
gx_);
514 this->
gx_->set(*this->
tG_); this->
gx_->plus(*this->
yG_);
515 this->
eG_->set(*this->
tG_); this->
eG_->axpy(-1.0,*this->
gx_); this->
eG_->plus(*this->
yG_);
517 Real Hx = 0.0, val = 0.0, err = 0.0, tmp = 0.0, y = 0.0;
518 for (
unsigned i = 0; i < this->
size(); i++) {
523 y = x[i]*(0.5*Hx + this->
alpha(i)/t) + err;
525 err = (tmp - val) + y;
527 g[i] = Hx + this->
alpha(i)/t;
532 unsigned solveDual_dim1(
const Real t,
const unsigned maxit = 1000,
const Real tol = 1.e-8) {
537 unsigned solveDual_dim2(
const Real t,
const unsigned maxit = 1000,
const Real tol = 1.e-8) {
541 Real diffg = this->
gx_->dot(*this->
gx_);
543 Real diffa = (this->
alpha(0)-this->
alpha(1))/t;
545 this->
setDualVariables(0,std::min(1.0,std::max(0.0,-(gdiffg+diffa)/diffg)));
565 void solveSystem(
int size,
char tran, Teuchos::SerialDenseMatrix<int,Real> &L, Teuchos::SerialDenseVector<int,Real> &v){
567 if( L.numRows()!=
size )
568 std::cout <<
"Error: Wrong size matrix!" << std::endl;
569 else if( v.numRows()!=
size )
570 std::cout <<
"Error: Wrong size vector!" << std::endl;
575 lapack_.TRTRS(
'L', tran,
'N', size, 1, L.values(), L.stride(), v.values(), v.stride(), &info );
580 bool isFeasible(Teuchos::SerialDenseVector<int,Real> &v,
const Real &tol){
582 for(
int i=0;i<v.numRows();i++){
590 unsigned solveDual_TT(
const Real t,
const unsigned maxit = 1000,
const Real tol = 1.e-8) {
592 std::cout <<
"Calling solveDual_TT" << std::endl;
593 std::cout <<
"t = " << t << std::endl;
594 std::cout <<
"maxit = " << maxit << std::endl;
595 std::cout <<
"tol = " << tol << std::endl;
608 L(0,0) = std::sqrt(
GiTGj(0,0));
616 z1.resize(1); z2.resize(1);
618 z2[0] = this->
alpha(0)/
L(0,0);
627 for (iter=0;iter<maxit;iter++){
631 switch( dependent_ ){
640 rho_ = ( 1 + z1Tz2/t )/z1Tz1;
641 tempv =
z1; tempv.scale(rho_);
642 tempw1 =
z2; tempw1.scale(1/t);
647 std::cout <<
"In case 0" << std::endl;
648 std::cout <<
"rho_ = " << rho_ << std::endl;
649 std::cout <<
"Solution tempv = \n" << tempv << std::endl;
659 Teuchos::SerialDenseMatrix<int,Real> LBprime( Teuchos::Copy,L,currSize_-1,currSize_-1);
660 lh.size(currSize_-1);
663 for(
unsigned i=0;i<currSize_-1;i++){
664 Real tmp =
L(currSize_-1,i);
670 bool singular =
false;
673 if (std::abs(lhTz1-1.0) <= tol*kappa){
680 tempv.resize(currSize_);
681 tempv[currSize_-1] = -1.0;
686 rho_ = ( (this->
alpha(base[currSize_-1]) -
lhTz2)/t ) / ( 1.0 -
lhTz1 );
689 Real tmp = ( 1.0 + z1Tz2 / t - rho_ *
z1Tz1 )/( 1.0 - lhTz1 );
690 tempw1 =
z1; tempw1.scale(rho_);
691 tempw2 =
z2; tempw2.scale(1.0/t);
693 tempw2 =
lh; tempw2.scale(tmp);
697 tempv.resize(currSize_);
698 tempv[currSize_-1] = tmp;
702 std::cout <<
"In case 1" << std::endl;
704 std::cout <<
"rho_ = " << rho_ << std::endl;
705 std::cout <<
"Solution tempv = \n" << tempv << std::endl;
708 std::cout <<
"Direction tempv = \n" << tempv << std::endl;
718 Teuchos::SerialDenseMatrix<int,Real> LBprime( Teuchos::Copy,L,currSize_-2,currSize_-2 );
719 lj.size(currSize_-2);
720 lh.size(currSize_-2);
723 for(
unsigned i=0;i<currSize_-2;i++){
724 Real tmp1 =
L(currSize_-1,i);
725 Real tmp2 =
L(currSize_-2,i);
731 if(std::abs(ljTz1-1.0) <= tol*kappa){
735 tempv.resize(currSize_);
736 tempv[currSize_-2] = 0.0;
737 tempv[currSize_-1] = -1.0;
741 Real mu = ( 1.0 -
lhTz1 )/( 1.0 - ljTz1 );
742 tempw1 =
lj; tempw1.scale(-mu);
746 tempv.resize(currSize_);
747 tempv[currSize_-2] = -1.0;
748 tempv[currSize_-1] = mu;
752 std::cout <<
"In case 2" << std::endl;
753 std::cout <<
"Direction tempv = \n" << tempv << std::endl;
763 std::cout <<
"Solution tempv is feasible" << std::endl;
773 std::cout <<
"Solution tempv is NOT feasible" << std::endl;
783 if ( ( entering_ < maxind ) && ( this->
getDualVariables(entering_) == 0.0 ) ){
784 if ( tempv[currSize_-1] < 0 )
790 sp += tempv[kk]*this->
alpha(base[kk]);
801 if ( (tempv[h] < -myeps) && (-this->
getDualVariables(base[h])/tempv[h] < step) )
805 std::cout <<
"h = " << h <<
" tempv[h] = " << tempv[h] <<
" dualV[base[h]] = " << this->
getDualVariables(base[h]) << std::endl;
810 std::cout <<
"Taking step of size " << step << std::endl;
813 if (step <= 0 || step == INF){
815 std::cout <<
"Invalid step!" << std::endl;
829 for (
unsigned baseitem=0;baseitem<
currSize_;baseitem++){
835 if( base[baseitem] == entering_ ){
837 std::cout <<
"Blocking " << entering_ <<
" because it just entered" << std::endl;
839 taboo_.push_back(entering_);
852 std::cout <<
"Returning because nothing deleted and not optimal" << std::endl;
867 newobjval = - Quad/2;
871 newobjval = (rho_ + Lin/t)/2;
875 std::cout <<
"New Obj value = " << newobjval << std::endl;
881 if( ( entering_ < maxind ) && ( objval <
INF ) ){
882 if( newobjval >= objval - std::max( tol*std::abs(objval), 1.e-16 ) ){
884 std::cout <<
"Blocking " << entering_ <<
" because it didn't provide decrease" << std::endl;
886 taboo_.push_back(entering_);
894 if ( objval + std::max( tol*std::abs(objval), 1.e-16 ) <= minobjval ){
897 std::cout <<
"Taboo list has been reset." << std::endl;
903 if ( (rho_ >= -INF) && (objval <= -INF) )
907 Real minro = - std::max( tol*currSize_*std::abs(objval),1.e-16 );
908 #if ( ! FIRST_VIOLATED ) 909 Real diff = -
INF, olddiff = -
INF;
912 for (
unsigned bundleitem=0;bundleitem<this->
size();bundleitem++){
914 if ( std::find(taboo_.begin(),taboo_.end(),bundleitem) != taboo_.end() ){
916 std::cout <<
"Item " << bundleitem <<
" is blocked." << std::endl;
921 if ( std::find(base.begin(),base.end(),bundleitem) == base.end() ){
924 std::cout <<
"Base does not contain index " << bundleitem << std::endl;
930 ro += this->
alpha(bundleitem)/t;
933 std::cout <<
"RO = " << ro << std::endl;
946 #if ( FIRST_VIOLATED ) 947 entering_ = bundleitem;
951 if ( diff > olddiff ){
952 entering_ = bundleitem;
962 std::cout <<
"entering_ = " << entering_ << std::endl;
966 if (entering_ < maxind){
968 std::cout <<
"Inserting " << entering_ << std::endl;
972 base.push_back(entering_);
974 std::cout <<
"Current base = " << std::endl;
975 for (
unsigned k=0;k<base.size();k++){
976 std::cout << base[k] << std::endl;
978 std::cout <<
"dependent_ = " << dependent_ << std::endl;
986 for (
unsigned i=0;i<zsize;i++){
987 lh[i] =
GiTGj(entering_,base[i]);
989 Teuchos::SerialDenseMatrix<int,Real> LBprime( Teuchos::Copy,L,zsize,zsize);
991 for (
unsigned i=0;i<zsize;i++){
992 lhTz1 += lh[i]*z1[i];
993 lhTz2 += lh[i]*z2[i];
996 Real nrm = lh.dot(lh);
997 Real delta =
GiTGj(entering_,entering_) - nrm;
999 std::cout <<
"GiTGj = " <<
GiTGj(entering_,entering_) << std::endl;
1000 std::cout <<
"lh_dot_lh = " << nrm << std::endl;
1001 std::cout <<
"delta = " << delta << std::endl;
1013 L.reshape(currSize_,currSize_);
1015 for (
unsigned i=0;i<zsize-1;i++){
1016 L(currSize_-1,i) = lh[i];
1019 Real deltaeps = tol*kappa*std::max(1.0,lh.dot(lh));
1021 std::cout <<
"kappa = " << kappa << std::endl;
1022 std::cout <<
"deltaeps = " << deltaeps << std::endl;
1024 if ( delta > deltaeps ){
1026 unsigned ind = currSize_-1;
1027 delta = std::sqrt(delta);
1030 else if(delta < -deltaeps){
1032 std::cout <<
"NEGATIVE delta!" << std::endl;
1043 std::cout <<
"Current L = " << std::endl;
1044 std::cout << L << std::endl;
1045 std::cout <<
"Current z1 = " << std::endl;
1046 std::cout << z1 << std::endl;
1047 std::cout <<
"Current z2 = " << std::endl;
1048 std::cout << z2 << std::endl;
1053 if( objval - std::max( tol*std::abs( objval ), 1.e-16 ) > minobjval ){
1055 std::cout <<
"Returning because no dual constraint violated and f cannot reach min value " << minobjval << std::endl;
1072 unsigned solveDual(
const Real t,
const unsigned maxit = 1000,
const Real tol = 1.e-8) {
1074 if (this->
size() == 1) {
1077 else if (this->
size() == 2) {
1082 unsigned outermaxit = 20;
1083 bool increase =
false, decrease =
false;
1085 for (
unsigned it=0; it < outermaxit; it++ ){
1087 if ( QPStatus_ == 1 )
1089 else if ( QPStatus_ == -2 || QPStatus_ == -3 ){
1097 if ( (mytol > 1.e-4) || (mytol < 1.e-16) ){
1100 if ( increase && decrease ){
1105 std::cout <<
"SolveDual returned after " << iter <<
" iterations with status " << QPStatus_ <<
" and solution" << std::endl;
1106 std::vector<Real> sol;
1107 for(
unsigned i=0;i<this->
size();i++){
1109 std::cout <<
"x[" << i <<
"] = " << sol[i] << std::endl;
1111 std::cout << std::endl;
1112 std::vector<Real> g(this->
size(),0.0);
1114 std::cout <<
"and objective value = " << val << std::endl;
1115 std::cout <<
"Checking constraints" << std::endl;
1117 std::cout << success << std::endl;
1123 bool checkPrimary(std::vector<Real> &g,
const std::vector<Real> &x,
const Real t)
const {
1124 bool success =
true;
1125 this->
gx_->zero(); this->
eG_->zero();
1126 for (
unsigned i = 0; i < this->
size(); i++) {
1128 this->
tG_->set(*this->
gx_);
1130 this->
gx_->set(*this->
tG_); this->
gx_->plus(*this->
yG_);
1131 this->
eG_->set(*this->
tG_); this->
eG_->axpy(-1.0,*this->
gx_); this->
eG_->plus(*this->
yG_);
1133 Real Hx = 0.0, v = 0.0, err = 0.0, tmp = 0.0, y = 0.0;
1134 for (
unsigned i = 0; i < this->
size(); i++) {
1139 y = x[i]*(t*Hx + this->
alpha(i)) + err;
1141 err = (tmp - v) + y;
1143 g[i] = - t*Hx - this->
alpha(i);
1147 for (
unsigned i = 0; i < this->
size(); i++) {
1148 if ( g[i] > v + myeps ){
1149 std::cout <<
"Constraint " << i <<
" is violated!: g[" << i <<
"] = " << g[i] <<
", v = " << v << std::endl;
1152 else if ( g[i] < v - myeps ){
1153 std::cout <<
"Constraint " << i <<
" is inactive" << std::endl;
1156 std::cout <<
"Constraint " << i <<
" is active" << std::endl;
unsigned solveDual_dim1(const Real t, const unsigned maxit=1000, const Real tol=1.e-8)
Real GiTGj(const int i, const int j)
unsigned solveDual_TT(const Real t, const unsigned maxit=1000, const Real tol=1.e-8)
void resetDualVariables(void)
unsigned solveDual(const Real t, const unsigned maxit=1000, const Real tol=1.e-8)
unsigned solveDual_dim2(const Real t, const unsigned maxit=1000, const Real tol=1.e-8)
Teuchos::SerialDenseVector< int, Real > tempw2
Teuchos::RCP< Vector< Real > > yG_
Teuchos::SerialDenseVector< int, Real > tempv
Contains definitions of custom data types in ROL.
Teuchos::SerialDenseMatrix< int, Real > Id
bool isFeasible(Teuchos::SerialDenseVector< int, Real > &v, const Real &tol)
Teuchos::LAPACK< int, Real > lapack_
Teuchos::SerialDenseVector< int, Real > tempw1
const Vector< Real > & subgradient(const unsigned i) const
Teuchos::SerialDenseVector< int, Real > z1
Teuchos::SerialDenseVector< int, Real > lj
const Real alpha(const unsigned i) const
bool checkPrimary(std::vector< Real > &g, const std::vector< Real > &x, const Real t) const
Bundle_TT(const unsigned maxSize=10, const Real coeff=0.0, const unsigned remSize=2)
Real getDualVariables(const unsigned i)
Teuchos::SerialDenseMatrix< int, Real > L
Real evaluateObjective(std::vector< Real > &g, const std::vector< Real > &x, const Real t) const
Teuchos::RCP< Vector< Real > > eG_
unsigned size(void) const
void solveSystem(int size, char tran, Teuchos::SerialDenseMatrix< int, Real > &L, Teuchos::SerialDenseVector< int, Real > &v)
Teuchos::RCP< Vector< Real > > tG_
void deleteSubgradFromBase(unsigned ind, Real tol)
Provides the interface for and implements a bundle. The semidefinite quadratic subproblem is solved u...
void setDualVariables(const unsigned i, const Real val)
void addSubgradToBase(unsigned ind, Real delta)
std::vector< int > taboo_
Teuchos::SerialDenseVector< int, Real > lh
void swapRowsL(unsigned ind1, unsigned ind2, bool trans=false)
Teuchos::SerialDenseVector< int, Real > z2
Teuchos::RCP< Vector< Real > > gx_
Provides the interface for and implments a bundle.
static const double ROL_EPSILON
Platform-dependent machine epsilon.