DoubleDouble.hpp File Reference

A double-double scalar: the unevaluated sum hi + lo of two doubles, carrying about 32 significant digits. More...

Classes
class	mundy::DoubleDouble
	A double-double scalar hi + lo, normalized so that |lo| <= ulp(hi) / 2. More...
struct	mundy::is_passive_scalar< DoubleDouble >
struct	mundy::NumTraits< DoubleDouble >
	Numeric traits of DoubleDouble, in DoubleDouble. More...
struct	mundy::ScalarBinaryOpTraits< DoubleDouble, A, Op >
	DoubleDouble with an arithmetic primitive yields DoubleDouble. More...
struct	mundy::ScalarBinaryOpTraits< A, DoubleDouble, Op >
	An arithmetic primitive with DoubleDouble yields DoubleDouble. More...
struct	Kokkos::reduction_identity< mundy::DoubleDouble >
	Identities for Kokkos reductions (Sum, Prod, Max, Min) over DoubleDouble. More...

Namespaces
namespace	mundy

Functions
Arithmetic
constexpr DoubleDouble	mundy::operator+ (const DoubleDouble &a)
constexpr DoubleDouble	mundy::operator- (const DoubleDouble &a)
constexpr DoubleDouble	mundy::operator+ (const DoubleDouble &a, const DoubleDouble &b)
	a + b, within 2 units of 2^-104 relative even under cancellation (QD's IEEE-style addition).
constexpr DoubleDouble	mundy::operator- (const DoubleDouble &a, const DoubleDouble &b)
constexpr DoubleDouble	mundy::operator* (const DoubleDouble &a, const DoubleDouble &b)
constexpr DoubleDouble	mundy::operator/ (const DoubleDouble &a, const DoubleDouble &b)
	a / b: the double quotient, corrected twice by the remainder it leaves (QD's accurate division).
constexpr DoubleDouble	mundy::operator/ (const DoubleDouble &a, double b)
	a / b for a double divisor: one correction of the double quotient, cheaper than the general division.
Comparison (lexicographic on (hi, lo), which orders normalized values; NaN compares false)
constexpr bool	mundy::operator== (const DoubleDouble &a, const DoubleDouble &b)
constexpr bool	mundy::operator!= (const DoubleDouble &a, const DoubleDouble &b)
constexpr bool	mundy::operator< (const DoubleDouble &a, const DoubleDouble &b)
constexpr bool	mundy::operator<= (const DoubleDouble &a, const DoubleDouble &b)
constexpr bool	mundy::operator> (const DoubleDouble &a, const DoubleDouble &b)
constexpr bool	mundy::operator>= (const DoubleDouble &a, const DoubleDouble &b)
Sign, rounding, and selection
constexpr DoubleDouble	mundy::abs (const DoubleDouble &a)
constexpr DoubleDouble	mundy::copysign (const DoubleDouble &a, const DoubleDouble &sign)
constexpr DoubleDouble	mundy::min (const DoubleDouble &a, const DoubleDouble &b)
constexpr DoubleDouble	mundy::max (const DoubleDouble &a, const DoubleDouble &b)
DoubleDouble	mundy::floor (const DoubleDouble &a)
	Round down. If hi is already an integer the fraction lives in lo, so lo is rounded instead.
DoubleDouble	mundy::ceil (const DoubleDouble &a)
	Round up. If hi is already an integer the fraction lives in lo, so lo is rounded instead.
DoubleDouble	mundy::round (const DoubleDouble &a)
	Round to nearest, halfway cases away from zero (like std::round).
Roots
constexpr DoubleDouble	mundy::sqrt (const DoubleDouble &a)
	Constexpr square root: one Newton step x + (a - x^2) / (2x) from the double root x.
DoubleDouble	mundy::cbrt (const DoubleDouble &a)
	Cube root: one Newton step x + (a - x^3) / (3x^2) from the double root x.
Exponentials and logarithms
DoubleDouble	mundy::exp (const DoubleDouble &a)
	e^a.
DoubleDouble	mundy::log (const DoubleDouble &a)
	Natural logarithm.
DoubleDouble	mundy::log10 (const DoubleDouble &a)
DoubleDouble	mundy::pow (const DoubleDouble &a, const DoubleDouble &b)
	a^b. Integer exponents use repeated squaring (any sign of a); others use exp(b log a), which needs a > 0 and whose relative error grows like |b log a| units of 2^-104.
Trigonometric functions
DoubleDouble	mundy::sin (const DoubleDouble &a)
DoubleDouble	mundy::cos (const DoubleDouble &a)
DoubleDouble	mundy::tan (const DoubleDouble &a)
DoubleDouble	mundy::atan2 (const DoubleDouble &y, const DoubleDouble &x)
	The angle of the point (x, y), in (-pi, pi].
DoubleDouble	mundy::atan (const DoubleDouble &a)
DoubleDouble	mundy::asin (const DoubleDouble &a)
	Arcsine, through atan2(a, sqrt((1 - a)(1 + a))); NaN for |a| > 1.
DoubleDouble	mundy::acos (const DoubleDouble &a)
	Arccosine, through atan2(sqrt((1 - a)(1 + a)), a); NaN for |a| > 1.
The MundyMath custom-scalar contract
std::ostream &	mundy::operator<< (std::ostream &os, const DoubleDouble &a)
	Write the exact pair as "hi + lo", at the stream's precision.

Detailed Description

DoubleDouble satisfies the MundyMath custom-scalar contract (NumTraits, ScalarBinaryOpTraits, and an ADL overload for every function in cmath.hpp's dispatch set), so it stands in for double in Scalar, Vector, Matrix, Quaternion, AutoDiffScalar, and Kokkos reductions, on the host and the device. Its arithmetic is constexpr.

Accuracy: arithmetic, the roots, exp, log, and the trigonometric functions are within a few units of 2^-104 relative; pow's error grows like |b log a| such units. sin, cos, and tan fall back to double accuracy beyond |x| = 2^50, where their argument reduction ends. Values must stay above NumTraits::norm_min() (about 2e-292) for lo to be normal.

zmort.hpp and the KokkosBlas/Belos solver backends do not support DoubleDouble.

The algorithms follow Hida, Li & Bailey's QD library.