3#ifndef DUNE_ISTL_MULTITYPEBLOCKVECTOR_HH
4#define DUNE_ISTL_MULTITYPEBLOCKVECTOR_HH
10#include <dune/common/dotproduct.hh>
11#include <dune/common/ftraits.hh>
12#include <dune/common/hybridutilities.hh>
13#include <dune/common/typetraits.hh>
19 template <
typename... Args >
20 class MultiTypeBlockVector;
38 template <
typename Arg0,
typename... Args>
41 typedef typename FieldTraits<Arg0>::field_type
field_type;
42 typedef typename FieldTraits<Arg0>::real_type
real_type;
53 template <
typename... Args >
55 :
public std::tuple<Args...>
58 typedef std::tuple<Args...> TupleType;
67 using TupleType::TupleType;
90 return sizeof...(Args);
97 return sizeof...(Args);
106 [[deprecated(
"Use method 'N' instead")]]
109 return sizeof...(Args);
116 Hybrid::forEach(std::make_index_sequence<
N()>{},
117 [&](
auto i){result += std::get<i>(*this).dim();});
140 template<
size_type index >
141 typename std::tuple_element<index,TupleType>::type&
142 operator[] ([[maybe_unused]]
const std::integral_constant< size_type, index > indexVariable)
144 return std::get<index>(*
this);
152 template<
size_type index >
153 const typename std::tuple_element<index,TupleType>::type&
154 operator[] ([[maybe_unused]]
const std::integral_constant< size_type, index > indexVariable)
const
156 return std::get<index>(*
this);
163 Dune::Hybrid::forEach(*
this, [&](
auto&& entry) {
172 using namespace Dune::Hybrid;
173 forEach(integralRange(Hybrid::size(*
this)), [&](
auto&& i) {
174 (*this)[i] += newv[i];
182 using namespace Dune::Hybrid;
183 forEach(integralRange(Hybrid::size(*
this)), [&](
auto&& i) {
184 (*this)[i] -= newv[i];
190 std::enable_if_t< IsNumber<T>::value,
int> = 0>
192 Hybrid::forEach(*
this, [&](
auto&& entry) {
199 std::enable_if_t< IsNumber<T>::value,
int> = 0>
201 Hybrid::forEach(*
this, [&](
auto&& entry) {
207 using namespace Dune::Hybrid;
208 return accumulate(integralRange(Hybrid::size(*
this)),
field_type(0), [&](
auto&& a,
auto&& i) {
209 return a + (*this)[i]*newv[i];
214 using namespace Dune::Hybrid;
215 return accumulate(integralRange(Hybrid::size(*
this)),
field_type(0), [&](
auto&& a,
auto&& i) {
216 return a + (*this)[i].dot(newv[i]);
223 using namespace Dune::Hybrid;
224 return accumulate(*
this,
typename FieldTraits<field_type>::real_type(0), [&](
auto&& a,
auto&& entry) {
225 return a + entry.one_norm();
232 using namespace Dune::Hybrid;
233 return accumulate(*
this,
typename FieldTraits<field_type>::real_type(0), [&](
auto&& a,
auto&& entry) {
234 return a + entry.one_norm_real();
240 typename FieldTraits<field_type>::real_type
two_norm2()
const {
241 using namespace Dune::Hybrid;
242 return accumulate(*
this,
typename FieldTraits<field_type>::real_type(0), [&](
auto&& a,
auto&& entry) {
243 return a + entry.two_norm2();
255 using namespace Dune::Hybrid;
257 using real_type =
typename FieldTraits<field_type>::real_type;
259 real_type result = 0.0;
262 if constexpr (HasNaN<field_type>()) {
264 real_type nanTracker = 1.0;
265 using namespace Dune::Hybrid;
266 forEach(*
this, [&](
auto&& entry) {
267 real_type entryNorm = entry.infinity_norm();
268 result = max(entryNorm, result);
269 nanTracker += entryNorm;
272 result *= (nanTracker / nanTracker);
274 using namespace Dune::Hybrid;
275 forEach(*
this, [&](
auto&& entry) {
276 result = max(entry.infinity_norm(), result);
286 using namespace Dune::Hybrid;
288 using real_type =
typename FieldTraits<field_type>::real_type;
290 real_type result = 0.0;
293 if constexpr (HasNaN<field_type>()) {
295 real_type nanTracker = 1.0;
296 using namespace Dune::Hybrid;
297 forEach(*
this, [&](
auto&& entry) {
298 real_type entryNorm = entry.infinity_norm_real();
299 result = max(entryNorm, result);
300 nanTracker += entryNorm;
303 result *= (nanTracker / nanTracker);
305 using namespace Dune::Hybrid;
306 forEach(*
this, [&](
auto&& entry) {
307 result = max(entry.infinity_norm_real(), result);
317 template<
typename Ta>
319 using namespace Dune::Hybrid;
320 forEach(integralRange(Hybrid::size(*
this)), [&](
auto&& i) {
321 (*this)[i].axpy(a, y[i]);
331 template <
typename... Args>
333 using namespace Dune::Hybrid;
334 forEach(integralRange(Dune::Hybrid::size(v)), [&](
auto&& i) {
335 s <<
"\t(" << i <<
"):\n" << v[i] <<
"\n";
Simple iterative methods like Jacobi, Gauss-Seidel, SOR, SSOR, etc. in a generic way.
FieldTraits< Arg0 >::field_type field_type
Definition: multitypeblockvector.hh:41
FieldTraits< Arg0 >::real_type real_type
Definition: multitypeblockvector.hh:42
void operator=(const T &newval)
Assignment operator.
Definition: multitypeblockvector.hh:162
std::size_t size_type
Type used for vector sizes.
Definition: multitypeblockvector.hh:62
int count() const
Definition: multitypeblockvector.hh:107
static constexpr size_type N()
Number of elements.
Definition: multitypeblockvector.hh:95
static constexpr size_type size()
Return the number of non-zero vector entries.
Definition: multitypeblockvector.hh:88
std::tuple_element< index, TupleType >::type & operator[](const std::integral_constant< size_type, index > indexVariable)
Random-access operator.
Definition: multitypeblockvector.hh:142
FieldTraits< field_type >::real_type two_norm() const
Compute the Euclidean norm.
Definition: multitypeblockvector.hh:249
size_type dim() const
Number of scalar elements.
Definition: multitypeblockvector.hh:113
field_type dot(const type &newv) const
Definition: multitypeblockvector.hh:213
void operator*=(const T &w)
Multiplication with a scalar.
Definition: multitypeblockvector.hh:191
void axpy(const Ta &a, const type &y)
Axpy operation on this vector (*this += a * y)
Definition: multitypeblockvector.hh:318
void operator/=(const T &w)
Division by a scalar.
Definition: multitypeblockvector.hh:200
MultiTypeBlockVector< Args... > type
Definition: multitypeblockvector.hh:72
FieldTraits< field_type >::real_type infinity_norm_real() const
Compute the simplified maximum norm (uses 1-norm for complex values)
Definition: multitypeblockvector.hh:284
auto one_norm() const
Compute the 1-norm.
Definition: multitypeblockvector.hh:222
void operator-=(const type &newv)
Definition: multitypeblockvector.hh:181
field_type operator*(const type &newv) const
Definition: multitypeblockvector.hh:206
void operator+=(const type &newv)
Definition: multitypeblockvector.hh:171
FieldTraits< field_type >::real_type infinity_norm() const
Compute the maximum norm.
Definition: multitypeblockvector.hh:253
auto one_norm_real() const
Compute the simplified 1-norm (uses 1-norm also for complex values)
Definition: multitypeblockvector.hh:231
double field_type
The type used for scalars.
Definition: multitypeblockvector.hh:81
FieldTraits< field_type >::real_type two_norm2() const
Compute the squared Euclidean norm.
Definition: multitypeblockvector.hh:240
Definition: allocator.hh:9
std::ostream & operator<<(std::ostream &s, const BlockVector< K, A > &v)
Send BlockVector to an output stream.
Definition: bvector.hh:588
A Vector class to support different block types.
Definition: multitypeblockvector.hh:56