refining numeric_limits for quad-double

Author: Ravenwater

include/universal/common/number_traits_reports.hpp

@@ -21,7 +21,7 @@ void numberTraits(std::ostream& ostr) {
 	ostr << "min                " << setw(ColumnWidth) << numeric_limits<Scalar>::min() << '\n';
 	ostr << "max                " << setw(ColumnWidth) << numeric_limits<Scalar>::max() << '\n';
 	ostr << "lowest             " << setw(ColumnWidth) << numeric_limits<Scalar>::lowest() << '\n';
-	ostr << "epsilon (1+1ULP-1) " << setw(ColumnWidth) << numeric_limits<Scalar>::epsilon() << '\n';
+	ostr << "epsilon ==ulp(1.0) " << setw(ColumnWidth) << numeric_limits<Scalar>::epsilon() << '\n';
 	ostr << "round_error        " << setw(ColumnWidth) << numeric_limits<Scalar>::round_error() << '\n';
 	ostr << "smallest value     " << setw(ColumnWidth) << numeric_limits<Scalar>::denorm_min() << '\n';
 	ostr << "infinity           " << setw(ColumnWidth) << numeric_limits<Scalar>::infinity() << '\n';

include/universal/number/qd/manipulators.hpp

@@ -19,23 +19,13 @@ namespace sw { namespace universal {
 	}
 
 	// generate a binary, color-coded representation of the quad-double
-	inline std::string color_print(qd const& r, bool nibbleMarker = false) {
-		//constexpr unsigned es = 11;
-		//constexpr unsigned fbits = 106;
+	inline std::string color_print(const qd& r, bool nibbleMarker = false) {
 		std::stringstream s;
-
-		/*
-		Color red(ColorCode::FG_RED);
-		Color yellow(ColorCode::FG_YELLOW);
-		Color blue(ColorCode::FG_BLUE);
-		Color magenta(ColorCode::FG_MAGENTA);
-		Color cyan(ColorCode::FG_CYAN);
-		Color white(ColorCode::FG_WHITE);
-		Color def(ColorCode::FG_DEFAULT);
-		*/
 		for (int i = 0; i < 4; ++i) {
+			std::string label = "x[" + std::to_string(i) + "]";
+			s << std::setw(20) << label << " : ";
 			s << color_print(r[i], nibbleMarker);
-			if (i < 3) s << ", ";
+			if (i < 3) s << '\n';
 		}
 		return s.str();
 	}

include/universal/number/qd/numeric_limits.hpp

@@ -15,38 +15,32 @@ class numeric_limits< sw::universal::qd > {
 	using QuadDouble = sw::universal::qd;
 	static constexpr bool is_specialized = true;
 	static constexpr QuadDouble min() { // return minimum value
-		// return QuadDouble(sw::universal::SpecificValue::minpos);
 		return QuadDouble(radix * (numeric_limits< double >::min() / numeric_limits< double >::epsilon()));
 	} 
 	static constexpr QuadDouble max() { // return maximum value
-		//return QuadDouble(sw::universal::SpecificValue::maxpos);
-		return QuadDouble(numeric_limits< double >::max());
+		return QuadDouble(sw::universal::SpecificValue::maxpos);
 	} 
 	static constexpr QuadDouble lowest() { // return most negative value
-		//return QuadDouble(sw::universal::SpecificValue::maxneg);
-		return (-(max)());
+		return QuadDouble(sw::universal::SpecificValue::maxneg);
 	} 
 	static constexpr QuadDouble epsilon() { // return smallest effective increment from 1.0
-		return numeric_limits< double >::epsilon() * numeric_limits< double >::epsilon() / radix;
+		constexpr double epsilon{ std::numeric_limits< double >::epsilon() };
+		return  (((epsilon * epsilon) * epsilon) * epsilon) * 0.5; // / QuadDouble(radix);
 	}
 	static constexpr QuadDouble round_error() { // return largest rounding error
 		return QuadDouble(1.0 / radix);
 	}
 	static constexpr QuadDouble denorm_min() {  // return minimum denormalized value
-		// return QuadDouble(sw::universal::SpecificValue::minpos);
-		return 0.0;
+		return QuadDouble(std::numeric_limits<double>::denorm_min());
 	}
 	static constexpr QuadDouble infinity() { // return positive infinity
 		return QuadDouble(sw::universal::SpecificValue::infpos);
-		//return numeric_limits< double >::infinity();
 	}
 	static constexpr QuadDouble quiet_NaN() { // return non-signaling NaN
 		return QuadDouble(sw::universal::SpecificValue::qnan);
-		//return numeric_limits< double >::quiet_NaN();
 	}
 	static constexpr QuadDouble signaling_NaN() { // return signaling NaN
 		return QuadDouble(sw::universal::SpecificValue::snan);
-		//return numeric_limits< double >::signaling_NaN();
 	}
 
 	static constexpr int  digits                   = 2 * std::numeric_limits<double>::digits;

include/universal/number/qd/qd_impl.hpp

@@ -34,7 +34,9 @@ qd operator+(const qd&, const qd&);
 qd operator-(const qd&, const qd&);
 qd operator*(const qd&, const qd&);
 qd operator/(const qd&, const qd&);
+std::ostream& operator<<(std::ostream&, const qd&);
 qd pown(const qd&, int);
+qd frexp(const qd&, int*);
 
 // qd is an unevaluated quadruple of IEEE-754 doubles that provides a (1,11,212) floating-point triple
 class qd {
@@ -300,14 +302,15 @@ class qd {
 
 	// create specific number system values of interest
 	constexpr qd& maxpos() noexcept {
-		x[0] = 0.0; 
+		x[0] = std::numeric_limits<double>::max(); 
 		x[1] = 0.0; 
 		x[2] = 0.0; 
 		x[3] = 0.0;
 		return *this;
 	}
+	// smallest positive normal number
 	constexpr qd& minpos() noexcept {
-		x[0] = 0.0;
+		x[0] = std::numeric_limits<double>::min();
 		x[1] = 0.0;
 		x[2] = 0.0;
 		x[3] = 0.0;
@@ -320,15 +323,16 @@ class qd {
 		x[3] = 0.0;
 		return *this;
 	}
+	// smallest negative normal number
 	constexpr qd& minneg() noexcept {
-		x[0] = 0.0;
+		x[0] = -std::numeric_limits<double>::min();
 		x[1] = 0.0;
 		x[2] = 0.0;
 		x[3] = 0.0;
 		return *this;
 	}
 	constexpr qd& maxneg() noexcept {
-		x[0] = 0.0;
+		x[0] = std::numeric_limits<double>::lowest();
 		x[1] = 0.0;
 		x[2] = 0.0;
 		x[3] = 0.0;
@@ -1105,6 +1109,16 @@ inline std::string to_quad(const qd& v, int precision = 17) {
 	return s.str();
 }
 
+inline std::string to_triple(const qd& v, int precision = 17) {
+	std::stringstream s;
+	bool isneg = v.isneg();
+	int scale = v.scale();
+	int exponent;
+	qd fraction = frexp(v, &exponent);
+	s << '(' << (isneg ? '1' : '0') << ", " << scale << ", " << std::setprecision(precision) << fraction << ')';
+	return s.str();
+}
+
 inline std::string to_binary(const qd& number, bool bNibbleMarker = false) {
 	std::stringstream s;
 	constexpr int nrLimbs = 4;