44 #ifndef ROL_BUNDLE_U_DEF_H
45 #define ROL_BUNDLE_U_DEF_H
51 template<
typename Real>
54 for (
unsigned j = ind.back()+1; j < size_; ++j) {
55 (subgradients_[j-1])->set(*(subgradients_[j]));
56 linearizationErrors_[j-1] = linearizationErrors_[j];
57 distanceMeasures_[j-1] = distanceMeasures_[j];
58 dualVariables_[j-1] = dualVariables_[j];
60 (subgradients_[size_-1])->
zero();
61 linearizationErrors_[size_-1] = ROL_OVERFLOW<Real>();
62 distanceMeasures_[size_-1] = ROL_OVERFLOW<Real>();
63 dualVariables_[size_-1] =
zero;
64 for (
unsigned i = ind.size()-1; i > 0; --i) {
65 for (
unsigned j = ind[i-1]+1; j < size_; ++j) {
66 (subgradients_[j-1])->set(*(subgradients_[j]));
67 linearizationErrors_[j-1] = linearizationErrors_[j];
68 distanceMeasures_[j-1] = distanceMeasures_[j];
69 dualVariables_[j-1] = dualVariables_[j];
75 template<
typename Real>
78 (subgradients_[size_])->set(g);
79 linearizationErrors_[size_] = le;
80 distanceMeasures_[size_] = dm;
81 dualVariables_[size_] =
zero;
85 template<
typename Real>
89 const unsigned remSize)
90 : size_(0), maxSize_(maxSize), remSize_(remSize),
91 coeff_(coeff), omega_(omega), isInitialized_(false) {
106 template<
typename Real>
108 if ( !isInitialized_ ) {
109 Real
zero(0), one(1);
110 for (
unsigned i = 0; i < maxSize_; ++i) {
111 subgradients_[i] = g.
clone();
113 (subgradients_[0])->set(g);
114 linearizationErrors_[0] =
zero;
115 distanceMeasures_[0] =
zero;
116 dualVariables_[0] = one;
118 isInitialized_ =
true;
127 template<
typename Real>
129 return linearizationErrors_[i];
132 template<
typename Real>
134 return distanceMeasures_[i];
137 template<
typename Real>
139 return *(subgradients_[i]);
142 template<
typename Real>
144 return dualVariables_[i];
147 template<
typename Real>
149 dualVariables_[i] = val;
152 template<
typename Real>
155 dualVariables_.assign(size_,zero);
158 template<
typename Real>
161 if ( coeff_ > ROL_EPSILON<Real>() ) {
162 alpha = std::max(coeff_*std::pow(dm,omega_),le);
167 template<
typename Real>
169 return computeAlpha(distanceMeasures_[i],linearizationErrors_[i]);
172 template<
typename Real>
177 template<
typename Real>
179 Real
zero(0), one(1);
181 Real eLE(0), eDM(0), yLE(0), yDM(0), tLE(0), tDM(0);
182 for (
unsigned i = 0; i < size_; ++i) {
185 yG_->set(*subgradients_[i]); yG_->scale(dualVariables_[i]); yG_->axpy(-one,*eG_);
186 tG_->set(aggSubGrad); tG_->plus(*yG_);
187 eG_->set(*tG_); eG_->axpy(-one,aggSubGrad); eG_->axpy(-one,*yG_);
188 aggSubGrad.
set(*tG_);
191 yLE = dualVariables_[i]*linearizationErrors_[i] - eLE;
192 tLE = aggLinErr + yLE;
193 eLE = (tLE - aggLinErr) - yLE;
197 yDM = dualVariables_[i]*distanceMeasures_[i] - eDM;
198 tDM = aggDistMeas + yDM;
199 eDM = (tDM - aggDistMeas) - yDM;
204 template<
typename Real>
206 if (size_ == maxSize_) {
208 unsigned loc = size_, cnt = 0;
209 std::vector<unsigned> ind(remSize_,0);
210 for (
unsigned i = size_; i > 0; --i) {
211 if ( std::abs(linearizationErrors_[i-1]) < ROL_EPSILON<Real>() ) {
216 for (
unsigned i = 0; i < size_; ++i) {
221 if (cnt == remSize_) {
232 template<
typename Real>
238 for (
unsigned i = 0; i < size_; ++i) {
240 linearizationErrors_[i] += linErr - subgradients_[i]->apply(s);
241 distanceMeasures_[i] += distMeas;
243 linearizationErrors_[size_] =
zero;
244 distanceMeasures_[size_] =
zero;
248 linearizationErrors_[size_] = linErr;
249 distanceMeasures_[size_] = distMeas;
252 (subgradients_[size_])->set(g);
254 dualVariables_[size_] =
zero;
259 template<
typename Real>
261 return subgradient(i).dot(subgradient(j));
264 template<
typename Real>
266 return x.
dot(subgradient(i));
269 template<
typename Real>
271 x.
axpy(a,subgradient(i));
274 template<
typename Real>
276 Real one(1), half(0.5);
277 gx_->zero(); eG_->zero();
278 for (
unsigned i = 0; i < size(); ++i) {
281 yG_->set(subgradient(i)); yG_->scale(x[i]); yG_->axpy(-one,*eG_);
282 tG_->set(*gx_); tG_->plus(*yG_);
283 eG_->set(*tG_); eG_->axpy(-one,*gx_); eG_->axpy(-one,*yG_);
286 Real Hx(0), val(0), err(0), tmp(0), y(0);
287 for (
unsigned i = 0; i < size(); ++i) {
292 y = x[i]*(half*Hx + alpha(i)/t) - err;
294 err = (tmp - val) - y;
297 g[i] = Hx + alpha(i)/t;
302 template<
typename Real>
304 setDualVariable(0,static_cast<Real>(1));
309 template<
typename Real>
311 Real diffg = gx_->dot(*gx_),
zero(0), one(1), half(0.5);
312 gx_->set(subgradient(0)); addGi(1,-one,*gx_);
313 if ( std::abs(diffg) > ROL_EPSILON<Real>() ) {
314 Real diffa = (alpha(0)-alpha(1))/t;
315 Real gdiffg = dotGi(1,*gx_);
316 setDualVariable(0,std::min(one,std::max(
zero,-(gdiffg+diffa)/diffg)));
317 setDualVariable(1,one-getDualVariable(0));
320 if ( std::abs(alpha(0)-alpha(1)) > ROL_EPSILON<Real>() ) {
321 if ( alpha(0) < alpha(1) ) {
322 setDualVariable(0,one); setDualVariable(1,
zero);
324 else if ( alpha(0) > alpha(1) ) {
325 setDualVariable(0,
zero); setDualVariable(1,one);
329 setDualVariable(0,half); setDualVariable(1,half);
void reset(const Vector< Real > &g, const Real le, const Real dm)
void remove(const std::vector< unsigned > &ind)
unsigned solveDual_dim1(const Real t, const unsigned maxit=1000, const Real tol=1.e-8)
virtual ROL::Ptr< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
void setDualVariable(const unsigned i, const Real val)
const Vector< Real > & subgradient(const unsigned i) const
virtual void axpy(const Real alpha, const Vector &x)
Compute where .
const Real dotGi(const unsigned i, const Vector< Real > &x) const
void update(const bool flag, const Real linErr, const Real distMeas, const Vector< Real > &g, const Vector< Real > &s)
Contains definitions of custom data types in ROL.
const Real GiGj(const unsigned i, const unsigned j) const
std::vector< Real > dualVariables_
const Real alpha(const unsigned i) const
virtual void zero()
Set to zero vector.
std::vector< Real > distanceMeasures_
Defines the linear algebra or vector space interface.
virtual void initialize(const Vector< Real > &g)
virtual Real dot(const Vector &x) const =0
Compute where .
void resetDualVariables(void)
Objective_SerialSimOpt(const Ptr< Obj > &obj, const V &ui) z0_ zero()
unsigned size(void) const
void aggregate(Vector< Real > &aggSubGrad, Real &aggLinErr, Real &aggDistMeas) const
Objective_SerialSimOpt(const Ptr< Obj > &obj, const V &ui) z0_ zero()
Bundle_U(const unsigned maxSize=10, const Real coeff=0.0, const Real omega=2.0, const unsigned remSize=2)
Real evaluateObjective(std::vector< Real > &g, const std::vector< Real > &x, const Real t) const
void add(const Vector< Real > &g, const Real le, const Real dm)
const Real distanceMeasure(const unsigned i) const
unsigned solveDual_dim2(const Real t, const unsigned maxit=1000, const Real tol=1.e-8)
std::vector< Ptr< Vector< Real > > > subgradients_
void addGi(const unsigned i, const Real a, Vector< Real > &x) const
const Real linearizationError(const unsigned i) const
std::vector< Real > linearizationErrors_
virtual void set(const Vector &x)
Set where .
const Real getDualVariable(const unsigned i) const
const Real computeAlpha(const Real dm, const Real le) const