Deals with empty intervals

interval/intervalAbs.cpp

@@ -28,6 +28,9 @@ namespace itv {
 
 interval interval_algebra::Abs(const interval& x)
 {
+    if (x.isEmpty())
+        return empty();
+
     // precision stays the same
     if (x.lo() >= 0) {
         return x;

interval/intervalAcos.cpp

@@ -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

interval/intervalAcosh.cpp

@@ -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

interval/intervalAdd.cpp

@@ -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

interval/intervalAnd.cpp

@@ -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());

interval/intervalAsinh.cpp

@@ -29,6 +29,9 @@ static const interval domain(-HUGE_VAL, HUGE_VAL);
 
 interval interval_algebra::Asinh(const interval& x)
 {
+    if (x.isEmpty())
+        return empty();
+
     double v    = maxValAbs(x);  // value at which the min slope is attained, here the bound of highest absolute value
     int    sign = signMaxValAbs(x);  // whether we compute the difference between f(v) and f(v+ε) or f(v-ε)
 

interval/intervalAtan.cpp

@@ -32,7 +32,7 @@ namespace itv {
 interval interval_algebra::Atan(const interval& x)
 {
     if (x.isEmpty()) {
-        return x;
+        return empty();
     }
 
     double v    = maxValAbs(x);  // value at which the min slope is attained, here the bound of highest absolute value

interval/intervalAtan2.cpp

@@ -31,6 +31,9 @@ namespace itv {
 // (where (x,y) are the cartesian coordinates of the point we wish to retrieve the angle of)
 interval interval_algebra::Atan2(const interval& y, const interval& x)
 {
+    if (x.isEmpty() || y.isEmpty())
+        return empty();
+
     double lo = -M_PI;
     double hi = M_PI;
 

interval/intervalAtanh.cpp

@@ -33,7 +33,7 @@ interval interval_algebra::Atanh(const interval& x)
 {
     interval i = intersection(domain, x);
     if (i.isEmpty()) {
-        return i;
+        return empty();
     }
 
     double v = minValAbs(x);

interval/intervalCeil.cpp

@@ -29,7 +29,7 @@ namespace itv {
 interval interval_algebra::Ceil(const interval& x)
 {
     if (x.isEmpty()) {
-        return x;
+        return empty();
     }
     return {ceil(x.lo()), ceil(x.hi()), -1 };  // even though the output of floor are mathematical integers,
                                             // they are implemented as floats and thus should not be given precision 0,

interval/intervalCos.cpp

@@ -38,7 +38,7 @@ interval interval_algebra::Cos(const interval& x)
     if (precision == INT_MIN or taylor_lsb) precision = 2*x.lsb() - 1; // if x.lsb() is so small that the automatic computation doesn't work
 
     if (x.isEmpty()) {
-        return {};
+        return empty();
     }
     if (x.size() >= 2*M_PI) {
         return {-1, 1, precision};

interval/intervalCosh.cpp

@@ -31,7 +31,7 @@ interval interval_algebra::Cosh(const interval& x)
     double v = 0;  // absolute lowest slope is at zero
 
     if (x.isEmpty()) {
-        return x;
+        return empty();
     }
 
     if (x.hasZero()) {

interval/intervalDelay.cpp

@@ -29,7 +29,7 @@ namespace itv {
 interval interval_algebra::Delay(const interval& x, const interval& y)
 {
     if (x.isEmpty() || y.isEmpty()) {
-        return {};
+        return empty();
     }
     if (y.isZero()) {
         return x;

interval/intervalEq.cpp

@@ -34,7 +34,7 @@ static double myEq(double x, double y)
 interval interval_algebra::Eq(const interval& x, const interval& y)
 {
     // boolean value => precision 0
-    if (x.isEmpty() || y.isEmpty()) return interval{};
+    if (x.isEmpty() || y.isEmpty()) return empty();
     if (x.lo() == x.hi() && x.lo() == y.lo() && x.lo() == y.hi()) return interval{1,1,0};
     if (x.hi() < y.lo() || x.lo() > y.hi()) return interval{0, 0, 0};
     return interval{0, 1, 0};

interval/intervalFloor.cpp

@@ -29,7 +29,7 @@ namespace itv {
 interval interval_algebra::Floor(const interval& x)
 {
     if (x.isEmpty()) {
-        return x;
+        return empty();
     }
     return {floor(x.lo()), floor(x.hi()), -1};  // even though the output of floor are mathematical integers,
                                                 // they are implemented as floats and thus should not be given precision 0,

interval/intervalGe.cpp

@@ -34,7 +34,7 @@ static double myGe(double x, double y)
 interval interval_algebra::Ge(const interval& x, const interval& y)
 {
     // boolean value => precision 0
-    if (x.isEmpty() || y.isEmpty()) return interval{};
+    if (x.isEmpty() || y.isEmpty()) return empty();
     if (x.lo() >= y.hi()) return interval{1,1,0};
     if (x.hi() < y.lo()) return interval{0,0,0};
     return interval{0, 1, 0};

interval/intervalGt.cpp

@@ -34,7 +34,7 @@ static double myGt(double x, double y)
 interval interval_algebra::Gt(const interval& x, const interval& y)
 {
     if (x.isEmpty() || y.isEmpty()) {
-        return interval{};
+        return empty();
     }
     if (x.lo() > y.hi()) {
         return interval{1,1,0};

interval/intervalHSlider.cpp

@@ -24,13 +24,16 @@ namespace itv {
 interval interval_algebra::HSlider(const interval& name, const interval& init, const interval& lo, const interval& hi,
                                    const interval& step)
 {
-    // elements of a slider with range [lo; hi] and step step are of the form lo + k·step <= hi with k an integer
+    if (init.isEmpty() || lo.isEmpty() || hi.isEmpty() || step.isEmpty())
+        return empty();
+
+    // elements of a slider with range [lo; hi] and step s are of the form lo + k·s <= hi with k an integer
     // the precision needed to represent such elements is the minimum between
     int lsb = std::min(step.lsb(), lo.lsb());       // the precision of the lower bound and that of the step
     if (step.lo() > 0) {                            // if we don't have negative or zero steps
         lsb = std::min(lsb, (int)log2(step.lo()));  // and that associated to the smallest value the step can take
     }
 
-    return {lo.lo(), hi.hi(), lsb};  // TODO: step, init
+    return {lo.lo(), hi.hi(), lsb}; 
 }
 }  // namespace itv

interval/intervalIntCast.cpp

@@ -31,7 +31,7 @@ namespace itv {
 interval interval_algebra::IntCast(const interval& x)
 {
     if (x.isEmpty()) {
-        return {};
+        return empty();
     }
     return {double(saturatedIntCast(x.lo())), double(saturatedIntCast(x.hi())),
             0};  // integer intervals have 0 bits of precision

interval/intervalInv.cpp

@@ -35,7 +35,7 @@ static double inv(double x)
 interval interval_algebra::Inv(const interval& x)
 {
     if (x.isEmpty()) {
-        return {};
+        return empty();
     }
 
     int    sign = signMaxValAbs(x);

interval/intervalLog.cpp

@@ -29,7 +29,7 @@ namespace itv {
 interval interval_algebra::Log(const interval& x)
 {
     if (x.isEmpty()) {
-        return {};
+        return empty();
     }
 
     // lowest slope is at the highest bound of the interval

interval/intervalLog10.cpp

@@ -29,7 +29,7 @@ namespace itv {
 interval interval_algebra::Log10(const interval& x)
 {
     if (x.isEmpty()) {
-        return {};
+        return empty();
     }
 
     // lowest slope is at the highest bound of the interval

interval/intervalLsh.cpp

@@ -34,6 +34,9 @@ static double lsh(double x, double y)
 
 interval interval_algebra::Lsh(const interval& x, const interval& y)
 {
+    if (x.isEmpty() || y.isEmpty())
+        return empty();
+
     interval j{pow(2, y.lo()), pow(2, y.hi())};
     interval z = Mul(x, j);
 

interval/intervalMax.cpp

@@ -29,7 +29,7 @@ namespace itv {
 interval interval_algebra::Max(const interval& x, const interval& y)
 {
     if (x.isEmpty() || y.isEmpty()) {
-        return {};
+        return empty();
     }
 
     return {

interval/intervalMem.cpp

@@ -29,7 +29,7 @@ namespace itv {
 interval interval_algebra::Mem(const interval& x)
 {
     if (x.isEmpty()) {
-        return {};
+        return empty();
     }
     interval z = reunion(x, interval{0});
     return {z.lo(), z.hi(), x.lsb()};

interval/intervalMin.cpp

@@ -28,7 +28,7 @@ namespace itv {
 
 interval interval_algebra::Min(const interval& x, const interval& y)
 {
-    if (x.isEmpty() || y.isEmpty()) return {};
+    if (x.isEmpty() || y.isEmpty()) return empty();
 
     return {std::min(x.lo(), y.lo()), std::min(x.hi(), y.hi()), std::min(x.lsb(), y.lsb())}; // resulting interval should be as precise as the most precise of the operands
 }

interval/intervalMul.cpp

@@ -43,7 +43,7 @@ static double specialmultint(double a, double b)
 interval interval_algebra::Mul(const interval& x, const interval& y)
 {
     if (x.isEmpty() || y.isEmpty()) {
-        return {};
+        return empty();
     }
 
     if (x.lsb() >= 0 and y.lsb() >= 0) // operation between integers

interval/intervalNe.cpp

@@ -29,7 +29,7 @@ namespace itv {
 interval interval_algebra::Ne(const interval& x, const interval& y)
 {
     if (x.isEmpty() || y.isEmpty()) {
-        return {};
+        return empty();
     }
     if ((x.hi() < y.lo()) || x.lo() > y.hi()) {
         return interval{1, 1, 0};

interval/intervalNeg.cpp

@@ -27,7 +27,7 @@ namespace itv {
 interval interval_algebra::Neg(const interval& x)
 {
     if (x.isEmpty()) {
-        return {};
+        return empty();
     }
 
     return {-x.hi(), -x.lo(), x.lsb()};

interval/intervalNot.cpp

@@ -29,7 +29,7 @@ namespace itv {
 interval interval_algebra::Not(const interval& x)
 {
     if (x.isEmpty()) {
-        return x;
+        return empty();
     }
     int x0 = saturatedIntCast(x.lo());
     int x1 = saturatedIntCast(x.hi());

interval/intervalNumEntry.cpp

@@ -22,6 +22,9 @@ namespace itv {
 interval interval_algebra::NumEntry(const interval& name, const interval& init, const interval& lo, const interval& hi,
                                     const interval& step)
 {
+    if (init.isEmpty() || lo.isEmpty() || hi.isEmpty() || step.isEmpty())
+        return empty();
+
     // elements of a slider with range [lo; hi] and step step are of the form lo + k·step <= hi with k an integer
     // the precision needed to represent such elements is the minimum between
     int lsb = std::min(step.lsb(), lo.lsb());       // the precision of the lower bound and that of the step

interval/intervalOr.cpp

@@ -38,7 +38,7 @@ static double myOr(double x, double y)
 interval interval_algebra::Or(const interval& x, const interval& y)
 {
     if (x.isEmpty() || y.isEmpty()) {
-        return {};
+        return empty();
     }
     int x0 = saturatedIntCast(x.lo());
     int x1 = saturatedIntCast(x.hi());

interval/intervalPow.cpp

@@ -74,13 +74,19 @@ static interval ipow(const interval& x, int y)
  */
 interval interval_algebra::fPow(const interval& x, const interval& y)
 {
+    if (x.isEmpty() || y.isEmpty())
+        return empty();
+
     assert(x.lo() > 0);
     // x all positive
     return Exp(Mul(y, Log(x)));
 }
 
 interval interval_algebra::iPow(const interval& x, const interval& y)
 {
+    if (x.isEmpty() || y.isEmpty())
+        return empty();
+
     int      y0 = std::max(0, saturatedIntCast(y.lo()));
     int      y1 = std::max(0, saturatedIntCast(y.hi()));
     interval z  = ipow(x, y0);

interval/intervalRint.cpp

@@ -28,6 +28,9 @@ namespace itv {
 
 interval interval_algebra::Rint(const interval& x)
 {
+    if (x.isEmpty())
+        return empty();
+
     return {rint(x.lo()), rint(x.hi()), 0}; // round to nearest integral value => integer => precision 0
 }
 

interval/intervalRound.cpp

@@ -28,6 +28,9 @@ namespace itv {
 
 interval interval_algebra::Round(const interval& x)
 {
+    if (x.isEmpty())
+        return empty();
+
     return {round(x.lo()), round(x.hi()), 0}; // round to integral value (regardless of rounding direction) => integer => precision 0
 }
 

interval/intervalRsh.cpp

@@ -33,6 +33,9 @@ static double rsh(double x, double y)
 
 interval interval_algebra::Rsh(const interval& x, const interval& y)
 {
+    if (x.isEmpty() || y.isEmpty())
+        return empty();
+
     interval j{pow(2, -y.hi()), pow(2, -y.lo())};
     interval z = Mul(x, j);
 

interval/intervalSin.cpp

@@ -34,6 +34,9 @@ static double sinPi(double x)
 
 interval interval_algebra::Sin(const interval& x)
 {
+    if (x.isEmpty())
+        return empty();
+
     int precision = exactPrecisionUnary(sin, 0.5, pow(2, x.lsb()));
     if (precision == INT_MIN or taylor_lsb) precision = 2*x.lsb() - 1; // if x.lsb() is so small that the automatic computation doesn't work
 

interval/intervalSinh.cpp

@@ -29,7 +29,7 @@ namespace itv {
 interval interval_algebra::Sinh(const interval& x)
 {
     if (x.isEmpty()) {
-        return x;
+        return empty();
     }
 
     double v    = 0;  // absolute lowest slope is at zero

interval/intervalSqrt.cpp

@@ -33,11 +33,11 @@ interval interval_algebra::Sqrt(const interval& x)
     interval i = intersection(SqrtDomain, x);
 
     if (i.isEmpty()) {
-        return i;
+        return empty();
     }
-    if (i.lo() < 0) {
+    /* if (i.lo() < 0) {
         return {};  // sqrt of negative numbers
-    }
+    }*/
 
     // lowest slope at the highest bound of the interval
     int precision = exactPrecisionUnary(sqrt, i.hi(), -pow(2, i.lsb()));

interval/intervalSub.cpp

@@ -32,7 +32,7 @@ static double sub(double a, double b)
 interval interval_algebra::Sub(const interval& x, const interval& y)
 {
     if (x.isEmpty() || y.isEmpty()) {
-        return {};
+        return empty();
     }
 
     return {x.lo() - y.hi(), x.hi() - y.lo(), std::min(x.lsb(), y.lsb())};