dune-istl  2.7.0
multitypeblockvector.hh
Go to the documentation of this file.
1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3 #ifndef DUNE_ISTL_MULTITYPEBLOCKVECTOR_HH
4 #define DUNE_ISTL_MULTITYPEBLOCKVECTOR_HH
5 
6 #include <cmath>
7 #include <iostream>
8 #include <tuple>
9 
10 #include <dune/common/dotproduct.hh>
11 #include <dune/common/ftraits.hh>
12 #include <dune/common/hybridutilities.hh>
13 #include <dune/common/typetraits.hh>
14 
15 #include "istlexception.hh"
16 
17 // forward declaration
18 namespace Dune {
19  template < typename... Args >
21 }
22 
23 #include "gsetc.hh"
24 
25 namespace Dune {
26 
38  template <typename Arg0, typename... Args>
39  struct FieldTraits< MultiTypeBlockVector<Arg0, Args...> >
40  {
41  typedef typename FieldTraits<Arg0>::field_type field_type;
42  typedef typename FieldTraits<Arg0>::real_type real_type;
43  };
53  template < typename... Args >
55  : public std::tuple<Args...>
56  {
58  typedef std::tuple<Args...> TupleType;
59  public:
60 
62  using size_type = std::size_t;
63 
67  using TupleType::TupleType;
68 
72  typedef MultiTypeBlockVector<Args...> type;
73 
81  typedef double field_type;
82 
88  static constexpr size_type size()
89  {
90  return sizeof...(Args);
91  }
92 
95  static constexpr size_type N()
96  {
97  return sizeof...(Args);
98  }
99 
103  int count() const DUNE_DEPRECATED_MSG("Use method 'N' instead")
104  {
105  return sizeof...(Args);
106  }
107 
109  size_type dim() const
110  {
111  size_type result = 0;
112  Hybrid::forEach(std::make_index_sequence<N()>{},
113  [&](auto i){result += std::get<i>(*this).dim();});
114 
115  return result;
116  }
117 
136  template< size_type index >
137  typename std::tuple_element<index,TupleType>::type&
138  operator[] ( const std::integral_constant< size_type, index > indexVariable )
139  {
140  DUNE_UNUSED_PARAMETER(indexVariable);
141  return std::get<index>(*this);
142  }
143 
149  template< size_type index >
150  const typename std::tuple_element<index,TupleType>::type&
151  operator[] ( const std::integral_constant< size_type, index > indexVariable ) const
152  {
153  DUNE_UNUSED_PARAMETER(indexVariable);
154  return std::get<index>(*this);
155  }
156 
159  template<typename T>
160  void operator= (const T& newval) {
161  Dune::Hybrid::forEach(*this, [&](auto&& entry) {
162  entry = newval;
163  });
164  }
165 
169  void operator+= (const type& newv) {
170  using namespace Dune::Hybrid;
171  forEach(integralRange(Hybrid::size(*this)), [&](auto&& i) {
172  (*this)[i] += newv[i];
173  });
174  }
175 
179  void operator-= (const type& newv) {
180  using namespace Dune::Hybrid;
181  forEach(integralRange(Hybrid::size(*this)), [&](auto&& i) {
182  (*this)[i] -= newv[i];
183  });
184  }
185 
187  template<class T,
188  std::enable_if_t< IsNumber<T>::value, int> = 0>
189  void operator*= (const T& w) {
190  Hybrid::forEach(*this, [&](auto&& entry) {
191  entry *= w;
192  });
193  }
194 
196  template<class T,
197  std::enable_if_t< IsNumber<T>::value, int> = 0>
198  void operator/= (const T& w) {
199  Hybrid::forEach(*this, [&](auto&& entry) {
200  entry /= w;
201  });
202  }
203 
204  field_type operator* (const type& newv) const {
205  using namespace Dune::Hybrid;
206  return accumulate(integralRange(Hybrid::size(*this)), field_type(0), [&](auto&& a, auto&& i) {
207  return a + (*this)[i]*newv[i];
208  });
209  }
210 
211  field_type dot (const type& newv) const {
212  using namespace Dune::Hybrid;
213  return accumulate(integralRange(Hybrid::size(*this)), field_type(0), [&](auto&& a, auto&& i) {
214  return a + (*this)[i].dot(newv[i]);
215  });
216  }
217 
220  auto one_norm() const {
221  using namespace Dune::Hybrid;
222  return accumulate(*this, typename FieldTraits<field_type>::real_type(0), [&](auto&& a, auto&& entry) {
223  return a + entry.one_norm();
224  });
225  }
226 
229  auto one_norm_real() const {
230  using namespace Dune::Hybrid;
231  return accumulate(*this, typename FieldTraits<field_type>::real_type(0), [&](auto&& a, auto&& entry) {
232  return a + entry.one_norm_real();
233  });
234  }
235 
238  typename FieldTraits<field_type>::real_type two_norm2() const {
239  using namespace Dune::Hybrid;
240  return accumulate(*this, typename FieldTraits<field_type>::real_type(0), [&](auto&& a, auto&& entry) {
241  return a + entry.two_norm2();
242  });
243  }
244 
247  typename FieldTraits<field_type>::real_type two_norm() const {return sqrt(this->two_norm2());}
248 
251  typename FieldTraits<field_type>::real_type infinity_norm() const
252  {
253  using namespace Dune::Hybrid;
254  using std::max;
255  using real_type = typename FieldTraits<field_type>::real_type;
256 
257  real_type result = 0.0;
258  // Compute max norm tracking appearing nan values
259  // if the field type supports nan.
260  ifElse(HasNaN<field_type>(), [&](auto&& id) {
261  // This variable will preserve any nan value
262  real_type nanTracker = 1.0;
263  using namespace Dune::Hybrid; // needed for icc, see issue #31
264  forEach(*this, [&](auto&& entry) {
265  real_type entryNorm = entry.infinity_norm();
266  result = max(entryNorm, result);
267  nanTracker += entryNorm;
268  });
269  // Incorporate possible nan value into result
270  result *= (nanTracker / nanTracker);
271  }, [&](auto&& id) {
272  using namespace Dune::Hybrid; // needed for icc, see issue #31
273  forEach(*this, [&](auto&& entry) {
274  result = max(entry.infinity_norm(), result);
275  });
276  });
277  return result;
278  }
279 
282  typename FieldTraits<field_type>::real_type infinity_norm_real() const
283  {
284  using namespace Dune::Hybrid;
285  using std::max;
286  using real_type = typename FieldTraits<field_type>::real_type;
287 
288  real_type result = 0.0;
289  // Compute max norm tracking appearing nan values
290  // if the field type supports nan.
291  ifElse(HasNaN<field_type>(), [&](auto&& id) {
292  // This variable will preserve any nan value
293  real_type nanTracker = 1.0;
294  using namespace Dune::Hybrid; // needed for icc, see issue #31
295  forEach(*this, [&](auto&& entry) {
296  real_type entryNorm = entry.infinity_norm_real();
297  result = max(entryNorm, result);
298  nanTracker += entryNorm;
299  });
300  // Incorporate possible nan value into result
301  result *= (nanTracker / nanTracker);
302  }, [&](auto&& id) {
303  using namespace Dune::Hybrid; // needed for icc, see issue #31
304  forEach(*this, [&](auto&& entry) {
305  result = max(entry.infinity_norm_real(), result);
306  });
307  });
308  return result;
309  }
310 
315  template<typename Ta>
316  void axpy (const Ta& a, const type& y) {
317  using namespace Dune::Hybrid;
318  forEach(integralRange(Hybrid::size(*this)), [&](auto&& i) {
319  (*this)[i].axpy(a, y[i]);
320  });
321  }
322 
323  };
324 
325 
326 
329  template <typename... Args>
330  std::ostream& operator<< (std::ostream& s, const MultiTypeBlockVector<Args...>& v) {
331  using namespace Dune::Hybrid;
332  forEach(integralRange(Dune::Hybrid::size(v)), [&](auto&& i) {
333  s << "\t(" << i << "):\n" << v[i] << "\n";
334  });
335  return s;
336  }
337 
338 
339 
340 } // end namespace
341 
342 #endif
gsetc.hh
Simple iterative methods like Jacobi, Gauss-Seidel, SOR, SSOR, etc. in a generic way.
Dune::MultiTypeBlockVector::operator/=
void operator/=(const T &w)
Division by a scalar.
Definition: multitypeblockvector.hh:198
Dune::MultiTypeBlockVector::two_norm
FieldTraits< field_type >::real_type two_norm() const
Compute the Euclidean norm.
Definition: multitypeblockvector.hh:247
Dune::MultiTypeBlockVector::operator+=
void operator+=(const type &newv)
Definition: multitypeblockvector.hh:169
Dune::MultiTypeBlockVector::type
MultiTypeBlockVector< Args... > type
Definition: multitypeblockvector.hh:72
Dune::MultiTypeBlockVector::dim
size_type dim() const
Number of scalar elements.
Definition: multitypeblockvector.hh:109
Dune::MultiTypeBlockVector::operator*
field_type operator*(const type &newv) const
Definition: multitypeblockvector.hh:204
istlexception.hh
Dune::MultiTypeBlockVector::operator-=
void operator-=(const type &newv)
Definition: multitypeblockvector.hh:179
Dune::MultiTypeBlockVector::operator*=
void operator*=(const T &w)
Multiplication with a scalar.
Definition: multitypeblockvector.hh:189
Dune::MultiTypeBlockVector::two_norm2
FieldTraits< field_type >::real_type two_norm2() const
Compute the squared Euclidean norm.
Definition: multitypeblockvector.hh:238
Dune::MultiTypeBlockVector::count
int count() const DUNE_DEPRECATED_MSG("Use method 'N' instead")
Definition: multitypeblockvector.hh:103
Dune::MultiTypeBlockVector::field_type
double field_type
The type used for scalars.
Definition: multitypeblockvector.hh:81
Dune::MultiTypeBlockVector::one_norm
auto one_norm() const
Compute the 1-norm.
Definition: multitypeblockvector.hh:220
Dune::MultiTypeBlockVector::axpy
void axpy(const Ta &a, const type &y)
Axpy operation on this vector (*this += a * y)
Definition: multitypeblockvector.hh:316
Dune::MultiTypeBlockVector
A Vector class to support different block types.
Definition: multitypeblockvector.hh:20
Dune
Definition: allocator.hh:7
Dune::FieldTraits< MultiTypeBlockVector< Arg0, Args... > >::real_type
FieldTraits< Arg0 >::real_type real_type
Definition: multitypeblockvector.hh:42
Dune::MultiTypeBlockVector::operator=
void operator=(const T &newval)
Assignment operator.
Definition: multitypeblockvector.hh:160
Dune::MultiTypeBlockVector::infinity_norm
FieldTraits< field_type >::real_type infinity_norm() const
Compute the maximum norm.
Definition: multitypeblockvector.hh:251
Dune::MultiTypeBlockVector::infinity_norm_real
FieldTraits< field_type >::real_type infinity_norm_real() const
Compute the simplified maximum norm (uses 1-norm for complex values)
Definition: multitypeblockvector.hh:282
Dune::MultiTypeBlockVector::N
static constexpr size_type N()
Number of elements.
Definition: multitypeblockvector.hh:95
Dune::operator<<
std::ostream & operator<<(std::ostream &s, const BlockVector< K, A > &v)
Send BlockVector to an output stream.
Definition: bvector.hh:650
Dune::MultiTypeBlockVector::one_norm_real
auto one_norm_real() const
Compute the simplified 1-norm (uses 1-norm also for complex values)
Definition: multitypeblockvector.hh:229
Dune::MultiTypeBlockVector::operator[]
std::tuple_element< index, TupleType >::type & operator[](const std::integral_constant< size_type, index > indexVariable)
Random-access operator.
Definition: multitypeblockvector.hh:138
Dune::MultiTypeBlockVector::size_type
std::size_t size_type
Type used for vector sizes.
Definition: multitypeblockvector.hh:62
Dune::MultiTypeBlockVector::dot
field_type dot(const type &newv) const
Definition: multitypeblockvector.hh:211
Dune::FieldTraits< MultiTypeBlockVector< Arg0, Args... > >::field_type
FieldTraits< Arg0 >::field_type field_type
Definition: multitypeblockvector.hh:41
Dune::MultiTypeBlockVector::size
static constexpr size_type size()
Return the number of non-zero vector entries.
Definition: multitypeblockvector.hh:88