Skip to content

Commit

Permalink
dd code hygiene, qd api test progression
Browse files Browse the repository at this point in the history
  • Loading branch information
@Ravenwater
Ravenwater committed Aug 26, 2024
1 parent 1313975 commit 1d9821c
Show file tree
Hide file tree
Showing 5 changed files with 136 additions and 56 deletions.
3 changes: 2 additions & 1 deletion include/universal/internal/bitblock/bitblock.hpp
View file
Original file line number Diff line number Diff line change
@@ -1,7 +1,8 @@
#pragma once
// bitblock.hpp : bitblock class
//
// Copyright (C) 2017-2023 Stillwater Supercomputing, Inc.
// Copyright (C) 2017 Stillwater Supercomputing, Inc.
// SPDX-License-Identifier: MIT
//
// This file is part of the universal numbers project, which is released under an MIT Open Source license.
#include <cstdint>
Expand Down
69 changes: 36 additions & 33 deletions include/universal/number/dd/dd_impl.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -494,11 +494,11 @@ class dd {
std::vector<char> t;

if (fixed) {
t.resize(nrDigitsForFixedFormat+1);
t.resize(static_cast<size_t>(nrDigitsForFixedFormat+1));
to_digits(t, e, nrDigitsForFixedFormat);
}
else {
t.resize(nrDigits+1);
t.resize(static_cast<size_t>(nrDigits+1));
to_digits(t, e, nrDigits);
}

Expand All @@ -508,30 +508,32 @@ class dd {

if (integerDigits > 0) {
int i;
for (i = 0; i < integerDigits; ++i) s += t[i];
for (i = 0; i < integerDigits; ++i) s += t[static_cast<unsigned>(i)];
if (precision > 0) {
s += '.';
for (int j = 0; j < precision; ++j, ++i) s += t[i];
for (int j = 0; j < precision; ++j, ++i) s += t[static_cast<unsigned>(i)];
}
}
else {
s += "0.";
if (integerDigits < 0) s.append(static_cast<size_t>(-integerDigits), '0');
for (int i = 0; i < nrDigits; ++i) s += t[i];
for (int i = 0; i < nrDigits; ++i) s += t[static_cast<unsigned>(i)];
}
}
else {
s += t[0];
s += t[0ull];
if (precision > 0) s += '.';

for (int i = 1; i <= precision; ++i)
s += t[i];
s += t[static_cast<unsigned>(i)];

}
}
}

// trap for improper offset with large values
// TBD: this is seriously broken and needs a redesign
//
// fix for improper offset with large values and small values
// without this trap, output of values of the for 10^j - 1 fail for j > 28
// and are output with the point in the wrong place, leading to a significant error
if (fixed && (precision > 0)) {
Expand All @@ -546,17 +548,18 @@ class dd {
for (std::string::size_type i = 1; i < s.length(); ++i) {
if (s[i] == '.') {
s[i] = s[i - 1];
s[i - 1] = '.';
s[i - 1] = '.'; // this will destroy the leading 0 when s[i==1] == '.';
break;
}
}
// BUG: the loop above, in particular s[i-1] = '.', destroys the leading 0
// in the fixed point representation if the point is located at i = 1;
// it also breaks the precision request as it adds a new digit to the fixed representation

from_string = atof(s.c_str());
// if this ratio is large, then the string has not been fixed
if (std::fabs(from_string / hi) > 3.0) {
std::cerr << "re-rounding unsuccessful in large number fixed point trap\n";
std::cerr << "re-rounding unsuccessful in fixed point fix\n";
}
}
}
Expand Down Expand Up @@ -670,22 +673,22 @@ class dd {
int nrDigits = precision;
// round decimal string and propagate carry
int lastDigit = nrDigits - 1;
if (s[lastDigit] >= '5') {
if (s[static_cast<unsigned>(lastDigit)] >= '5') {
if constexpr(bTraceDecimalRounding) std::cout << "need to round\n";
int i = nrDigits - 2;
s[i]++;
while (i > 0 && s[i] > '9') {
s[i] -= 10;
s[--i]++;
s[static_cast<unsigned>(i)]++;
while (i > 0 && s[static_cast<unsigned>(i)] > '9') {
s[static_cast<unsigned>(i)] -= 10;
s[static_cast<unsigned>(--i)]++;
}
}

// if first digit is 10, shift everything.
if (s[0] > '9') {
if constexpr(bTraceDecimalRounding) std::cout << "shift right to handle overflow\n";
for (int i = precision; i >= 2; --i) s[i] = s[i - 1];
s[0] = '1';
s[1] = '0';
for (int i = precision; i >= 2; --i) s[static_cast<unsigned>(i)] = s[static_cast<unsigned>(i - 1)];
s[0u] = '1';
s[1u] = '0';

(*decimalPoint)++; // increment decimal point
++precision;
Expand Down Expand Up @@ -722,7 +725,7 @@ class dd {

if (iszero()) {
exponent = 0;
for (int i = 0; i < precision; ++i) s[i] = '0';
for (int i = 0; i < precision; ++i) s[static_cast<unsigned>(i)] = '0';
return;
}

Expand Down Expand Up @@ -781,20 +784,20 @@ class dd {
r -= mostSignificantDigit;
r *= 10.0;

s[i] = static_cast<char>(mostSignificantDigit + '0');
s[static_cast<unsigned>(i)] = static_cast<char>(mostSignificantDigit + '0');
if constexpr (bTraceDecimalConversion) std::cout << "to_digits digit[" << i << "] : " << s << '\n';
}

// Fix out of range digits
for (int i = nrDigits - 1; i > 0; --i) {
if (s[i] < '0') {
s[i - 1]--;
s[i] += 10;
if (s[static_cast<unsigned>(i)] < '0') {
s[static_cast<unsigned>(i - 1)]--;
s[static_cast<unsigned>(i)] += 10;
}
else {
if (s[i] > '9') {
s[i - 1]++;
s[i] -= 10;
if (s[static_cast<unsigned>(i)] > '9') {
s[static_cast<unsigned>(i - 1)]++;
s[static_cast<unsigned>(i)] -= 10;
}
}
}
Expand All @@ -806,26 +809,26 @@ class dd {

// Round and propagate carry
int lastDigit = nrDigits - 1;
if (s[lastDigit] >= '5') {
if (s[static_cast<unsigned>(lastDigit)] >= '5') {
int i = nrDigits - 2;
s[i]++;
while (i > 0 && s[i] > '9') {
s[i] -= 10;
s[--i]++;
s[static_cast<unsigned>(i)]++;
while (i > 0 && s[static_cast<unsigned>(i)] > '9') {
s[static_cast<unsigned>(i)] -= 10;
s[static_cast<unsigned>(--i)]++;
}
}

// If first digit is 10, shift left and increment exponent
if (s[0] > '9') {
++e;
for (int i = precision; i >= 2; --i) {
s[i] = s[i - 1];
s[static_cast<unsigned>(i)] = s[static_cast<unsigned>(i - 1)];
}
s[0] = '1';
s[1] = '0';
}

s[precision] = 0; // termination null
s[static_cast<unsigned>(precision)] = 0; // termination null
exponent = e;
}

Expand Down
1 change: 0 additions & 1 deletion include/universal/number/dd/manipulators.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,6 @@
// This file is part of the universal numbers project, which is released under an MIT Open Source license.
#include <string>
#include <iomanip>
#include <universal/internal/blockbinary/blockbinary.hpp>
#include <universal/number/dd/dd_fwd.hpp>
#include <universal/native/manipulators.hpp>
// pull in the color printing for shells utility
Expand Down
23 changes: 12 additions & 11 deletions static/dd/api/api.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -78,7 +78,7 @@ namespace sw {
int expOfFirstLimb = scale(a);
std::cout << to_binary(expOfFirstLimb) << " : " << expOfFirstLimb << '\n';
// second limb exponent
int expOfSecondLimb = expOfFirstLimb - std::log10(1ull << 53);
int expOfSecondLimb = expOfFirstLimb - 15; // floor(log10(2^53)) = floor(15.9245) = 15
std::cout << "exponent of the first limb : " << expOfFirstLimb << '\n';
std::cout << "exponent of the second limb : " << expOfSecondLimb << '\n';
// construct the second limb
Expand Down Expand Up @@ -115,7 +115,8 @@ namespace sw {
dd u = ulp(from);
std::cout << "ulp(" << std::scientific << std::setprecision(0) << from
<< ") gives "
<< std::fixed << std::setprecision(6) << u
// << std::fixed << std::setprecision(6) << u
<< " : " << std::setprecision(6) << std::defaultfloat << u
<< '\n';
}
}
Expand Down Expand Up @@ -211,8 +212,6 @@ try {
++d;
if (isdenorm(d)) std::cout << d << " is a subnormal number\n";
}


}

// helper api
Expand Down Expand Up @@ -268,11 +267,11 @@ try {
double x1{ std::pow(2.0, -53.0) };

// now let's walk that bit down to the ULP
int precisionForRange = 16;
unsigned precisionForRange = 16;
std::cout << std::setprecision(precisionForRange);
x0 = 1.0;
dd a(x0, x1);
std::cout << centered("double-double", precisionForRange + 6) << " : ";
std::cout << centered("double-double", precisionForRange + 6u) << " : ";
std::cout << centered("binary form of x0", 68) << " : ";
std::cout << centered("real value of x0", 15) << '\n';
std::cout << a << " : " << to_binary(x0) << " : " << x0 << '\n';
Expand All @@ -286,7 +285,7 @@ try {
x0 = 1.0;
precisionForRange = 32;
std::cout << std::setprecision(precisionForRange);
std::cout << centered("double-double", precisionForRange + 6) << " : ";
std::cout << centered("double-double", precisionForRange + 6u) << " : ";
std::cout << centered("binary form of x1", 68) << " : ";
std::cout << centered("real value of x1", 15) << '\n';
for (int i = 0; i < 54; ++i) {
Expand All @@ -298,7 +297,7 @@ try {
std::cout << "\nvalue and binary pattern of the double-double\n";
precisionForRange = 32;
std::cout << std::setprecision(precisionForRange);
std::cout << centered("double-double", precisionForRange + 6) << " : ";
std::cout << centered("double-double", precisionForRange + 6u) << " : ";
std::cout << centered("binary form of double-double", 110) << '\n';
for (int i = 0; i < 54; ++i) {
x1 = (std::pow(2.0, -53.0 - double(i)));
Expand All @@ -323,16 +322,18 @@ try {
std::cout << std::setprecision(32);
a.maxpos();
std::cout << "maxpos double-double : " << to_binary(a) << " : " << a << " : " << scale(a) << '\n';
a.setbits(0x0080); // positive min normal
std::cout << "minnorm double-double : " << to_binary(a) << " : " << a << " : " << scale(a) << '\n';
a.minpos();
std::cout << "minpos double-double : " << to_binary(a) << " : " << a << " : " << scale(a) << '\n';
a = std::numeric_limits<dd>::denorm_min();
std::cout << "smallest double-double: " << to_binary(a) << " : " << a << " : " << scale(a) << '\n';
a.zero();
std::cout << "zero : " << to_binary(a) << " : " << a << " : " << scale(a) << '\n';
a.minneg();
std::cout << "minneg double-double : " << to_binary(a) << " : " << a << " : " << scale(a) << '\n';
a.maxneg();
std::cout << "maxneg double-double : " << to_binary(a) << " : " << a << " : " << scale(a) << '\n';

std::cout << "Notice that minpos is the smallest normal number, not the smallest number, which is a denorm\n";
std::cout << std::setprecision(defaultPrecision);
std::cout << "---\n";
}
Expand Down Expand Up @@ -486,7 +487,7 @@ try {

std::cout << "---------- Unit in the Last Place --------+\n";
{
ulp_progression("\nULP progression for dd:\n", dd(10.0));
ulp_progression("\nULP progression for dd:\n", dd(10.0e01));

for (int i = -5; i < 6; ++i) {
dd a(std::pow(2.0, double(i)));
Expand Down
Loading

0 comments on commit 1d9821c

Please sign in to comment.