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> 18 template <
typename... Args >
37 template <
typename Arg0,
typename... Args>
40 typedef typename FieldTraits<Arg0>::field_type
field_type;
41 typedef typename FieldTraits<Arg0>::real_type
real_type;
52 template <
typename... Args >
54 :
public std::tuple<Args...>
57 typedef std::tuple<Args...> tupleType;
75 static constexpr std::size_t
size()
77 return sizeof...(Args);
85 return sizeof...(Args);
106 template< std::
size_t index >
107 typename std::tuple_element<index,tupleType>::type&
108 operator[] (
const std::integral_constant< std::size_t, index > indexVariable )
110 DUNE_UNUSED_PARAMETER(indexVariable);
111 return std::get<index>(*this);
119 template< std::
size_t index >
120 const typename std::tuple_element<index,tupleType>::type&
121 operator[] (
const std::integral_constant< std::size_t, index > indexVariable )
const 123 DUNE_UNUSED_PARAMETER(indexVariable);
124 return std::get<index>(*this);
130 void operator= (
const T& newval) {
131 Dune::Hybrid::forEach(*
this, [&](
auto&& entry) {
139 void operator+= (
const type& newv) {
140 using namespace Dune::Hybrid;
141 forEach(integralRange(Hybrid::size(*
this)), [&](
auto&& i) {
142 (*this)[i] += newv[i];
149 void operator-= (
const type& newv) {
150 using namespace Dune::Hybrid;
151 forEach(integralRange(Hybrid::size(*
this)), [&](
auto&& i) {
152 (*this)[i] -= newv[i];
167 void operator*= (
const int& w) {
168 Dune::Hybrid::forEach(*
this, [&](
auto&& entry) {
173 void operator*= (
const float& w) {
174 Dune::Hybrid::forEach(*
this, [&](
auto&& entry) {
179 void operator*= (
const double& w) {
180 Dune::Hybrid::forEach(*
this, [&](
auto&& entry) {
185 field_type operator* (
const type& newv)
const {
186 using namespace Dune::Hybrid;
187 return accumulate(integralRange(Hybrid::size(*
this)), field_type(0), [&](
auto&& a,
auto&& i) {
188 return a + (*this)[i]*newv[i];
193 using namespace Dune::Hybrid;
194 return accumulate(integralRange(Hybrid::size(*
this)), field_type(0), [&](
auto&& a,
auto&& i) {
195 return a + (*this)[i].dot(newv[i]);
201 typename FieldTraits<field_type>::real_type
two_norm2()
const {
202 using namespace Dune::Hybrid;
203 return accumulate(*
this,
typename FieldTraits<field_type>::real_type(0), [&](
auto&& a,
auto&& entry) {
204 return a + entry.two_norm2();
210 typename FieldTraits<field_type>::real_type
two_norm()
const {
return sqrt(this->two_norm2());}
216 using namespace Dune::Hybrid;
218 using real_type =
typename FieldTraits<field_type>::real_type;
220 real_type result = 0.0;
223 ifElse(has_nan<field_type>(), [&](
auto&&
id) {
225 real_type nanTracker = 1.0;
226 forEach(*
this, [&](
auto&& entry) {
227 real_type entryNorm = entry.infinity_norm();
228 result = max(entryNorm, result);
229 nanTracker += entryNorm;
232 result *= (nanTracker / nanTracker);
234 forEach(*
this, [&](
auto&& entry) {
235 result = max(entry.infinity_norm(), result);
245 template<
typename Ta>
247 using namespace Dune::Hybrid;
248 forEach(integralRange(Hybrid::size(*
this)), [&](
auto&& i) {
249 (*this)[i].axpy(a, y[i]);
259 template <
typename... Args>
261 using namespace Dune::Hybrid;
262 forEach(integralRange(Dune::Hybrid::size(v)), [&](
auto&& i) {
263 s <<
"\t(" << i <<
"):\n" << v[i] <<
"\n";
Definition: allocator.hh:7
void axpy(const Ta &a, const type &y)
Axpy operation on this vector (*this += a * y)
Definition: multitypeblockvector.hh:246
A Vector class to support different block types.
Definition: multitypeblockvector.hh:19
FieldTraits< field_type >::real_type two_norm() const
Compute the Euclidean norm.
Definition: multitypeblockvector.hh:210
Simple iterative methods like Jacobi, Gauss-Seidel, SOR, SSOR, etc. in a generic way.
std::ostream & operator<<(std::ostream &s, const BlockVector< K, A > &v)
Send BlockVector to an output stream.
Definition: bvector.hh:620
field_type dot(const type &newv) const
Definition: multitypeblockvector.hh:192
int count()
Definition: multitypeblockvector.hh:83
FieldTraits< Arg0 >::field_type field_type
Definition: multitypeblockvector.hh:40
FieldTraits< field_type >::real_type two_norm2() const
Compute the squared Euclidean norm.
Definition: multitypeblockvector.hh:201
double field_type
The type used for scalars.
Definition: multitypeblockvector.hh:72
static constexpr std::size_t size()
Return the number of vector entries.
Definition: multitypeblockvector.hh:75
FieldTraits< Arg0 >::real_type real_type
Definition: multitypeblockvector.hh:41
MultiTypeBlockVector< Args... > type
Definition: multitypeblockvector.hh:63
FieldTraits< field_type >::real_type infinity_norm() const
Compute the maximum norm.
Definition: multitypeblockvector.hh:214