Reference documentation for deal.II version 9.5.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
vector_adaptor.h
Go to the documentation of this file.
1//-----------------------------------------------------------
2//
3// Copyright (C) 2017 - 2022 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.md at
12// the top level directory of deal.II.
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
23
24# include <deal.II/lac/vector.h>
25
26# include <ROL_Vector.hpp>
27
28# include <limits>
29# include <type_traits>
30
31
33
39namespace Rol
40{
110 template <typename VectorType>
111 class VectorAdaptor : public ROL::Vector<typename VectorType::value_type>
112 {
116 using size_type = typename VectorType::size_type;
117
121 using value_type = typename VectorType::value_type;
122
126 using real_type = typename VectorType::real_type;
127
128 static_assert(std::is_convertible<real_type, value_type>::value,
129 "The real_type of the current VectorType is not "
130 "convertible to the value_type.");
131
132 private:
137 Teuchos::RCP<VectorType> vector_ptr;
138
139 public:
143 VectorAdaptor(const Teuchos::RCP<VectorType> &vector_ptr);
144
149 Teuchos::RCP<VectorType>
151
155 Teuchos::RCP<const VectorType>
156 getVector() const;
157
161 int
162 dimension() const;
163
172 void
173 set(const ROL::Vector<value_type> &rol_vector);
174
178 void
179 plus(const ROL::Vector<value_type> &rol_vector);
180
185 void
186 axpy(const value_type alpha, const ROL::Vector<value_type> &rol_vector);
187
191 void
192 scale(const value_type alpha);
193
198 dot(const ROL::Vector<value_type> &rol_vector) const;
199
210 norm() const;
211
215 Teuchos::RCP<ROL::Vector<value_type>>
216 clone() const;
217
223 Teuchos::RCP<ROL::Vector<value_type>>
224 basis(const int i) const;
225
229 void
230 applyUnary(const ROL::Elementwise::UnaryFunction<value_type> &f);
231
236 void
237 applyBinary(const ROL::Elementwise::BinaryFunction<value_type> &f,
238 const ROL::Vector<value_type> & rol_vector);
239
245 reduce(const ROL::Elementwise::ReductionOp<value_type> &r) const;
246
250 void
251 print(std::ostream &outStream) const;
252 };
253
254
255 /*------------------------------member definitions--------------------------*/
256# ifndef DOXYGEN
257
258
259 template <typename VectorType>
261 const Teuchos::RCP<VectorType> &vector_ptr)
262 : vector_ptr(vector_ptr)
263 {}
264
265
266
267 template <typename VectorType>
268 Teuchos::RCP<VectorType>
270 {
271 return vector_ptr;
272 }
273
274
275
276 template <typename VectorType>
277 Teuchos::RCP<const VectorType>
279 {
280 return vector_ptr;
281 }
282
283
284
285 template <typename VectorType>
286 void
287 VectorAdaptor<VectorType>::set(const ROL::Vector<value_type> &rol_vector)
288 {
289 const VectorAdaptor &vector_adaptor =
290 Teuchos::dyn_cast<const VectorAdaptor>(rol_vector);
291
292 (*vector_ptr) = *(vector_adaptor.getVector());
293 }
294
295
296
297 template <typename VectorType>
298 void
299 VectorAdaptor<VectorType>::plus(const ROL::Vector<value_type> &rol_vector)
300 {
301 Assert(this->dimension() == rol_vector.dimension(),
302 ExcDimensionMismatch(this->dimension(), rol_vector.dimension()));
303
304 const VectorAdaptor &vector_adaptor =
305 Teuchos::dyn_cast<const VectorAdaptor>(rol_vector);
306
307 *vector_ptr += *(vector_adaptor.getVector());
308 }
309
310
311
312 template <typename VectorType>
313 void
314 VectorAdaptor<VectorType>::axpy(const value_type alpha,
315 const ROL::Vector<value_type> &rol_vector)
316 {
317 Assert(this->dimension() == rol_vector.dimension(),
318 ExcDimensionMismatch(this->dimension(), rol_vector.dimension()));
319
320 const VectorAdaptor &vector_adaptor =
321 Teuchos::dyn_cast<const VectorAdaptor>(rol_vector);
322
323 vector_ptr->add(alpha, *(vector_adaptor.getVector()));
324 }
325
326
327
328 template <typename VectorType>
329 int
331 {
332 Assert(vector_ptr->size() < std::numeric_limits<int>::max(),
333 ExcMessage("The size of the vector being used is greater than "
334 "largest value of type int."));
335 return static_cast<int>(vector_ptr->size());
336 }
337
338
339
340 template <typename VectorType>
341 void
342 VectorAdaptor<VectorType>::scale(const value_type alpha)
343 {
344 (*vector_ptr) *= alpha;
345 }
346
347
348
349 template <typename VectorType>
350 typename VectorType::value_type
352 const ROL::Vector<value_type> &rol_vector) const
353 {
354 Assert(this->dimension() == rol_vector.dimension(),
355 ExcDimensionMismatch(this->dimension(), rol_vector.dimension()));
356
357 const VectorAdaptor &vector_adaptor =
358 Teuchos::dyn_cast<const VectorAdaptor>(rol_vector);
359
360 return (*vector_ptr) * (*vector_adaptor.getVector());
361 }
362
363
364
365 template <typename VectorType>
366 typename VectorType::value_type
368 {
369 return vector_ptr->l2_norm();
370 }
371
372
373
374 template <typename VectorType>
375 Teuchos::RCP<ROL::Vector<typename VectorType::value_type>>
377 {
378 Teuchos::RCP<VectorType> vec_ptr = Teuchos::rcp(new VectorType);
379 (*vec_ptr) = (*vector_ptr);
380
381 return Teuchos::rcp(new VectorAdaptor(vec_ptr));
382 }
383
384
385
386 template <typename VectorType>
387 Teuchos::RCP<ROL::Vector<typename VectorType::value_type>>
388 VectorAdaptor<VectorType>::basis(const int i) const
389 {
390 Teuchos::RCP<VectorType> vec_ptr = Teuchos::rcp(new VectorType);
391
392 // Zero all the entries in dealii vector.
393 vec_ptr->reinit(*vector_ptr, false);
394
395 if (vector_ptr->locally_owned_elements().is_element(i))
396 vec_ptr->operator[](i) = 1.;
397
398 if (vec_ptr->has_ghost_elements())
399 {
400 vec_ptr->update_ghost_values();
401 }
402 else
403 {
404 vec_ptr->compress(VectorOperation::insert);
405 }
406
407 Teuchos::RCP<VectorAdaptor> e = Teuchos::rcp(new VectorAdaptor(vec_ptr));
408
409 return e;
410 }
411
412
413
414 template <typename VectorType>
415 void
417 const ROL::Elementwise::UnaryFunction<value_type> &f)
418 {
419 const typename VectorType::iterator vend = vector_ptr->end();
420
421 for (typename VectorType::iterator iterator = vector_ptr->begin();
422 iterator != vend;
423 iterator++)
424 *iterator = f.apply(*iterator);
425
426 if (vector_ptr->has_ghost_elements())
427 {
428 vector_ptr->update_ghost_values();
429 }
430 else
431 {
432 vector_ptr->compress(VectorOperation::insert);
433 }
434 }
435
436
437
438 template <typename VectorType>
439 void
441 const ROL::Elementwise::BinaryFunction<value_type> &f,
442 const ROL::Vector<value_type> & rol_vector)
443 {
444 Assert(this->dimension() == rol_vector.dimension(),
445 ExcDimensionMismatch(this->dimension(), rol_vector.dimension()));
446
447 const VectorAdaptor &vector_adaptor =
448 Teuchos::dyn_cast<const VectorAdaptor>(rol_vector);
449
450 const VectorType &given_rol_vector = *(vector_adaptor.getVector());
451
452 const typename VectorType::iterator vend = vector_ptr->end();
453 const typename VectorType::const_iterator rolend = given_rol_vector.end();
454
455 typename VectorType::const_iterator r_iterator = given_rol_vector.begin();
456 for (typename VectorType::iterator l_iterator = vector_ptr->begin();
457 l_iterator != vend && r_iterator != rolend;
458 l_iterator++, r_iterator++)
459 *l_iterator = f.apply(*l_iterator, *r_iterator);
460
461 if (vector_ptr->has_ghost_elements())
462 {
463 vector_ptr->update_ghost_values();
464 }
465 else
466 {
467 vector_ptr->compress(VectorOperation::insert);
468 }
469 }
470
471
472
473 template <typename VectorType>
474 typename VectorType::value_type
476 const ROL::Elementwise::ReductionOp<value_type> &r) const
477 {
478 typename VectorType::value_type result = r.initialValue();
479
480 const typename VectorType::iterator vend = vector_ptr->end();
481
482 for (typename VectorType::iterator iterator = vector_ptr->begin();
483 iterator != vend;
484 iterator++)
485 r.reduce(*iterator, result);
486 // Parallel reduction among processes is handled internally.
487
488 return result;
489 }
490
491
492
493 template <typename VectorType>
494 void
495 VectorAdaptor<VectorType>::print(std::ostream &outStream) const
496 {
497 vector_ptr->print(outStream);
498 }
499
500
501# endif // DOXYGEN
502
503
504} // namespace Rol
505
506
508
509
510#endif // DEAL_II_TRILINOS_WITH_ROL
511
512#endif // dealii_optimization_rol_vector_adaptor_h
void applyBinary(const ROL::Elementwise::BinaryFunction< value_type > &f, const ROL::Vector< value_type > &rol_vector)
Teuchos::RCP< ROL::Vector< value_type > > clone() const
void scale(const value_type alpha)
value_type dot(const ROL::Vector< value_type > &rol_vector) const
void applyUnary(const ROL::Elementwise::UnaryFunction< value_type > &f)
Teuchos::RCP< VectorType > vector_ptr
Teuchos::RCP< VectorType > getVector()
int dimension() const
typename VectorType::value_type value_type
VectorAdaptor(const Teuchos::RCP< VectorType > &vector_ptr)
typename VectorType::real_type real_type
void print(std::ostream &outStream) const
typename VectorType::size_type size_type
void set(const ROL::Vector< value_type > &rol_vector)
Teuchos::RCP< const VectorType > getVector() const
void plus(const ROL::Vector< value_type > &rol_vector)
Teuchos::RCP< ROL::Vector< value_type > > basis(const int i) const
void axpy(const value_type alpha, const ROL::Vector< value_type > &rol_vector)
value_type reduce(const ROL::Elementwise::ReductionOp< value_type > &r) const
value_type norm() const
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
#define Assert(cond, exc)
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
static ::ExceptionBase & ExcMessage(std::string arg1)
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)