Reference documentation for deal.II version 9.0.0
vector_adaptor.h
1 //-----------------------------------------------------------
2 //
3 // Copyright (C) 2017 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE at
12 // the top level of the deal.II distribution.
13 //
14 //---------------------------------------------------------------
15 
16 #ifndef dealii_optimization_rol_vector_adaptor_h
17 #define dealii_optimization_rol_vector_adaptor_h
18 
19 #include <deal.II/base/config.h>
20 
21 #ifdef DEAL_II_TRILINOS_WITH_ROL
22 #include <ROL_Vector.hpp>
23 
24 #include <deal.II/base/exceptions.h>
25 #include <deal.II/lac/vector.h>
26 
27 #include <type_traits>
28 
29 
30 DEAL_II_NAMESPACE_OPEN
31 
32 
33 using namespace dealii;
34 
35 
41 namespace Rol
42 {
43 
114  template<typename VectorType>
115  class VectorAdaptor : public ROL::Vector<typename VectorType::value_type>
116  {
117 
121  using size_type = typename VectorType::size_type;
122 
126  using value_type = typename VectorType::value_type;
127 
131  using real_type = typename VectorType::real_type;
132 
133  static_assert (std::is_convertible<real_type, value_type>::value,
134  "The real_type of the current VectorType is not "
135  "convertible to the value_type.");
136 
137  private:
138 
143  Teuchos::RCP<VectorType> vector_ptr;
144 
145  public:
146 
150  VectorAdaptor (const Teuchos::RCP<VectorType> &vector_ptr);
151 
156  Teuchos::RCP<VectorType> getVector ();
157 
161  Teuchos::RCP<const VectorType> getVector () const;
162 
166  int dimension () const;
167 
176  void set (const ROL::Vector<value_type> &rol_vector);
177 
181  void plus (const ROL::Vector<value_type> &rol_vector);
182 
187  void axpy (const value_type alpha,
188  const ROL::Vector<value_type> &rol_vector);
189 
193  void scale (const value_type alpha);
194 
198  value_type dot (const ROL::Vector<value_type> &rol_vector) const;
199 
209  value_type norm() const;
210 
214  Teuchos::RCP<ROL::Vector<value_type>> clone() const;
215 
221  Teuchos::RCP<ROL::Vector<value_type>> basis (const int i) const;
222 
226  void applyUnary (const ROL::Elementwise::UnaryFunction<value_type> &f);
227 
232  void applyBinary (const ROL::Elementwise::BinaryFunction<value_type> &f,
233  const ROL::Vector<value_type> &rol_vector);
234 
239  value_type reduce (const ROL::Elementwise::ReductionOp<value_type> &r) const;
240 
244  void print (std::ostream &outStream) const;
245 
246  };
247 
248 
249  /*------------------------------member definitions--------------------------*/
250 #ifndef DOXYGEN
251 
252 
253  template<typename VectorType>
255  VectorAdaptor (const Teuchos::RCP<VectorType> &vector_ptr)
256  :
257  vector_ptr (vector_ptr)
258  {}
259 
260 
261 
262  template<typename VectorType>
263  Teuchos::RCP<VectorType>
265  {
266  return vector_ptr;
267  }
268 
269 
270 
271  template<typename VectorType>
272  Teuchos::RCP<const VectorType>
274  {
275  return vector_ptr;
276  }
277 
278 
279 
280  template<typename VectorType>
281  void
282  VectorAdaptor<VectorType>::set (const ROL::Vector<value_type> &rol_vector)
283  {
284  const VectorAdaptor &vector_adaptor =
285  Teuchos::dyn_cast<const VectorAdaptor>(rol_vector);
286 
287  (*vector_ptr) = *(vector_adaptor.getVector());
288  }
289 
290 
291 
292  template<typename VectorType>
293  void
294  VectorAdaptor<VectorType>::plus (const ROL::Vector<value_type> &rol_vector)
295  {
296  Assert (this->dimension() == rol_vector.dimension(),
297  ExcDimensionMismatch(this->dimension(), rol_vector.dimension()));
298 
299  const VectorAdaptor &vector_adaptor =
300  Teuchos::dyn_cast<const VectorAdaptor>(rol_vector);
301 
302  *vector_ptr += *(vector_adaptor.getVector());
303  }
304 
305 
306 
307  template<typename VectorType>
308  void
309  VectorAdaptor<VectorType>::axpy (const value_type alpha,
310  const ROL::Vector<value_type> &rol_vector)
311  {
312  Assert (this->dimension() == rol_vector.dimension(),
313  ExcDimensionMismatch(this->dimension(), rol_vector.dimension()));
314 
315  const VectorAdaptor &vector_adaptor =
316  Teuchos::dyn_cast<const VectorAdaptor>(rol_vector);
317 
318  vector_ptr->add (alpha, *(vector_adaptor.getVector()));
319  }
320 
321 
322 
323  template<typename VectorType>
324  int
326  {
327  Assert (vector_ptr->size() < std::numeric_limits<int>::max(),
328  ExcMessage("The size of the vector being used is greater than "
329  "largest value of type int."));
330  return static_cast<int>(vector_ptr->size());
331  }
332 
333 
334 
335  template<typename VectorType>
336  void
337  VectorAdaptor<VectorType>::scale (const value_type alpha)
338  {
339  (*vector_ptr) *= alpha;
340  }
341 
342 
343 
344  template<typename VectorType>
345  typename VectorType::value_type
347  dot (const ROL::Vector<value_type> &rol_vector) const
348  {
349  Assert (this->dimension() == rol_vector.dimension(),
350  ExcDimensionMismatch(this->dimension(), rol_vector.dimension()));
351 
352  const VectorAdaptor &vector_adaptor =
353  Teuchos::dyn_cast< const VectorAdaptor>(rol_vector);
354 
355  return (*vector_ptr) * (*vector_adaptor.getVector());
356  }
357 
358 
359 
360  template<typename VectorType>
361  typename VectorType::value_type
363  {
364  return vector_ptr->l2_norm();
365  }
366 
367 
368 
369  template<typename VectorType>
370  Teuchos::RCP<ROL::Vector<typename VectorType::value_type> >
372  {
373  Teuchos::RCP< VectorType> vec_ptr = Teuchos::rcp (new VectorType);
374  (*vec_ptr) = (*vector_ptr);
375 
376  return Teuchos::rcp (new VectorAdaptor (vec_ptr));
377  }
378 
379 
380 
381  template<typename VectorType>
382  Teuchos::RCP<ROL::Vector<typename VectorType::value_type> >
383  VectorAdaptor<VectorType>::basis (const int i) const
384  {
385  Teuchos::RCP<VectorType> vec_ptr = Teuchos::rcp (new VectorType);
386 
387  // Zero all the entries in dealii vector.
388  vec_ptr->reinit(*vector_ptr, false);
389 
390  if (vector_ptr->locally_owned_elements().is_element(i))
391  vec_ptr->operator[](i) = 1.;
392 
393  vec_ptr->compress(VectorOperation::insert);
394 
395  Teuchos::RCP<VectorAdaptor> e = Teuchos::rcp (new VectorAdaptor(vec_ptr));
396 
397  return e;
398  }
399 
400 
401 
402  template<typename VectorType>
403  void
405  applyUnary (const ROL::Elementwise::UnaryFunction<value_type> &f)
406  {
407  const typename VectorType::iterator vend = vector_ptr->end();
408 
409  for (typename VectorType::iterator
410  iterator = vector_ptr->begin();
411  iterator != vend;
412  iterator++)
413  *iterator = f.apply(*iterator);
414 
415  vector_ptr->compress (VectorOperation::insert);
416  }
417 
418 
419 
420  template<typename VectorType>
421  void
423  applyBinary (const ROL::Elementwise::BinaryFunction<value_type> &f,
424  const ROL::Vector<value_type> &rol_vector)
425  {
426  Assert (this->dimension() == rol_vector.dimension(),
427  ExcDimensionMismatch(this->dimension(), rol_vector.dimension()));
428 
429  const VectorAdaptor &vector_adaptor =
430  Teuchos::dyn_cast<const VectorAdaptor>(rol_vector);
431 
432  const VectorType &given_rol_vector = *(vector_adaptor.getVector());
433 
434  const typename VectorType::iterator vend = vector_ptr->end();
435  const typename VectorType::const_iterator rolend = given_rol_vector.end();
436 
437  typename VectorType::const_iterator r_iterator = given_rol_vector.begin();
438  for (typename VectorType::iterator l_iterator = vector_ptr->begin();
439  l_iterator != vend && r_iterator != rolend;
440  l_iterator++, r_iterator++)
441  *l_iterator = f.apply(*l_iterator, *r_iterator);
442 
443  vector_ptr->compress (VectorOperation::insert);
444  }
445 
446 
447 
448  template<typename VectorType>
449  typename VectorType::value_type
451  reduce (const ROL::Elementwise::ReductionOp<value_type> &r) const
452  {
453  typename VectorType::value_type result = r.initialValue();
454 
455  const typename VectorType::iterator vend = vector_ptr->end();
456 
457  for (typename VectorType::iterator
458  iterator = vector_ptr->begin();
459  iterator != vend;
460  iterator++)
461  r.reduce(*iterator, result);
462  // Parallel reduction among processes is handled internally.
463 
464  return result;
465  }
466 
467 
468 
469  template<typename VectorType>
470  void
471  VectorAdaptor<VectorType>::print (std::ostream &outStream) const
472  {
473  vector_ptr->print(outStream);
474  }
475 
476 
477 #endif // DOXYGEN
478 
479 
480 } // namespace Rol
481 
482 
483 DEAL_II_NAMESPACE_CLOSE
484 
485 
486 #endif // DEAL_II_TRILINOS_WITH_ROL
487 
488 #endif // dealii_optimization_rol_vector_adaptor_h
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
void applyBinary(const ROL::Elementwise::BinaryFunction< value_type > &f, const ROL::Vector< value_type > &rol_vector)
Teuchos::RCP< ROL::Vector< value_type > > clone() const
Teuchos::RCP< VectorType > getVector()
void applyUnary(const ROL::Elementwise::UnaryFunction< value_type > &f)
value_type norm() const
Teuchos::RCP< VectorType > vector_ptr
typename VectorType::value_type value_type
static ::ExceptionBase & ExcMessage(std::string arg1)
void plus(const ROL::Vector< value_type > &rol_vector)
#define Assert(cond, exc)
Definition: exceptions.h:1142
value_type dot(const ROL::Vector< value_type > &rol_vector) const
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
void print(std::ostream &outStream) const
int dimension() const
value_type reduce(const ROL::Elementwise::ReductionOp< value_type > &r) const
Teuchos::RCP< ROL::Vector< value_type > > basis(const int i) const
typename VectorType::real_type real_type
typename VectorType::size_type size_type
void axpy(const value_type alpha, const ROL::Vector< value_type > &rol_vector)
void set(const ROL::Vector< value_type > &rol_vector)
void scale(const value_type alpha)