Reference documentation for deal.II version 9.2.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\}}\)
packaged_operation.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2014 - 2020 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_packaged_operation_h
17 #define dealii_packaged_operation_h
18 
19 #include <deal.II/base/config.h>
20 
22 
24 
25 #include <functional>
26 
28 
29 // Forward declarations:
30 #ifndef DOXYGEN
31 template <typename Number>
32 class Vector;
33 template <typename Range, typename Domain, typename Payload>
34 class LinearOperator;
35 template <typename Range = Vector<double>>
36 class PackagedOperation;
37 #endif
38 
39 
105 template <typename Range>
107 {
108 public:
115  {
116  apply = [](Range &) {
117  Assert(false,
118  ExcMessage(
119  "Uninitialized PackagedOperation<Range>::apply called"));
120  };
121 
122  apply_add = [](Range &) {
123  Assert(false,
124  ExcMessage(
125  "Uninitialized PackagedOperation<Range>::apply_add called"));
126  };
127 
128  reinit_vector = [](Range &, bool) {
129  Assert(false,
130  ExcMessage("Uninitialized PackagedOperation<Range>::reinit_vector "
131  "method called"));
132  };
133  }
134 
138  PackagedOperation(const PackagedOperation<Range> &) = default;
139 
149  PackagedOperation(const Range &u)
150  {
151  *this = u;
152  }
153 
158  operator=(const PackagedOperation<Range> &) = default;
159 
170  operator=(const Range &u)
171  {
172  apply = [&u](Range &v) { v = u; };
173 
174  apply_add = [&u](Range &v) { v += u; };
175 
176  reinit_vector = [&u](Range &v, bool omit_zeroing_entries) {
177  v.reinit(u, omit_zeroing_entries);
178  };
179 
180  return *this;
181  }
182 
189  operator Range() const
190  {
191  Range result_vector;
192 
193  reinit_vector(result_vector, /*bool omit_zeroing_entries=*/true);
194  apply(result_vector);
195 
196  return result_vector;
197  }
198 
203 
209  {
210  *this = *this + second_comp;
211  return *this;
212  }
213 
220  {
221  *this = *this - second_comp;
222  return *this;
223  }
224 
230  operator+=(const Range &offset)
231  {
232  *this = *this + PackagedOperation<Range>(offset);
233  return *this;
234  }
235 
241  operator-=(const Range &offset)
242  {
243  *this = *this - PackagedOperation<Range>(offset);
244  return *this;
245  }
246 
251  operator*=(typename Range::value_type number)
252  {
253  *this = *this * number;
254  return *this;
255  }
257 
262  std::function<void(Range &v)> apply;
263 
268  std::function<void(Range &v)> apply_add;
269 
277  std::function<void(Range &v, bool omit_zeroing_entries)> reinit_vector;
278 };
279 
280 
285 
294 template <typename Range>
297  const PackagedOperation<Range> &second_comp)
298 {
299  PackagedOperation<Range> return_comp;
300 
301  return_comp.reinit_vector = first_comp.reinit_vector;
302 
303  // ensure to have valid PackagedOperation objects by catching first_comp and
304  // second_comp by value
305 
306  return_comp.apply = [first_comp, second_comp](Range &v) {
307  first_comp.apply(v);
308  second_comp.apply_add(v);
309  };
310 
311  return_comp.apply_add = [first_comp, second_comp](Range &v) {
312  first_comp.apply_add(v);
313  second_comp.apply_add(v);
314  };
315 
316  return return_comp;
317 }
318 
327 template <typename Range>
330  const PackagedOperation<Range> &second_comp)
331 {
332  PackagedOperation<Range> return_comp;
333 
334  return_comp.reinit_vector = first_comp.reinit_vector;
335 
336  // ensure to have valid PackagedOperation objects by catching first_comp and
337  // second_comp by value
338 
339  return_comp.apply = [first_comp, second_comp](Range &v) {
340  second_comp.apply(v);
341  v *= -1.;
342  first_comp.apply_add(v);
343  };
344 
345  return_comp.apply_add = [first_comp, second_comp](Range &v) {
346  first_comp.apply_add(v);
347  v *= -1.;
348  second_comp.apply_add(v);
349  v *= -1.;
350  };
351 
352  return return_comp;
353 }
354 
363 template <typename Range>
365  typename Range::value_type number)
366 {
367  PackagedOperation<Range> return_comp;
368 
369  return_comp.reinit_vector = comp.reinit_vector;
370 
371  // the trivial case: number is zero
372  if (number == 0.)
373  {
374  return_comp.apply = [](Range &v) { v = 0.; };
375 
376  return_comp.apply_add = [](Range &) {};
377  }
378  else
379  {
380  return_comp.apply = [comp, number](Range &v) {
381  comp.apply(v);
382  v *= number;
383  };
384 
385  return_comp.apply_add = [comp, number](Range &v) {
386  v /= number;
387  comp.apply_add(v);
388  v *= number;
389  };
390  }
391 
392  return return_comp;
393 }
394 
403 template <typename Range>
404 PackagedOperation<Range> operator*(typename Range::value_type number,
405  const PackagedOperation<Range> &comp)
406 {
407  return comp * number;
408 }
409 
418 template <typename Range>
420 operator+(const PackagedOperation<Range> &comp, const Range &offset)
421 {
422  return comp + PackagedOperation<Range>(offset);
423 }
424 
433 template <typename Range>
435 operator+(const Range &offset, const PackagedOperation<Range> &comp)
436 {
437  return PackagedOperation<Range>(offset) + comp;
438 }
439 
448 template <typename Range>
450 operator-(const PackagedOperation<Range> &comp, const Range &offset)
451 {
452  return comp - PackagedOperation<Range>(offset);
453 }
454 
455 
465 template <typename Range>
467 operator-(const Range &offset, const PackagedOperation<Range> &comp)
468 {
469  return PackagedOperation<Range>(offset) - comp;
470 }
471 
473 
474 
479 
480 namespace internal
481 {
482  namespace PackagedOperationImplementation
483  {
484  // Poor man's trait class that determines whether type T is a vector:
485  // FIXME: Implement this as a proper type trait - similar to
486  // isBlockVector
487 
488  template <typename T>
490  {
491  template <typename C>
492  static std::false_type
493  test(...);
494 
495  template <typename C>
496  static std::true_type
497  test(decltype(&C::operator+=),
498  decltype(&C::operator-=),
499  decltype(&C::l2_norm));
500 
501  public:
502  // type is std::true_type if Matrix provides vmult_add and Tvmult_add,
503  // otherwise it is std::false_type
504 
505  using type = decltype(test<T>(nullptr, nullptr, nullptr));
506  }; // namespace
507  } // namespace PackagedOperationImplementation
508 } // namespace internal
509 
510 
525 template <typename Range,
526  typename = typename std::enable_if<
528  Range>::type::value>::type>
530 operator+(const Range &u, const Range &v)
531 {
532  PackagedOperation<Range> return_comp;
533 
534  // ensure to have valid PackagedOperation objects by catching op by value
535  // u is caught by reference
536 
537  return_comp.reinit_vector = [&u](Range &x, bool omit_zeroing_entries) {
538  x.reinit(u, omit_zeroing_entries);
539  };
540 
541  return_comp.apply = [&u, &v](Range &x) {
542  x = u;
543  x += v;
544  };
545 
546  return_comp.apply_add = [&u, &v](Range &x) {
547  x += u;
548  x += v;
549  };
550 
551  return return_comp;
552 }
553 
554 
570 template <typename Range,
571  typename = typename std::enable_if<
573  Range>::type::value>::type>
575 operator-(const Range &u, const Range &v)
576 {
577  PackagedOperation<Range> return_comp;
578 
579  // ensure to have valid PackagedOperation objects by catching op by value
580  // u is caught by reference
581 
582  return_comp.reinit_vector = [&u](Range &x, bool omit_zeroing_entries) {
583  x.reinit(u, omit_zeroing_entries);
584  };
585 
586  return_comp.apply = [&u, &v](Range &x) {
587  x = u;
588  x -= v;
589  };
590 
591  return_comp.apply_add = [&u, &v](Range &x) {
592  x += u;
593  x -= v;
594  };
595 
596  return return_comp;
597 }
598 
599 
614 template <typename Range,
615  typename = typename std::enable_if<
617  Range>::type::value>::type>
619  typename Range::value_type number)
620 {
621  return PackagedOperation<Range>(u) * number;
622 }
623 
624 
639 template <typename Range,
640  typename = typename std::enable_if<
642  Range>::type::value>::type>
643 PackagedOperation<Range> operator*(typename Range::value_type number,
644  const Range & u)
645 {
646  return number * PackagedOperation<Range>(u);
647 }
648 
649 
666 template <typename Range, typename Domain, typename Payload>
669 {
670  PackagedOperation<Range> return_comp;
671 
672  return_comp.reinit_vector = op.reinit_range_vector;
673 
674  // ensure to have valid PackagedOperation objects by catching op by value
675  // u is caught by reference
676 
677  return_comp.apply = [op, &u](Range &v) { op.vmult(v, u); };
678 
679  return_comp.apply_add = [op, &u](Range &v) { op.vmult_add(v, u); };
680 
681  return return_comp;
682 }
683 
684 
701 template <typename Range, typename Domain, typename Payload>
704 {
705  PackagedOperation<Range> return_comp;
706 
707  return_comp.reinit_vector = op.reinit_domain_vector;
708 
709  // ensure to have valid PackagedOperation objects by catching op by value
710  // u is caught by reference
711 
712  return_comp.apply = [op, &u](Domain &v) { op.Tvmult(v, u); };
713 
714  return_comp.apply_add = [op, &u](Domain &v) { op.Tvmult_add(v, u); };
715 
716  return return_comp;
717 }
718 
719 
728 template <typename Range, typename Domain, typename Payload>
731  const PackagedOperation<Domain> & comp)
732 {
733  PackagedOperation<Range> return_comp;
734 
735  return_comp.reinit_vector = op.reinit_range_vector;
736 
737  // ensure to have valid PackagedOperation objects by catching op by value
738  // u is caught by reference
739 
740  return_comp.apply = [op, comp](Domain &v) {
741  GrowingVectorMemory<Range> vector_memory;
742 
743  typename VectorMemory<Range>::Pointer i(vector_memory);
744  op.reinit_domain_vector(*i, /*bool omit_zeroing_entries =*/true);
745 
746  comp.apply(*i);
747  op.vmult(v, *i);
748  };
749 
750  return_comp.apply_add = [op, comp](Domain &v) {
751  GrowingVectorMemory<Range> vector_memory;
752 
753  typename VectorMemory<Range>::Pointer i(vector_memory);
754  op.reinit_range_vector(*i, /*bool omit_zeroing_entries =*/true);
755 
756  comp.apply(*i);
757  op.vmult_add(v, *i);
758  };
759 
760  return return_comp;
761 }
762 
763 
772 template <typename Range, typename Domain, typename Payload>
776 {
777  PackagedOperation<Range> return_comp;
778 
779  return_comp.reinit_vector = op.reinit_domain_vector;
780 
781  // ensure to have valid PackagedOperation objects by catching op by value
782  // u is caught by reference
783 
784  return_comp.apply = [op, comp](Domain &v) {
785  GrowingVectorMemory<Range> vector_memory;
786 
787  typename VectorMemory<Range>::Pointer i(vector_memory);
788  op.reinit_range_vector(*i, /*bool omit_zeroing_entries =*/true);
789 
790  comp.apply(*i);
791  op.Tvmult(v, *i);
792  };
793 
794  return_comp.apply_add = [op, comp](Domain &v) {
795  GrowingVectorMemory<Range> vector_memory;
796 
797  typename VectorMemory<Range>::Pointer i(vector_memory);
798  op.reinit_range_vector(*i, /*bool omit_zeroing_entries =*/true);
799 
800  comp.apply(*i);
801  op.Tvmult_add(v, *i);
802  };
803 
804  return return_comp;
805 }
806 
808 
810 
811 #endif
PackagedOperation< Range > & operator+=(const PackagedOperation< Range > &second_comp)
PackagedOperation(const Range &u)
std::function< void(Range &v)> apply
decltype(test< T >(nullptr, nullptr, nullptr)) type
std::function< void(Range &v, bool omit_zeroing_entries)> reinit_vector
std::function< void(Range &v, bool omit_zeroing_entries)> reinit_range_vector
std::vector< value_type > l2_norm(const typename ::Triangulation< dim, spacedim >::cell_iterator &parent, const value_type parent_value)
std::function< void(Domain &v, const Range &u)> Tvmult_add
std::function< void(Domain &v, const Range &u)> Tvmult
PackagedOperation< Range > & operator-=(const PackagedOperation< Range > &second_comp)
PackagedOperation< Range > & operator=(const PackagedOperation< Range > &)=default
static ::ExceptionBase & ExcMessage(std::string arg1)
PackagedOperation< Range > & operator=(const Range &u)
#define Assert(cond, exc)
Definition: exceptions.h:1419
std::function< void(Range &v, const Domain &u)> vmult
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:359
PackagedOperation< Range > operator-(const PackagedOperation< Range > &first_comp, const PackagedOperation< Range > &second_comp)
std::function< void(Domain &v, bool omit_zeroing_entries)> reinit_domain_vector
PackagedOperation< Range > & operator*=(typename Range::value_type number)
std::function< void(Range &v)> apply_add
std::function< void(Range &v, const Domain &u)> vmult_add
PackagedOperation< Range > & operator+=(const Range &offset)
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:358
PackagedOperation< Range > operator*(const PackagedOperation< Range > &comp, typename Range::value_type number)
static const bool value
PackagedOperation< Range > operator+(const PackagedOperation< Range > &first_comp, const PackagedOperation< Range > &second_comp)
PackagedOperation< Range > & operator-=(const Range &offset)