ROL
ROL_PartitionedVector.hpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // Rapid Optimization Library (ROL) Package
4 //
5 // Copyright 2014 NTESS and the ROL contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
10 #include "ROL_Vector.hpp"
11 
12 #include <initializer_list>
13 
14 #ifndef ROL_PARTITIONED_VECTOR_H
15 #define ROL_PARTITIONED_VECTOR_H
16 
23 namespace ROL {
24 
25 template<class Real>
26 class PartitionedVector : public Vector<Real> {
27 
28  typedef Vector<Real> V;
29  typedef ROL::Ptr<V> Vp;
31 
32 private:
33  const std::vector<Vp> vecs_;
34  mutable std::vector<Vp> dual_vecs_;
35  mutable ROL::Ptr<PV> dual_pvec_;
36 public:
37 
39 
40  PartitionedVector( const std::vector<Vp> &vecs ) : vecs_(vecs) {
41  for( size_type i=0; i<vecs_.size(); ++i ) {
42  dual_vecs_.push_back((vecs_[i]->dual()).clone());
43  }
44  }
45 
46  void set( const V &x ) {
47  const PV &xs = dynamic_cast<const PV&>(x);
48  ROL_TEST_FOR_EXCEPTION( numVectors() != xs.numVectors(),
49  std::invalid_argument,
50  "Error: Vectors must have the same number of subvectors." );
51  for( size_type i=0; i<vecs_.size(); ++i ) {
52  vecs_[i]->set(*xs.get(i));
53  }
54  }
55 
56  void plus( const V &x ) {
57  const PV &xs = dynamic_cast<const PV&>(x);
58  ROL_TEST_FOR_EXCEPTION( numVectors() != xs.numVectors(),
59  std::invalid_argument,
60  "Error: Vectors must have the same number of subvectors." );
61  for( size_type i=0; i<vecs_.size(); ++i ) {
62  vecs_[i]->plus(*xs.get(i));
63  }
64  }
65 
66  void scale( const Real alpha ) {
67  for( size_type i=0; i<vecs_.size(); ++i ) {
68  vecs_[i]->scale(alpha);
69  }
70  }
71 
72  void axpy( const Real alpha, const V &x ) {
73  const PV &xs = dynamic_cast<const PV&>(x);
74  ROL_TEST_FOR_EXCEPTION( numVectors() != xs.numVectors(),
75  std::invalid_argument,
76  "Error: Vectors must have the same number of subvectors." );
77 
78  for( size_type i=0; i<vecs_.size(); ++i ) {
79  vecs_[i]->axpy(alpha,*xs.get(i));
80  }
81  }
82 
83  Real dot( const V &x ) const {
84  const PV &xs = dynamic_cast<const PV&>(x);
85  ROL_TEST_FOR_EXCEPTION( numVectors() != xs.numVectors(),
86  std::invalid_argument,
87  "Error: Vectors must have the same number of subvectors." );
88  Real result = 0;
89  for( size_type i=0; i<vecs_.size(); ++i ) {
90  result += vecs_[i]->dot(*xs.get(i));
91  }
92  return result;
93  }
94 
95  Real norm() const {
96  Real result = 0;
97  for( size_type i=0; i<vecs_.size(); ++i ) {
98  result += std::pow(vecs_[i]->norm(),2);
99  }
100  return std::sqrt(result);
101  }
102 
103  Vp clone() const {
104  std::vector<Vp> clonevec;
105  for( size_type i=0; i<vecs_.size(); ++i ) {
106  clonevec.push_back(vecs_[i]->clone());
107  }
108  return ROL::makePtr<PV>(clonevec);
109  }
110 
111  const V& dual(void) const {
112  for( size_type i=0; i<vecs_.size(); ++i ) {
113  dual_vecs_[i]->set(vecs_[i]->dual());
114  }
115  dual_pvec_ = ROL::makePtr<PV>( dual_vecs_ );
116  return *dual_pvec_;
117  }
118 
119  Real apply(const Vector<Real> &x) const {
120  const PV &xs = dynamic_cast<const PV&>(x);
121  ROL_TEST_FOR_EXCEPTION( numVectors() != xs.numVectors(),
122  std::invalid_argument,
123  "Error: Vectors must have the same number of subvectors." );
124  Real result = 0;
125  for( size_type i=0; i<vecs_.size(); ++i ) {
126  result += vecs_[i]->apply(*xs.get(i));
127  }
128  return result;
129  }
130 
131  Vp basis( const int i ) const {
132  ROL_TEST_FOR_EXCEPTION( i >= dimension() || i<0,
133  std::invalid_argument,
134  "Error: Basis index must be between 0 and vector dimension." );
135  Vp bvec = clone();
136  // Downcast
137  PV &eb = dynamic_cast<PV&>(*bvec);
138 
139  // Iterate over subvectors
140  int begin = 0, end = 0;
141  for( size_type j=0; j<vecs_.size(); ++j ) {
142  end += vecs_[j]->dimension();
143  if( begin<= i && i<end ) {
144  eb.set(j, *(vecs_[j]->basis(i-begin)) );
145  }
146  else {
147  eb.zero(j);
148  }
149  begin = end;
150  }
151  return bvec;
152  }
153 
154  int dimension() const {
155  int total_dim = 0;
156  for( size_type j=0; j<vecs_.size(); ++j ) {
157  total_dim += vecs_[j]->dimension();
158  }
159  return total_dim;
160  }
161 
162  void zero() {
163  for( size_type j=0; j<vecs_.size(); ++j ) {
164  vecs_[j]->zero();
165  }
166  }
167 
168  // Apply the same unary function to each subvector
169  void applyUnary( const Elementwise::UnaryFunction<Real> &f ) {
170  for( size_type i=0; i<vecs_.size(); ++i ) {
171  vecs_[i]->applyUnary(f);
172  }
173  }
174 
175  // Apply the same binary function to each pair of subvectors in this vector and x
176  void applyBinary( const Elementwise::BinaryFunction<Real> &f, const V &x ) {
177  const PV &xs = dynamic_cast<const PV&>(x);
178 
179  for( size_type i=0; i<vecs_.size(); ++i ) {
180  vecs_[i]->applyBinary(f,*xs.get(i));
181  }
182  }
183 
184  Real reduce( const Elementwise::ReductionOp<Real> &r ) const {
185  Real result = r.initialValue();
186 
187  for( size_type i=0; i<vecs_.size(); ++i ) {
188  r.reduce(vecs_[i]->reduce(r),result);
189  }
190  return result;
191  }
192 
193  void setScalar( const Real C ) {
194  for (size_type i=0; i<vecs_.size(); ++i) {
195  vecs_[i]->setScalar(C);
196  }
197  }
198 
199  void randomize( const Real l = 0.0, const Real u = 1.0 ) {
200  for (size_type i=0; i<vecs_.size(); ++i) {
201  vecs_[i]->randomize(l,u);
202  }
203  }
204 
205  void print( std::ostream &outStream ) const {
206  for( size_type i=0; i<vecs_.size(); ++i ) {
207  outStream << "V[" << i << "]: ";
208  vecs_[i]->print(outStream);
209  }
210  }
211 
212  // Methods that do not exist in the base class
213  const Vector<Real>& operator[]( size_type i ) const {
214  return *(vecs_[i]);
215  }
216 
218  return *(vecs_[i]);
219  }
220 
221  ROL::Ptr<const Vector<Real> > get(size_type i) const {
222  return vecs_[i];
223  }
224 
225  ROL::Ptr<Vector<Real> > get(size_type i) {
226  return vecs_[i];
227  }
228 
229  void set(size_type i, const V &x) {
230  vecs_[i]->set(x);
231  }
232 
233  void zero(size_type i) {
234  vecs_[i]->zero();
235  }
236 
238  return vecs_.size();
239  }
240 
241 public:
242 
243  // Make a new PartitionedVector from an initializer_list of pointers to vectors
244  static Ptr<PartitionedVector> create( std::initializer_list<Vp> vs ) {
245  std::vector<Vp> subvecs{vs};
246  return ROL::makePtr<PartitionedVector>( subvecs );
247  }
248 
249  // Make a new PartitionedVector by cloning the given vector N times
250  static Ptr<PartitionedVector> create( const V& x, size_type N ) {
251  std::vector<Vp> subvecs(N);
252  for( size_type i=0; i<N; ++i ) subvecs.at(i) = x.clone();
253  return ROL::makePtr<PartitionedVector>( subvecs );
254  }
255 
256 };
257 
258 // Helper methods
259 
260 template<typename Real>
262  return static_cast<PartitionedVector<Real>&>(x);
263 }
264 
265 template<typename Real>
266 inline const PartitionedVector<Real>& partition( const Vector<Real>& x ) {
267  return static_cast<const PartitionedVector<Real>&>(x);
268 }
269 
270 template<typename Real>
271 inline Ptr<PartitionedVector<Real>> partition( const Ptr<Vector<Real>>& x ) {
272  return staticPtrCast<PartitionedVector<Real>>(x);
273 }
274 
275 template<typename Real>
276 inline Ptr<const PartitionedVector<Real>> partition( const Ptr<const Vector<Real>>& x ) {
277  return staticPtrCast<const PartitionedVector<Real>>(x);
278 }
279 
280 
281 
282 template<class Real>
283 ROL::Ptr<Vector<Real>>
284 CreatePartitionedVector( const ROL::Ptr<Vector<Real>> &a ) {
285 
286  using Vp = ROL::Ptr<Vector<Real>>;
287  using PV = PartitionedVector<Real>;
288 
289  Vp temp[] = {a};
290  return ROL::makePtr<PV>( std::vector<Vp>(temp, temp+1) );
291 }
292 
293 template<class Real>
294 ROL::Ptr<const Vector<Real> >
295 CreatePartitionedVector( const ROL::Ptr<const Vector<Real>> &a ) {
296 
297  using Vp = ROL::Ptr<const Vector<Real>>;
298  using PV = const PartitionedVector<Real>;
299 
300  Vp temp[] = {a};
301  return ROL::makePtr<PV>( std::vector<Vp>(temp, temp+1) );
302 }
303 
304 template<class Real>
305 ROL::Ptr<Vector<Real>>
307  const ROL::Ptr<Vector<Real>> &b ) {
308  using Vp = ROL::Ptr<Vector<Real>>;
309  using PV = PartitionedVector<Real>;
310 
311  Vp temp[] = {a,b};
312  return ROL::makePtr<PV>( std::vector<Vp>(temp, temp+2) );
313 }
314 
315 template<class Real>
316 ROL::Ptr<const Vector<Real> >
317 CreatePartitionedVector( const ROL::Ptr<const Vector<Real>> &a,
318  const ROL::Ptr<const Vector<Real>> &b ) {
319  using Vp = ROL::Ptr<const Vector<Real>>;
320  using PV = const PartitionedVector<Real>;
321 
322  Vp temp[] = {a,b};
323  return ROL::makePtr<PV>( std::vector<Vp>(temp, temp+2) );
324 }
325 
326 
327 template<class Real>
328 ROL::Ptr<Vector<Real>>
330  const ROL::Ptr<Vector<Real>> &b,
331  const ROL::Ptr<Vector<Real>> &c ) {
332 
333  using Vp = ROL::Ptr<Vector<Real>>;
334  using PV = PartitionedVector<Real>;
335 
336  Vp temp[] = {a,b,c};
337  return ROL::makePtr<PV>( std::vector<Vp>(temp, temp+3) );
338 }
339 
340 template<class Real>
341 ROL::Ptr<const Vector<Real> >
342 CreatePartitionedVector( const ROL::Ptr<const Vector<Real>> &a,
343  const ROL::Ptr<const Vector<Real>> &b,
344  const ROL::Ptr<const Vector<Real>> &c ) {
345 
346  using Vp = ROL::Ptr<const Vector<Real>>;
347  using PV = const PartitionedVector<Real>;
348 
349  Vp temp[] = {a,b,c};
350  return ROL::makePtr<PV>( std::vector<Vp>(temp, temp+3) );
351 }
352 
353 template<class Real>
354 ROL::Ptr<Vector<Real> >
356  const ROL::Ptr<Vector<Real>> &b,
357  const ROL::Ptr<Vector<Real>> &c,
358  const ROL::Ptr<Vector<Real>> &d ) {
359 
360 
361  typedef ROL::Ptr<Vector<Real> > Vp;
362  typedef PartitionedVector<Real> PV;
363 
364  Vp temp[] = {a,b,c,d};
365  return ROL::makePtr<PV>( std::vector<Vp>(temp, temp+4) );
366 }
367 
368 template<class Real>
369 ROL::Ptr<const Vector<Real> >
370 CreatePartitionedVector( const ROL::Ptr<const Vector<Real>> &a,
371  const ROL::Ptr<const Vector<Real>> &b,
372  const ROL::Ptr<const Vector<Real>> &c,
373  const ROL::Ptr<const Vector<Real>> &d ) {
374 
375 
376  using Vp = ROL::Ptr<const Vector<Real>>;
377  using PV = const PartitionedVector<Real>;
378 
379  Vp temp[] = {a,b,c,d};
380  return ROL::makePtr<PV>( std::vector<Vp>(temp, temp+4) );
381 }
382 
383 } // namespace ROL
384 
385 #endif // ROL_PARTITIONED_VECTOR_H
386 
PartitionedVector< Real > & partition(Vector< Real > &x)
PartitionedVector(const std::vector< Vp > &vecs)
static Ptr< PartitionedVector > create(std::initializer_list< Vp > vs)
PartitionedVector< Real > PV
typename PV< Real >::size_type size_type
virtual ROL::Ptr< Vector > clone() const =0
Clone to make a new (uninitialized) vector.
ROL::Ptr< const Vector< Real > > get(size_type i) const
Real norm() const
Returns where .
Defines the linear algebra of vector space on a generic partitioned vector.
ROL::Ptr< Vector< Real > > CreatePartitionedVector(const ROL::Ptr< Vector< Real >> &a)
const V & dual(void) const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
Vector< Real > & operator[](size_type i)
void applyBinary(const Elementwise::BinaryFunction< Real > &f, const V &x)
void print(std::ostream &outStream) const
static Ptr< PartitionedVector > create(const V &x, size_type N)
Vp clone() const
Clone to make a new (uninitialized) vector.
Defines the linear algebra or vector space interface.
Definition: ROL_Vector.hpp:46
Real dot(const V &x) const
Compute where .
Real apply(const Vector< Real > &x) const
Apply to a dual vector. This is equivalent to the call .
void zero()
Set to zero vector.
void applyUnary(const Elementwise::UnaryFunction< Real > &f)
void setScalar(const Real C)
Set where .
PartitionedVector< Real > PV
void set(size_type i, const V &x)
void scale(const Real alpha)
Compute where .
int dimension() const
Return dimension of the vector space.
void set(const V &x)
Set where .
Vp basis(const int i) const
Return i-th basis vector.
const Vector< Real > & operator[](size_type i) const
std::vector< PV >::size_type size_type
void randomize(const Real l=0.0, const Real u=1.0)
Set vector to be uniform random between [l,u].
const std::vector< Vp > vecs_
void plus(const V &x)
Compute , where .
Real reduce(const Elementwise::ReductionOp< Real > &r) const
void axpy(const Real alpha, const V &x)
Compute where .