ROL
ROL_PartitionedVector.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ************************************************************************
3 //
4 // Rapid Optimization Library (ROL) Package
5 // Copyright (2014) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact lead developers:
38 // Drew Kouri (dpkouri@sandia.gov) and
39 // Denis Ridzal (dridzal@sandia.gov)
40 //
41 // ************************************************************************
42 // @HEADER
43 
44 #include "ROL_Vector.hpp"
45 
46 #ifndef ROL_PARTITIONED_VECTOR_H
47 #define ROL_PARTITIONED_VECTOR_H
48 
55 namespace ROL {
56 
57 template<class Real>
58 class PartitionedVector : public Vector<Real> {
59 
60  typedef Vector<Real> V;
61  typedef Teuchos::RCP<V> RCPV;
63 
64 private:
65  Teuchos::RCP<std::vector<RCPV> > vecs_;
66  mutable std::vector<RCPV> dual_vecs_;
67  mutable Teuchos::RCP<PV> dual_pvec_;
68 public:
69 
70  typedef typename std::vector<PV>::size_type size_type;
71 
72  PartitionedVector( const Teuchos::RCP<std::vector<RCPV> > &vecs ) :
73  vecs_(vecs) {
74  for( size_type i=0; i<vecs_->size(); ++i ) {
75  dual_vecs_.push_back(((*vecs_)[i]->dual()).clone());
76  }
77  }
78 
79  void set( const V &x ) {
80 
81 
82  using Teuchos::dyn_cast;
83  const PV &xs = dyn_cast<const PV>(dyn_cast<const V>(x));
84 
85  TEUCHOS_TEST_FOR_EXCEPTION( numVectors() != xs.numVectors(),
86  std::invalid_argument,
87  "Error: Vectors must have the same number of subvectors." );
88 
89  for( size_type i=0; i<vecs_->size(); ++i ) {
90  (*vecs_)[i]->set(*xs.get(i));
91  }
92  }
93 
94  void plus( const V &x ) {
95  using Teuchos::dyn_cast;
96  const PV &xs = dyn_cast<const PV>(dyn_cast<const V>(x));
97 
98  TEUCHOS_TEST_FOR_EXCEPTION( numVectors() != xs.numVectors(),
99  std::invalid_argument,
100  "Error: Vectors must have the same number of subvectors." );
101 
102  for( size_type i=0; i<vecs_->size(); ++i ) {
103  (*vecs_)[i]->plus(*xs.get(i));
104  }
105  }
106 
107  void scale( const Real alpha ) {
108  for( size_type i=0; i<vecs_->size(); ++i ) {
109  (*vecs_)[i]->scale(alpha);
110  }
111  }
112 
113  void axpy( const Real alpha, const V &x ) {
114  using Teuchos::dyn_cast;
115  const PV &xs = dyn_cast<const PV>(x);
116 
117  TEUCHOS_TEST_FOR_EXCEPTION( numVectors() != xs.numVectors(),
118  std::invalid_argument,
119  "Error: Vectors must have the same number of subvectors." );
120 
121  for( size_type i=0; i<vecs_->size(); ++i ) {
122  (*vecs_)[i]->axpy(alpha,*xs.get(i));
123  }
124  }
125 
126  Real dot( const V &x ) const {
127  using Teuchos::dyn_cast;
128  const PV &xs = dyn_cast<const PV>(x);
129 
130  TEUCHOS_TEST_FOR_EXCEPTION( numVectors() != xs.numVectors(),
131  std::invalid_argument,
132  "Error: Vectors must have the same number of subvectors." );
133 
134  Real result = 0;
135  for( size_type i=0; i<vecs_->size(); ++i ) {
136  result += (*vecs_)[i]->dot(*xs.get(i));
137  }
138  return result;
139  }
140 
141  Real norm() const {
142  Real result = 0;
143  for( size_type i=0; i<vecs_->size(); ++i ) {
144  result += std::pow((*vecs_)[i]->norm(),2);
145  }
146  return std::sqrt(result);
147  }
148 
149  RCPV clone() const {
150  using Teuchos::RCP;
151  using Teuchos::rcp;
152 
153  RCP<std::vector<RCPV> > clonevec = rcp( new std::vector<RCPV> );
154 
155  for( size_type i=0; i<vecs_->size(); ++i ) {
156  clonevec->push_back((*vecs_)[i]->clone());
157  }
158  return rcp( new PV(clonevec) );
159  }
160 
161  const V& dual(void) const {
162 
163  using Teuchos::rcp;
164 
165  for( size_type i=0; i<vecs_->size(); ++i ) {
166  dual_vecs_[i]->set((*vecs_)[i]->dual());
167  }
168  dual_pvec_ = rcp( new PV( rcp( &dual_vecs_, false ) ) );
169  return *dual_pvec_;
170  }
171 
172  RCPV basis( const int i ) const {
173 
174  TEUCHOS_TEST_FOR_EXCEPTION( i >= dimension() || i<0,
175  std::invalid_argument,
176  "Error: Basis index must be between 0 and vector dimension." );
177 
178  using Teuchos::RCP;
179  using Teuchos::rcp;
180  using Teuchos::dyn_cast;
181 
182  RCPV bvec = clone();
183 
184  // Downcast
185  PV &eb = dyn_cast<PV>(*bvec);
186 
187  int begin = 0;
188  int end = 0;
189 
190  // Iterate over subvectors
191  for( size_type j=0; j<vecs_->size(); ++j ) {
192 
193  end += (*vecs_)[j]->dimension();
194 
195  if( begin<= i && i<end ) {
196  eb.set(j, *((*vecs_)[j]->basis(i-begin)) );
197  }
198  else {
199  eb.zero(j);
200  }
201 
202  begin = end;
203 
204  }
205  return bvec;
206  }
207 
208  int dimension() const {
209  int total_dim = 0;
210  for( size_type j=0; j<vecs_->size(); ++j ) {
211  total_dim += (*vecs_)[j]->dimension();
212  }
213  return total_dim;
214  }
215 
216  void zero() {
217  for( size_type j=0; j<vecs_->size(); ++j ) {
218  (*vecs_)[j]->zero();
219  }
220  }
221 
222  // Apply the same unary function to each subvector
223  void applyUnary( const Elementwise::UnaryFunction<Real> &f ) {
224  for( size_type i=0; i<vecs_->size(); ++i ) {
225  (*vecs_)[i]->applyUnary(f);
226  }
227  }
228 
229  // Apply the same binary function to each pair of subvectors in this vector and x
230  void applyBinary( const Elementwise::BinaryFunction<Real> &f, const V &x ) {
231  const PV &xs = Teuchos::dyn_cast<const PV>(x);
232 
233  for( size_type i=0; i<vecs_->size(); ++i ) {
234  (*vecs_)[i]->applyBinary(f,*xs.get(i));
235  }
236  }
237 
238  Real reduce( const Elementwise::ReductionOp<Real> &r ) const {
239 
240  Real result = r.initialValue();
241 
242  for( size_type i=0; i<vecs_->size(); ++i ) {
243  r.reduce((*vecs_)[i]->reduce(r),result);
244  }
245  return result;
246 
247  }
248 
249 
250 
251 
252 
253  // Methods that do not exist in the base class
254 
255  Teuchos::RCP<const Vector<Real> > get(size_type i) const {
256  return (*vecs_)[i];
257  }
258 
259  Teuchos::RCP<Vector<Real> > get(size_type i) {
260  return (*vecs_)[i];
261  }
262 
263  void set(size_type i, const V &x) {
264  (*vecs_)[i]->set(x);
265  }
266 
267  void zero(size_type i) {
268  (*vecs_)[i]->zero();
269  }
270 
271  size_type numVectors() const {
272  return vecs_->size();
273  }
274 
275 };
276 
277 // Helper methods
278 template<class Real>
279 Teuchos::RCP<Vector<Real> > CreatePartitionedVector( const Teuchos::RCP<Vector<Real> > &a ) {
280  using Teuchos::RCP;
281  using Teuchos::rcp;
282  typedef RCP<Vector<Real> > RCPV;
283  typedef PartitionedVector<Real> PV;
284 
285  RCPV temp[] = {a};
286  return rcp( new PV( rcp( new std::vector<RCPV>(temp, temp+1) ) ) );
287 }
288 
289 template<class Real>
290 Teuchos::RCP<const Vector<Real> > CreatePartitionedVector( const Teuchos::RCP<const Vector<Real> > &a ) {
291  using Teuchos::RCP;
292  using Teuchos::rcp;
293  typedef RCP<const Vector<Real> > RCPV;
294  typedef const PartitionedVector<Real> PV;
295 
296  RCPV temp[] = {a};
297  return rcp( new PV( rcp( new std::vector<RCPV>(temp, temp+1) ) ) );
298 }
299 
300 template<class Real>
301 Teuchos::RCP<Vector<Real> > CreatePartitionedVector( const Teuchos::RCP<Vector<Real> > &a,
302  const Teuchos::RCP<Vector<Real> > &b ) {
303  using Teuchos::RCP;
304  using Teuchos::rcp;
305  typedef RCP<Vector<Real> > RCPV;
306  typedef PartitionedVector<Real> PV;
307 
308  RCPV temp[] = {a,b};
309  return rcp( new PV( rcp( new std::vector<RCPV>(temp, temp+2) ) ) );
310 }
311 
312 template<class Real>
313 Teuchos::RCP<const Vector<Real> > CreatePartitionedVector( const Teuchos::RCP<const Vector<Real> > &a,
314  const Teuchos::RCP<const Vector<Real> > &b ) {
315  using Teuchos::RCP;
316  using Teuchos::rcp;
317  typedef RCP<const Vector<Real> > RCPV;
318  typedef const PartitionedVector<Real> PV;
319 
320  RCPV temp[] = {a,b};
321  return rcp( new PV( rcp( new std::vector<RCPV>(temp, temp+2) ) ) );
322 }
323 
324 template<class Real>
325 Teuchos::RCP<Vector<Real> > CreatePartitionedVector( const Teuchos::RCP<Vector<Real> > &a,
326  const Teuchos::RCP<Vector<Real> > &b,
327  const Teuchos::RCP<Vector<Real> > &c ) {
328  using Teuchos::RCP;
329  using Teuchos::rcp;
330  typedef RCP<Vector<Real> > RCPV;
331  typedef PartitionedVector<Real> PV;
332 
333  RCPV temp[] = {a,b,c};
334  return rcp( new PV( rcp( new std::vector<RCPV>(temp, temp+3) ) ) );
335 }
336 
337 template<class Real>
338 Teuchos::RCP<const Vector<Real> > CreatePartitionedVector( const Teuchos::RCP<const Vector<Real> > &a,
339  const Teuchos::RCP<const Vector<Real> > &b,
340  const Teuchos::RCP<const Vector<Real> > &c ) {
341  using Teuchos::RCP;
342  using Teuchos::rcp;
343  typedef RCP<const Vector<Real> > RCPV;
344  typedef const PartitionedVector<Real> PV;
345 
346  RCPV temp[] = {a,b,c};
347  return rcp( new PV( rcp( new std::vector<RCPV>(temp, temp+3) ) ) );
348 }
349 
350 template<class Real>
351 Teuchos::RCP<Vector<Real> > CreatePartitionedVector( const Teuchos::RCP<Vector<Real> > &a,
352  const Teuchos::RCP<Vector<Real> > &b,
353  const Teuchos::RCP<Vector<Real> > &c,
354  const Teuchos::RCP<Vector<Real> > &d ) {
355  using Teuchos::RCP;
356  using Teuchos::rcp;
357  typedef RCP<Vector<Real> > RCPV;
358  typedef PartitionedVector<Real> PV;
359 
360  RCPV temp[] = {a,b,c,d};
361  return rcp( new PV( rcp( new std::vector<RCPV>(temp, temp+4) ) ) );
362 }
363 
364 template<class Real>
365 Teuchos::RCP<const Vector<Real> > CreatePartitionedVector( const Teuchos::RCP<const Vector<Real> > &a,
366  const Teuchos::RCP<const Vector<Real> > &b,
367  const Teuchos::RCP<const Vector<Real> > &c,
368  const Teuchos::RCP<const Vector<Real> > &d ) {
369  using Teuchos::RCP;
370  using Teuchos::rcp;
371  typedef RCP<const Vector<Real> > RCPV;
372  typedef const PartitionedVector<Real> PV;
373 
374  RCPV temp[] = {a,b,c,d};
375  return rcp( new PV( rcp( new std::vector<RCPV>(temp, temp+4) ) ) );
376 }
377 
378 
379 
380 
381 } // namespace ROL
382 
383 #endif // ROL_PARTITIONED_VECTOR_H
384 
PartitionedVector(const Teuchos::RCP< std::vector< RCPV > > &vecs)
Real norm() const
Returns where .
Defines the linear algebra of vector space on a generic partitioned vector.
const V & dual(void) const
Return dual representation of , for example, the result of applying a Riesz map, or change of basis...
void applyBinary(const Elementwise::BinaryFunction< Real > &f, const V &x)
Teuchos::RCP< Vector< Real > > CreatePartitionedVector(const Teuchos::RCP< Vector< Real > > &a)
Defines the linear algebra or vector space interface.
Definition: ROL_Vector.hpp:74
Real dot(const V &x) const
Compute where .
void zero()
Set to zero vector.
void applyUnary(const Elementwise::UnaryFunction< Real > &f)
PartitionedVector< Real > PV
void scale(const Real alpha)
Compute where .
int dimension() const
Return dimension of the vector space.
Teuchos::RCP< const Vector< Real > > get(size_type i) const
RCPV basis(const int i) const
Return i-th basis vector.
Teuchos::RCP< std::vector< RCPV > > vecs_
std::vector< RCPV > dual_vecs_
void set(const V &x)
Set where .
RCPV clone() const
Clone to make a new (uninitialized) vector.
std::vector< PV >::size_type size_type
virtual void set(const Vector &x)
Set where .
Definition: ROL_Vector.hpp:196
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 .