5#ifndef DUNE_ISTL_MULTITYPEBLOCKVECTOR_HH
6#define DUNE_ISTL_MULTITYPEBLOCKVECTOR_HH
12#include <dune/common/dotproduct.hh>
13#include <dune/common/ftraits.hh>
14#include <dune/common/hybridutilities.hh>
15#include <dune/common/typetraits.hh>
16#include <dune/common/std/type_traits.hh>
22 template <
typename... Args >
23 class MultiTypeBlockVector;
41 template <
typename... Args>
45 using real_type =
typename FieldTraits<field_type>::real_type;
56 template <
typename... Args >
58 :
public std::tuple<Args...>
61 typedef std::tuple<Args...> TupleType;
70 using TupleType::TupleType;
86 static_assert (
sizeof...(Args) == 0 or not std::is_same_v<field_type, Std::nonesuch>,
87 "No std::common_type implemented for the given field_types of the Args. Please provide an implementation and check that a FieldTraits class is present for your type.");
97 return sizeof...(Args);
104 return sizeof...(Args);
113 [[deprecated(
"Use method 'N' instead")]]
116 return sizeof...(Args);
123 Hybrid::forEach(std::make_index_sequence<
N()>{},
124 [&](
auto i){result += std::get<i>(*this).dim();});
147 template<
size_type index >
148 typename std::tuple_element<index,TupleType>::type&
149 operator[] ([[maybe_unused]]
const std::integral_constant< size_type, index > indexVariable)
151 return std::get<index>(*
this);
159 template<
size_type index >
160 const typename std::tuple_element<index,TupleType>::type&
161 operator[] ([[maybe_unused]]
const std::integral_constant< size_type, index > indexVariable)
const
163 return std::get<index>(*
this);
170 Dune::Hybrid::forEach(*
this, [&](
auto&& entry) {
179 using namespace Dune::Hybrid;
180 forEach(integralRange(Hybrid::size(*
this)), [&](
auto&& i) {
181 (*this)[i] += newv[i];
189 using namespace Dune::Hybrid;
190 forEach(integralRange(Hybrid::size(*
this)), [&](
auto&& i) {
191 (*this)[i] -= newv[i];
197 std::enable_if_t< IsNumber<T>::value,
int> = 0>
199 Hybrid::forEach(*
this, [&](
auto&& entry) {
206 std::enable_if_t< IsNumber<T>::value,
int> = 0>
208 Hybrid::forEach(*
this, [&](
auto&& entry) {
214 using namespace Dune::Hybrid;
215 return accumulate(integralRange(Hybrid::size(*
this)),
field_type(0), [&](
auto&& a,
auto&& i) {
216 return a + (*this)[i]*newv[i];
221 using namespace Dune::Hybrid;
222 return accumulate(integralRange(Hybrid::size(*
this)),
field_type(0), [&](
auto&& a,
auto&& i) {
223 return a + (*this)[i].dot(newv[i]);
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();
239 using namespace Dune::Hybrid;
240 return accumulate(*
this,
typename FieldTraits<field_type>::real_type(0), [&](
auto&& a,
auto&& entry) {
241 return a + entry.one_norm_real();
247 typename FieldTraits<field_type>::real_type
two_norm2()
const {
248 using namespace Dune::Hybrid;
249 return accumulate(*
this,
typename FieldTraits<field_type>::real_type(0), [&](
auto&& a,
auto&& entry) {
250 return a + entry.two_norm2();
262 using namespace Dune::Hybrid;
264 using real_type =
typename FieldTraits<field_type>::real_type;
266 real_type result = 0.0;
269 if constexpr (HasNaN<field_type>()) {
271 real_type nanTracker = 1.0;
272 using namespace Dune::Hybrid;
273 forEach(*
this, [&](
auto&& entry) {
274 real_type entryNorm = entry.infinity_norm();
275 result = max(entryNorm, result);
276 nanTracker += entryNorm;
279 result *= (nanTracker / nanTracker);
281 using namespace Dune::Hybrid;
282 forEach(*
this, [&](
auto&& entry) {
283 result = max(entry.infinity_norm(), result);
293 using namespace Dune::Hybrid;
295 using real_type =
typename FieldTraits<field_type>::real_type;
297 real_type result = 0.0;
300 if constexpr (HasNaN<field_type>()) {
302 real_type nanTracker = 1.0;
303 using namespace Dune::Hybrid;
304 forEach(*
this, [&](
auto&& entry) {
305 real_type entryNorm = entry.infinity_norm_real();
306 result = max(entryNorm, result);
307 nanTracker += entryNorm;
310 result *= (nanTracker / nanTracker);
312 using namespace Dune::Hybrid;
313 forEach(*
this, [&](
auto&& entry) {
314 result = max(entry.infinity_norm_real(), result);
324 template<
typename Ta>
326 using namespace Dune::Hybrid;
327 forEach(integralRange(Hybrid::size(*
this)), [&](
auto&& i) {
328 (*this)[i].axpy(a, y[i]);
338 template <
typename... Args>
340 using namespace Dune::Hybrid;
341 forEach(integralRange(Dune::Hybrid::size(v)), [&](
auto&& i) {
342 s <<
"\t(" << i <<
"):\n" << v[i] <<
"\n";
355 template <
size_t i,
typename... Args>
356 struct tuple_element<i,
Dune::MultiTypeBlockVector<Args...> >
358 using type =
typename std::tuple_element<i, std::tuple<Args...> >::type;
Simple iterative methods like Jacobi, Gauss-Seidel, SOR, SSOR, etc. in a generic way.
void operator=(const T &newval)
Assignment operator.
Definition: multitypeblockvector.hh:169
std::size_t size_type
Type used for vector sizes.
Definition: multitypeblockvector.hh:65
typename FieldTraits< field_type >::real_type real_type
Definition: multitypeblockvector.hh:45
int count() const
Definition: multitypeblockvector.hh:114
static constexpr size_type N()
Number of elements.
Definition: multitypeblockvector.hh:102
Std::detected_t< std::common_type_t, typename FieldTraits< std::decay_t< Args > >::field_type... > field_type
The type used for scalars.
Definition: multitypeblockvector.hh:82
static constexpr size_type size()
Return the number of non-zero vector entries.
Definition: multitypeblockvector.hh:95
std::tuple_element< index, TupleType >::type & operator[](const std::integral_constant< size_type, index > indexVariable)
Random-access operator.
Definition: multitypeblockvector.hh:149
FieldTraits< field_type >::real_type two_norm() const
Compute the Euclidean norm.
Definition: multitypeblockvector.hh:256
typename MultiTypeBlockVector< Args... >::field_type field_type
Definition: multitypeblockvector.hh:44
size_type dim() const
Number of scalar elements.
Definition: multitypeblockvector.hh:120
field_type dot(const type &newv) const
Definition: multitypeblockvector.hh:220
void operator*=(const T &w)
Multiplication with a scalar.
Definition: multitypeblockvector.hh:198
void axpy(const Ta &a, const type &y)
Axpy operation on this vector (*this += a * y)
Definition: multitypeblockvector.hh:325
void operator/=(const T &w)
Division by a scalar.
Definition: multitypeblockvector.hh:207
MultiTypeBlockVector< Args... > type
Definition: multitypeblockvector.hh:75
FieldTraits< field_type >::real_type infinity_norm_real() const
Compute the simplified maximum norm (uses 1-norm for complex values)
Definition: multitypeblockvector.hh:291
auto one_norm() const
Compute the 1-norm.
Definition: multitypeblockvector.hh:229
void operator-=(const type &newv)
Definition: multitypeblockvector.hh:188
field_type operator*(const type &newv) const
Definition: multitypeblockvector.hh:213
void operator+=(const type &newv)
Definition: multitypeblockvector.hh:178
typename std::tuple_element< i, std::tuple< Args... > >::type type
Definition: multitypeblockvector.hh:358
FieldTraits< field_type >::real_type infinity_norm() const
Compute the maximum norm.
Definition: multitypeblockvector.hh:260
auto one_norm_real() const
Compute the simplified 1-norm (uses 1-norm also for complex values)
Definition: multitypeblockvector.hh:238
FieldTraits< field_type >::real_type two_norm2() const
Compute the squared Euclidean norm.
Definition: multitypeblockvector.hh:247
Definition: allocator.hh:11
std::ostream & operator<<(std::ostream &s, const BlockVector< K, A > &v)
Send BlockVector to an output stream.
Definition: bvector.hh:590
A Vector class to support different block types.
Definition: multitypeblockvector.hh:59