Skip to content

Commit

Permalink
Update interval library.
Browse files Browse the repository at this point in the history
  • Loading branch information
@sletz
sletz committed Apr 9, 2024
1 parent b41cfee commit 0e040ff
Show file tree
Hide file tree
Showing 51 changed files with 364 additions and 203 deletions.
96 changes: 64 additions & 32 deletions compiler/interval/check.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -34,9 +34,13 @@ void check(const std::string& expected, const itv::interval& exp)
std::stringstream ss;
ss << exp;
if (ss.str().compare(expected) == 0) {
std::cout << "OK: " << expected << std::endl;
std::cout << "\033[32m"
<< "OK: " << expected
<< "\033[0m" << std::endl;
} else {
std::cout << "ERR: We got " << ss.str() << " instead of " << expected << std::endl;
std::cout << "\033[31m"
<< "ERR: We got " << ss.str() << " instead of " << expected
<< "\033[0m" << std::endl;
}
}

Expand All @@ -50,9 +54,17 @@ void check(const std::string& expected, const itv::interval& exp)
void check(const std::string& testname, const itv::interval& exp, const itv::interval& res)
{
if (exp == res) {
std::cout << "OK: " << testname << " " << exp << " = " << res << std::endl;
std::cout << "\033[32m"
<< "OK: " << testname << " " << exp << " = " << res
<< "\033[0m" << std::endl;
if (exp.lsb() != res.lsb())
std::cout << "\033[33m"
<< "\t But precisions differ by " << exp.lsb() - res.lsb()
<< "\033[0m" << std::endl;
} else {
std::cout << "ERR:" << testname << " FAILED. We got " << exp << " instead of " << res << std::endl;
std::cout << "\033[31m"
<< "ERR:" << testname << " FAILED. We got " << exp << " instead of " << res
<< "\033[0m" << std::endl;
}
}

Expand All @@ -66,9 +78,13 @@ void check(const std::string& testname, const itv::interval& exp, const itv::int
void check(const std::string& testname, bool exp, bool res)
{
if (exp == res) {
std::cout << "OK: " << testname << std::endl;
std::cout << "\033[32m"
<< "OK: " << testname
<< "\033[0m" << std::endl;
} else {
std::cout << "ERR:" << testname << " FAILED. We got " << exp << " instead of " << res << std::endl;
std::cout << "\033[31m"
<< "ERR:" << testname << " FAILED. We got " << exp << " instead of " << res
<< "\033[0m" << std::endl;
}
}

Expand Down Expand Up @@ -232,19 +248,27 @@ void analyzeUnaryMethod(int E, int M, const char* title, const itv::interval& D,
}

for (int m = 0; m < M; m++) { // M measurements
double presample = rx(generator);
sample = truncate(presample, D.lsb());
double presample = rx(generator); // non-truncated sample
sample = truncate(presample, D.lsb()); // truncated sample
double pre_y = f(presample);
y = f(sample);
// y = truncate(y, -30); // workaround to avoid artefacts in trigonometric functions

measurements.insert(y);

if (!std::isnan(y)) {
if (y < y0) {
// interval bounds
if (!std::isnan(pre_y)) {
/* if (y < y0) {
y0 = y;
}
if (y > y1) {
y1 = y;
} */
if (pre_y < y0) {
y0 = pre_y;
}
if (pre_y > y1) {
y1 = pre_y;
}
}
}
Expand All @@ -264,19 +288,20 @@ void analyzeUnaryMethod(int E, int M, const char* title, const itv::interval& D,
}

itv::interval Y(y0, y1, lsb);
if (y0 > y1) Y = itv::empty(); // if we didn't manage to draw any samples
itv::interval Z = (A.*mp)(X);

if (Z >= Y and Z.lsb() <= Y.lsb()) {
double precision = (Z.size() == 0) ? 1 : Y.size() / Z.size();

std::cout << "\033[32m"
<< "OK " << e << ": " << title << "(" << X << ") = \t" << Z << "\t >= \t" << Y
<< "\t (precision " << precision << ", LSB diff = " << Y.lsb() - Z.lsb() << ")"
<< "OK " << e << ": " << title << "(" << X << ") = \n" << Z << "(c)\t >= \t" << Y
<< "(m)\t (precision " << precision << ", LSB diff = " << Y.lsb() - Z.lsb() << ")"
<< "\033[0m" << std::endl;
} else {
std::cout << "\033[31m"
<< "ERROR " << e << ": " << title << "(" << X << ") = \t" << Z << "\t INSTEAD OF \t" << Y
<< ", \t LSB diff = " << Y.lsb() - Z.lsb() << "\033[0m" << std::endl;
<< "ERROR " << e << ": " << title << "(" << X << ") = \n" << Z << "(c)\t INSTEAD OF \t" << Y
<< "(m), \t LSB diff = " << Y.lsb() - Z.lsb() << "\033[0m" << std::endl;
}
std::cout << std::endl;
}
Expand Down Expand Up @@ -310,7 +335,6 @@ void analyzeBinaryMethod(int E, int M, const char* title, const itv::interval& D

// store output values in order to measure the output precision
std::set<double> measurements;

if (Dx.lsb() < 0 or Dy.lsb() < 0) // if we're not doing an integer operation
{
// X: random input interval X < Dx
Expand Down Expand Up @@ -378,33 +402,34 @@ void analyzeBinaryMethod(int E, int M, const char* title, const itv::interval& D
}

itv::interval Zm(zlo, zhi, lsb); // the measured Z
if (zlo > zhi) Zm = itv::empty(); // if we didn't manage to draw any samples
itv::interval Zc = (A.*bm)(X, Y); // the computed Z
double precision = (Zm.size() == Zc.size()) ? 1 : Zm.size() / Zc.size();

if (Zc >= Zm and Zc.lsb() <= Zm.lsb()) {
std::string color = "\033[32m";
if (precision < 0.8 or Zm.lsb() - Zc.lsb() >= 10) color = "\033[36m"; // cyan instead of green if approximation is technically correct but of poor quality
std::cout << color
<< "OK " << e << ": " << title << "(" << X << ",\t" << Y << ")\t =c=> " << Zc << " >= " << Zm
<< "OK " << e << ": " << title << "(" << X << ",\t" << Y << ")\n =c=> " << Zc << "(c) >= " << Zm << "(m)"
<< "\t (precision " << precision << "), \t LSB diff = " << Zm.lsb() - Zc.lsb()
<< "\033[0m" << std::endl;
} else {
std::cout << "\033[31m"
<< "ERROR " << e << ": " << title << "(" << X << ",\t" << Y << ")\t =c=> " << Zc << " != " << Zm
<< "ERROR " << e << ": " << title << "(" << X << ",\t" << Y << ")\n =c=> " << Zc << "(c) != " << Zm << "(m)"
<< "\t LSB diff = " << Zm.lsb() - Zc.lsb()
<< "\033[0m" << std::endl;
}

} else { // integer operation
std::cout << "Testing integer version of " << title << std::endl;
// std::cout << "Testing integer version of " << title << std::endl;
// X: random input interval X < Dx
int x0 = truncate(rdx(generator), Dx.lsb());
int x1 = truncate(rdx(generator), Dx.lsb());
double x0 = truncate(rdx(generator), Dx.lsb());
double x1 = truncate(rdx(generator), Dx.lsb());
itv::interval X(std::min(x0, x1), std::max(x0, x1), Dx.lsb());

// Y: random input interval Y < Dy
int y0 = truncate(rdy(generator), Dy.lsb());
int y1 = truncate(rdy(generator), Dy.lsb());
double y0 = truncate(rdy(generator), Dy.lsb());
double y1 = truncate(rdy(generator), Dy.lsb());
itv::interval Y(std::min(y0, y1), std::max(y0, y1), Dy.lsb());

// boundaries of the resulting interval Z
Expand All @@ -431,15 +456,22 @@ void analyzeBinaryMethod(int E, int M, const char* title, const itv::interval& D

// measure the interval Z using the numerical function f
for (int m = 0; m < M; m++) { // M measurements
z = f(truncate(ivx(generator), Dx.lsb()), truncate(ivy(generator), Dy.lsb()));
int pre_x = ivx(generator);
int x = truncate(pre_x, Dx.lsb());
int pre_y = ivy(generator);
int y = truncate(pre_y, Dy.lsb());
z = f(x, y);
int pre_z = f(pre_x, pre_y);

measurements.insert(z);

if (z < zlo) {
zlo = z;
}
if (z > zhi) {
zhi = z;
if (!std::isnan(pre_z)){
if (z < zlo) {
zlo = pre_z;
}
if (z > zhi) {
zhi = pre_z;
}
}
}

Expand All @@ -465,12 +497,12 @@ void analyzeBinaryMethod(int E, int M, const char* title, const itv::interval& D
std::string color = "\033[32m";
if (precision < 0.8 or Zm.lsb() - Zc.lsb() >= 10) color = "\033[36m"; // cyan instead of green if approximation is technically correct but of poor quality
std::cout << color
<< "OK " << e << ": " << title << "(" << X << ",\t" << Y << ")\t =c=> " << Zc << " >= " << Zm
<< "OK " << e << ": " << title << "(" << X << ",\t" << Y << ")\n =c=> " << Zc << "(c) >= " << Zm << "(m)"
<< "\t (precision " << precision << "), \t LSB diff = " << Zm.lsb() - Zc.lsb()
<< "\033[0m" << std::endl;
} else {
std::cout << "\033[31m"
<< "ERROR " << e << ": " << title << "(" << X << ",\t" << Y << ")\t =c=> " << Zc << " != " << Zm
<< "ERROR " << e << ": " << title << "(" << X << ",\t" << Y << ")\n =c=> " << Zc << "(c) != " << Zm << "(m)"
<< "\t LSB diff = " << Zm.lsb() - Zc.lsb()
<< "\033[0m" << std::endl;
}
Expand Down Expand Up @@ -562,10 +594,10 @@ void propagateBackwardsComposition(std::vector<const char*> titles, std::vector<
return;
}

int n = titles.size();
unsigned long n = titles.size();

std::cout << "Shaving input " << X << " of ";
for(auto t: titles)
for(const auto *t: titles)
std::cout << t << " ○ ";
std::cout << "\b\b\b";
std::cout << " to achieve an output lsb of " << l << std::endl << std::endl;
Expand Down
22 changes: 22 additions & 0 deletions compiler/interval/intervalAbs.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -28,21 +28,43 @@ namespace itv {

interval interval_algebra::Abs(const interval& x)
{
if (x.isEmpty())
return empty();

// precision stays the same
if (x.lo() >= 0) {
return x;
}

// integer overflowing
if (x.lsb() >= 0 and x.lo() <= (double)INT_MIN){
double lo = (x.hi() >= 0) ? 0 : std::min(std::abs(x.hi()), (double)INT_MAX);
return {lo,
(double)INT_MAX,
x.lsb()};
}
if (x.hi() <= 0) {
return {-x.hi(), -x.lo(), x.lsb()};
}

return {0, std::max(std::abs(x.lo()), std::abs(x.hi())), x.lsb()};
}

void interval_algebra::testAbs()
{
std::cout << "Testing abs on finite intervals" << std::endl;
analyzeUnaryMethod(
10, 10000, "abs", interval(-10, 10, 0), [](double x) { return std::abs(x); }, &interval_algebra::Abs);
analyzeUnaryMethod(
10, 10000, "abs", interval(-10, 10, -15), [](double x) { return std::abs(x); }, &interval_algebra::Abs);

std::cout << "Testing abs on intervals with limit bounds" << std::endl;
check("abs", Abs(interval(INT_MIN, INT_MIN, 0)), interval(INT_MAX, INT_MAX, 0));
check("abs", Abs(interval(INT_MIN, -10, 0)), interval(10, INT_MAX, 0));
check("abs", Abs(interval(INT_MIN, 10, 0)), interval(0, INT_MAX, 0));
check("abs", Abs(interval(-10, INT_MAX, 0)), interval(0, INT_MAX, 0));

check("abs", Abs(interval(INT_MIN, INT_MIN, 5)), interval(INT_MAX, INT_MAX, 5));
check("abs", Abs(interval(INT_MIN, INT_MIN, -5)), interval(-(float)INT_MIN, -(float)INT_MIN, -5));
}
} // namespace itv
4 changes: 2 additions & 2 deletions compiler/interval/intervalAcos.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -32,7 +32,7 @@ interval interval_algebra::Acos(const interval& x)
{
interval i = intersection(AcosDomain, x); // TODO: warn about interval violations
if (i.isEmpty()) {
return i;
return empty();
}

double v = 0; // value at which the min slope is attained, zero if it is present
Expand All @@ -47,7 +47,7 @@ interval interval_algebra::Acos(const interval& x)
int precision = exactPrecisionUnary(acos, v, sign * pow(2, i.lsb()));

if (precision == INT_MIN or taylor_lsb)
precision = floor(x.lsb() - (double)log2(1 - v*v)/2);
precision = floor(i.lsb() - (double)log2(1 - v*v)/2);

return {acos(i.hi()), acos(i.lo()), precision};
}
Expand Down
2 changes: 1 addition & 1 deletion compiler/interval/intervalAcosh.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -32,7 +32,7 @@ interval interval_algebra::Acosh(const interval& x)
{
interval i = intersection(domain, x);
if (i.isEmpty()) {
return i;
return empty();
}

// the min slope is attained at the highest bound of the interval
Expand Down
49 changes: 37 additions & 12 deletions compiler/interval/intervalAdd.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -38,7 +38,7 @@ static double addint(double x, double y)
interval interval_algebra::Add(const interval& x, const interval& y)
{
if (x.isEmpty() || y.isEmpty()) {
return {};
return empty();
}

if (x.lsb() >= 0 and y.lsb() >= 0) { // if both intervals are integers
Expand All @@ -49,12 +49,28 @@ interval interval_algebra::Add(const interval& x, const interval& y)
const int yhi = (int)y.hi();

// detect wrapping
if (std::abs((double)xhi + (double)yhi) >= (double) INT_MAX
/* if (std::abs((double)xhi + (double)yhi) >= (double) INT_MAX
or std::abs((double)xhi + (double)yhi) <= (double) INT_MIN
or std::abs((double)xlo + (double)ylo) >= (double) INT_MAX
or std::abs((double)xlo + (double)ylo) <= (double) INT_MIN)
return {(double) INT_MIN, (double) INT_MAX, std::min(x.lsb(), y.lsb())};
return {(double) INT_MIN, (double) INT_MAX, std::min(x.lsb(), y.lsb())};*/

double lo = x.lo() + y.lo();
double hi = x.hi() + y.hi();

// if there is a discontinuity by the lower end of integers
if (lo <= (double)INT_MIN -1 and hi >= (double)INT_MIN)
{
return {(double)INT_MIN, (double)INT_MAX, std::min(x.lsb(), y.lsb())};
}

// if there is a discontinuity by the higher end of integers
if (lo <= (double)INT_MAX and hi >= (double)INT_MAX+1)
{
return {(double)INT_MIN, (double)INT_MAX, std::min(x.lsb(), y.lsb())};
}

// if there is potential wrapping but no discontinuity
return {(double)(xlo + ylo), (double)(xhi + yhi), std::min(x.lsb(), y.lsb())};
}

Expand All @@ -65,15 +81,24 @@ interval interval_algebra::Add(const interval& x, const interval& y)

void interval_algebra::testAdd()
{
// check("test algebra Add", Add(interval(0, 100), interval(10, 500)), interval(10, 600));
// analyzeBinaryMethod(10, 2000, "add", interval(0, 10, 0), interval(0, 10, 0), add, &interval_algebra::Add);
std::cout << "Testing add on regular intervals" << std::endl;
check("test algebra Add", Add(interval(0, 100), interval(10, 500)), interval(10, 600));

std::cout << "Testing add on positive intervals" << std::endl;
analyzeBinaryMethod(5, 2000, "add", interval(0, 10, 0), interval(0, 10, 0), add, &interval_algebra::Add);
analyzeBinaryMethod(5, 2000, "add", interval(0, 10, -10), interval(0, 10, 0), add, &interval_algebra::Add);
analyzeBinaryMethod(5, 2000, "add", interval(0, 10, 0), interval(0, 10, -10), add, &interval_algebra::Add);
analyzeBinaryMethod(5, 2000, "add", interval(0, 10, -10), interval(0, 10, -10), add, &interval_algebra::Add);

std::cout << "Testing add on negative intervals" << std::endl;
analyzeBinaryMethod(5, 2000, "add", interval(-10, 10, -5), interval(-10, 0, -5), add, &interval_algebra::Add);
analyzeBinaryMethod(5, 2000, "add", interval(-10, 0, -5), interval(-10, 0, -5), add, &interval_algebra::Add);

std::cout << "Testing add with wrapping" << std::endl;

analyzeBinaryMethod(10, 2000, "add", interval(0, pow(2, 31), 0), interval(0, pow(2, 31), 0), addint, &interval_algebra::Add);
/* analyzeBinaryMethod(10, 2000, "add", interval(0, 10, -5), interval(0, 10, 0), add, &interval_algebra::Add);
analyzeBinaryMethod(10, 2000, "add", interval(0, 10, -10), interval(0, 10, 0), add, &interval_algebra::Add);
analyzeBinaryMethod(10, 2000, "add", interval(0, 10, -15), interval(0, 10, 0), add, &interval_algebra::Add);
analyzeBinaryMethod(10, 2000, "add", interval(0, 10, 0), interval(0, 10, -10), add, &interval_algebra::Add);
analyzeBinaryMethod(10, 2000, "add", interval(0, 10, -5), interval(0, 10, -10), add, &interval_algebra::Add);
analyzeBinaryMethod(10, 2000, "add", interval(0, 10, -10), interval(0, 10, -10), add, &interval_algebra::Add);
analyzeBinaryMethod(10, 2000, "add", interval(0, 10, -15), interval(0, 10, -10), add, &interval_algebra::Add);*/
analyzeBinaryMethod(10, 2000, "add", interval((double)INT_MAX/2, (double)INT_MAX, 0), interval((double)INT_MAX/2, (double)INT_MAX, 0), addint, &interval_algebra::Add);
analyzeBinaryMethod(10, 2000, "add", interval((double)INT_MIN, (double)INT_MIN/2, 0), interval((double)INT_MAX/2, (double)INT_MAX, 0), addint, &interval_algebra::Add);
analyzeBinaryMethod(10, 2000, "add", interval((double)INT_MIN, (double)INT_MIN/2, 0), interval((double)INT_MIN, (double)INT_MIN/2, 0), addint, &interval_algebra::Add);
}
} // namespace itv
12 changes: 8 additions & 4 deletions compiler/interval/intervalAnd.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -130,7 +130,7 @@ interval interval_algebra::And(const interval& x, const interval& y)
interval interval_algebra::And(const interval& x, const interval& y)
{
if (x.isEmpty() || y.isEmpty()) {
return {};
return empty();
}
int x0 = saturatedIntCast(x.lo());
int x1 = saturatedIntCast(x.hi());
Expand All @@ -139,7 +139,9 @@ interval interval_algebra::And(const interval& x, const interval& y)

SInterval z = bitwiseSignedAnd({x0, x1}, {y0, y1});

int precision = std::max(x.lsb(), y.lsb()); // output precision cannot be finer than that of the input intervals
int precision = std::min(x.lsb(), y.lsb());

/* int precision = std::max(x.lsb(), y.lsb()); // output precision cannot be finer than that of the input intervals
// however, if one of the intervals is reduced to one element, the mask can make it so
int precisionx = 0;
Expand All @@ -162,9 +164,11 @@ interval interval_algebra::And(const interval& x, const interval& y)
v = v / 2;
precisiony++;
}
}
}*/

return {double(z.lo), double(z.hi), std::max(precision, std::max(precisionx, precisiony))};
return {double(z.lo), double(z.hi),
// std::max(precision, std::max(precisionx, precisiony))
precision};
}

void interval_algebra::testAnd()
Expand Down
Loading

0 comments on commit 0e040ff

Please sign in to comment.