Skip to content
Snippets Groups Projects
Commit 61b68a36 authored by Christian Engwer's avatar Christian Engwer
Browse files

Merge branch 'feature/FS1665'

This fixes bug FS#1665

* we introduce limit the constructor of GMPField to match only those
  types which can initialize the underlying mpf_class
* we extend the FieldVector<K,1> constructors to accept all scalar
  types which are compatible to K
* we update DenseVector to use ADL to find functions like max(...),
  otherwise it would fail for GMPField
parents ec3aa65a 9d46150e
No related branches found
No related tags found
No related merge requests found
......@@ -42,7 +42,8 @@ namespace Dune {
template<class K>
inline typename FieldTraits<K>::real_type absreal (const K& k)
{
return std::abs(k);
using std::abs;
return abs(k);
}
/**
......@@ -52,7 +53,8 @@ namespace Dune {
template<class K>
inline typename FieldTraits<K>::real_type absreal (const std::complex<K>& c)
{
return std::abs(c.real()) + std::abs(c.imag());
using std::abs;
return abs(c.real()) + abs(c.imag());
}
/**
......@@ -84,7 +86,8 @@ namespace Dune {
{
static inline typename FieldTraits<K>::real_type sqrt (const K& k)
{
return std::sqrt(k);
using std::sqrt;
return sqrt(k);
}
};
......@@ -97,7 +100,8 @@ namespace Dune {
{
static inline typename FieldTraits<K>::real_type sqrt (const K& k)
{
return typename FieldTraits<K>::real_type(std::sqrt(double(k)));
using std::sqrt;
return typename FieldTraits<K>::real_type(sqrt(double(k)));
}
};
......@@ -564,9 +568,10 @@ namespace Dune {
//! one norm (sum over absolute values of entries)
typename FieldTraits<value_type>::real_type one_norm() const {
using std::abs;
typename FieldTraits<value_type>::real_type result( 0 );
for (size_type i=0; i<size(); i++)
result += std::abs((*this)[i]);
result += abs((*this)[i]);
return result;
}
......@@ -601,29 +606,35 @@ namespace Dune {
//! infinity norm (maximum of absolute values of entries)
typename FieldTraits<value_type>::real_type infinity_norm () const
{
using std::abs;
using std::max;
typedef typename FieldTraits<value_type>::real_type real_type;
if (size() == 0)
return 0.0;
ConstIterator it = begin();
typename FieldTraits<value_type>::real_type max = std::abs(*it);
real_type max_val = abs(*it);
for (it = it + 1; it != end(); ++it)
max = std::max(max, std::abs(*it));
max_val = max(max_val, real_type(abs(*it)));
return max;
return max_val;
}
//! simplified infinity norm (uses Manhattan norm for complex values)
typename FieldTraits<value_type>::real_type infinity_norm_real () const
{
using std::max;
if (size() == 0)
return 0.0;
ConstIterator it = begin();
typename FieldTraits<value_type>::real_type max = fvmeta::absreal(*it);
typename FieldTraits<value_type>::real_type max_val = fvmeta::absreal(*it);
for (it = it + 1; it != end(); ++it)
max = std::max(max, fvmeta::absreal(*it));
max_val = max(max_val, fvmeta::absreal(*it));
return max;
return max_val;
}
//===== sizes
......
......@@ -17,6 +17,8 @@
#include "typetraits.hh"
#include "exceptions.hh"
#include "array.hh"
#include "ftraits.hh"
#include "densevector.hh"
#include "unused.hh"
......@@ -226,7 +228,14 @@ namespace Dune {
{}
/** \brief Constructor with a given scalar */
FieldVector (const K& k) : _data(k) {}
template<typename T,
typename EnableIf = typename std::enable_if<
std::is_convertible<T, K>::value &&
! std::is_same<K, DenseVector<typename FieldTraits<T>::field_type>
>::value
>::type
>
FieldVector (const T& k) : _data(k) {}
//! Constructor making vector with identical coordinates
template<class C>
......@@ -243,7 +252,14 @@ namespace Dune {
{}
//! Assignment operator for scalar
inline FieldVector& operator= (const K& k)
template<typename T,
typename EnableIf = typename std::enable_if<
std::is_convertible<T, K>::value &&
! std::is_same<K, DenseVector<typename FieldTraits<T>::field_type>
>::value
>::type
>
inline FieldVector& operator= (const T& k)
{
_data = k;
return *this;
......
......@@ -8,6 +8,7 @@
*/
#include <iostream>
#include <string>
#if HAVE_GMP
......@@ -23,22 +24,34 @@ namespace Dune
typedef mpf_class Base;
public:
/** default constructor, iitialize to zero */
GMPField ()
: Base(0,precision)
{}
template< class T >
GMPField ( const T &v )
: Base( v,precision )
/** \brief initialize from a string
\note this is the only reliable way to initialize with higher values
*/
GMPField ( const char* str )
: Base(str,precision)
{}
/*
GMPField &operator=(const GMPField &other)
{
Base(*this) = Base(other);
return *this;
}
/** \brief initialize from a string
\note this is the only reliable way to initialize with higher values
*/
GMPField ( const std::string& str )
: Base(str,precision)
{}
/** \brief initialize from a compatible scalar type
*/
template< class T,
typename EnableIf = typename std::enable_if<
std::is_convertible<T, mpf_class>::value>::type
>
GMPField ( const T &v )
: Base( v,precision )
{}
// type conversion operators
operator double () const
......@@ -46,82 +59,8 @@ namespace Dune
return this->get_d();
}
operator float () const
{
return this->get_d();
}
};
template< unsigned int precision >
inline GMPField< precision >
operator+ ( const GMPField< precision > &a, const GMPField< precision > &b )
{
typedef mpf_class F;
return ((const F &)a + (const F &)b);
}
template< unsigned int precision >
inline GMPField< precision >
operator- ( const GMPField< precision > &a, const GMPField< precision > &b )
{
typedef mpf_class F;
return ((const F &)a - (const F &)b);
}
template< unsigned int precision >
inline GMPField< precision >
operator- ( const GMPField< precision > &a )
{
typedef mpf_class F;
return -((const F &)a);
}
template< unsigned int precision >
inline GMPField< precision >
operator* ( const GMPField< precision > &a, const GMPField< precision > &b )
{
typedef mpf_class F;
return ((const F &)a * (const F &)b);
}
template< unsigned int precision >
inline GMPField< precision >
operator/ ( const GMPField< precision > &a, const GMPField< precision > &b )
{
typedef mpf_class F;
return ((const F &)a / (const F &)b);
}
template< unsigned int precision >
inline std::ostream &
operator<< ( std::ostream &out, const GMPField< precision > &value )
{
return out << static_cast<const mpf_class&>(value);
}
}
namespace std
{
template< unsigned int precision >
inline Dune::GMPField< precision >
sqrt ( const Dune::GMPField< precision > &a )
{
return Dune::GMPField< precision >(sqrt(static_cast<const mpf_class&>(a)));
}
template< unsigned int precision >
inline Dune::GMPField< precision >
abs ( const Dune::GMPField< precision > &a )
{
return Dune::GMPField< precision >( abs( static_cast< const mpf_class & >( a ) ) );
}
}
#endif // HAVE_GMP
......
......@@ -7,6 +7,7 @@
#include <dune/common/exceptions.hh>
#include <dune/common/typetraits.hh>
#include <dune/common/classname.hh>
#include <dune/common/gmpfield.hh>
#include <iostream>
#include <complex>
#include <typeinfo>
......@@ -380,6 +381,17 @@ int main()
FieldVectorTest<double, 3>();
FieldVectorTest<int, 1>();
FieldVectorTest<double, 1>();
#if HAVE_GMP
// we skip the complex test and the int test, as these will be very hard to implement with GMPField
typedef Dune::GMPField<128u> ft;
FieldVectorMainTest<ft,ft,3>();
FieldVectorMainTest<ft,ft,2>();
FieldVectorMainTest<ft,ft,1>();
FieldVectorMainTest<ft,ft,0>();
ScalarOperatorTest<ft>();
ScalarOrderingTest<ft>();
DotProductTest<ft,3>();
#endif // HAVE_GMP
test_nan();
test_infinity_norms();
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment