In Interval, now 'qinterval' is accurate to documentation;
some code is duplicated from AIntervalCategory, which was
moved to there before, but since Interval is for FPS which
is inaccurate, so the rounding is inevitable and the code
has to be duplicated.
Also slightly improves accuracy for 'abs' and '* for ZZ'.
- Qian
--
You received this message because you are subscribed to the Google Groups
"FriCAS - computer algebra system" group.
To unsubscribe from this group and stop receiving emails from it, send an email
to [email protected].
To view this discussion visit
https://groups.google.com/d/msgid/fricas-devel/f2a1ddd7-9271-40f5-9a7f-1c776bf9a46e%40gmail.com.
diff --git a/src/algebra/ainterval.spad b/src/algebra/ainterval.spad
index bdc67939..86d4c9c2 100644
--- a/src/algebra/ainterval.spad
+++ b/src/algebra/ainterval.spad
@@ -111,7 +111,9 @@ AIntervalCategory(R : OrderedAbelianSemiGroup) : Category ==
interval(i : R, s : R) : % ==
s < i => qinterval(s, i)
qinterval(i, s)
- ((z : PP) * (x : %)) : % == qinterval(z * inf x, z * sup x)
+ ((z : PP) * (x : %)) : % ==
+ one? z => x
+ qinterval(z * inf x, z * sup x)
coerce(x : %) : List R == [inf x, sup x]
coerce(x : %) : OutputForm == x :: List(R) :: OutputForm
((x : %) = (y : %)) : Boolean == inf x = inf y and sup x = sup y
diff --git a/src/algebra/interval.spad b/src/algebra/interval.spad
index 8bd84ec1..92e6825e 100644
--- a/src/algebra/interval.spad
+++ b/src/algebra/interval.spad
@@ -33,9 +33,6 @@ IntervalCategory(R : Join(FloatingPointSystem, TranscendentalFunctionCategory)):
Interval(R : Join(FloatingPointSystem, TranscendentalFunctionCategory)) :
IntervalCategory(R) == add
- -- Because of SPAD semantics, Interval deliberately does not inherit
- -- from AInterval, since the implementation of qinterval is different.
-
import from Integer
Rep := Record(Inf : R, Sup : R)
@@ -90,13 +87,58 @@ Interval(R : Join(FloatingPointSystem, TranscendentalFunctionCategory)) :
[roundDown f, roundUp f]
[f, f]
- qinterval(i : R, s : R) : % ==
+ rinterval(i : R, s : R) : % ==
[roundDown i, roundUp s]
- exactInterval(i : R, s : R) : % == [i, s]
+ qinterval(i : R, s : R) : % == [i, s]
exactSupInterval(i : R, s : R) : % == [roundDown i, s]
exactInfInterval(i : R, s : R) : % == [i, roundUp s]
+ -- definitions copied from AIntervalCategory,
+ -- with 'qinterval' replaced by 'rinterval'.
+ interval(i : R, s : R) : % ==
+ s < i => rinterval(s, i)
+ rinterval(i, s)
+
+ ((z : PositiveInteger) * (x : %)) : % ==
+ one? z => x
+ rinterval(z * inf x, z * sup x)
+
+ ((z : Integer) * (x : %)) : % ==
+ zero? z => 0
+ one? z or zero? x => x
+ z > 0 => rinterval(z * inf x, z * sup x)
+ rinterval(z * sup x, z * inf x)
+
+ ((x : %) + (y : %)) : % ==
+ zero? x => y
+ zero? y => x
+ x = y => rinterval(2 * inf x, 2 * sup x)
+ rinterval(inf x + inf y, sup x + sup y)
+
+ ((x : %) + (r : R)) : % ==
+ zero? r => x
+ rinterval(inf x + r, sup x + r)
+
+ ((x : %) - (y : %)) : % ==
+ zero? x => - y
+ zero? y => x
+ rinterval(inf x - sup y, sup x - inf y)
+
+ ((x : %) * (y : %)) : % ==
+ one? x => y
+ one? y => x
+ zero? x or zero? y => 0
+ ix := inf x; sx := sup x; iy := inf y; sy := sup y
+ products : List R := sort [ix * iy, ix * sy, sx * iy, sx * sy]
+ rinterval(first products, last products)
+ -- copy from AIntervalCategory ends here.
+
+ ((x : %) ^ (n : PositiveInteger)) : % ==
+ contains?(x, 0$R) and even? n =>
+ exactInfInterval(0, max(- inf(x), sup(x))^n)
+ interval(inf(x)^n, sup(x)^n) -- will sort so that lower < upper
+
inf(u : %) : R == u.Inf
sup(u : %) : R == u.Sup
@@ -108,7 +150,7 @@ Interval(R : Join(FloatingPointSystem, TranscendentalFunctionCategory)) :
else
sup(a) < sup(b)
- _- (a : %) : % == exactInterval(-sup(a), -inf(a))
+ _- (a : %) : % == qinterval(-sup(a), -inf(a))
_= (a : %, b : %) : Boolean == (inf(a)=inf(b)) and (sup(a)=sup(b))
_~_= (a : %, b : %) : Boolean == (inf(a)~=inf(b)) or (sup(a)~=sup(b))
@@ -122,7 +164,7 @@ Interval(R : Join(FloatingPointSystem, TranscendentalFunctionCategory)) :
recip(u:%):Union(%,"failed") ==
contains?(u,0) => "failed"
vals : List R := sort [1/inf(u), 1/sup(u)]$List(R)
- qinterval(first vals, last vals)
+ rinterval(first vals, last vals)
_exquo(u:%,v:%):Union(%,"failed") ==
@@ -132,14 +174,14 @@ Interval(R : Join(FloatingPointSystem, TranscendentalFunctionCategory)) :
u = v => 1
u=-v => -1
vals : List R := sort [inf(u)/inf(v), inf(u)/sup(v), sup(u)/inf(v), sup(u)/sup(v)]$List(R)
- qinterval(first vals, last vals)
+ rinterval(first vals, last vals)
gcd(u : %, v : %) : % == 1
coerce(u : Integer) : % ==
ur := normaliseFloat(u::R)
- exactInterval(ur, ur)
+ qinterval(ur, ur)
interval(u : Fraction Integer) : % ==
@@ -154,9 +196,9 @@ Interval(R : Join(FloatingPointSystem, TranscendentalFunctionCategory)) :
bin : Union(Integer,"failed") := retractIfCan(log2(denom(u)::Float))
bin case Integer and length(numer u)$Integer < bits() =>
flt := normaliseFloat flt
- exactInterval(flt, flt)
+ qinterval(flt, flt)
- qinterval(flt, flt)
+ rinterval(flt, flt)
retractIfCan(u:%):Union(Integer,"failed") ==
@@ -177,22 +219,22 @@ Interval(R : Join(FloatingPointSystem, TranscendentalFunctionCategory)) :
-- Explicit export from TranscendentalFunctionCategory
- pi() : % == qinterval(pi(), pi())
+ pi() : % == rinterval(pi(), pi())
-- From ElementaryFunctionCategory
log(u : %) : % ==
- positive?(u) => qinterval(log inf u, log sup u)
+ positive?(u) => rinterval(log inf u, log sup u)
error "negative logs in interval"
- exp(u : %) : % == qinterval(exp inf u, exp sup u)
+ exp(u : %) : % == rinterval(exp inf u, exp sup u)
_^ (u : %, v : %) : % ==
zero?(v) => if zero?(u) then error "0^0 is undefined" else 1
one?(u) => 1
expts : List R := sort [inf(u)^inf(v), sup(u)^sup(v),
inf(u)^sup(v), sup(u)^inf(v)]
- qinterval(first expts, last expts)
+ rinterval(first expts, last expts)
-- From TrigonometricFunctionCategory
@@ -216,7 +258,7 @@ Interval(R : Join(FloatingPointSystem, TranscendentalFunctionCategory)) :
hasMinusOne? : Boolean := hasTwoPiMultiple(3*ipi/(2::R), ipi, u)
if hasOne? and hasMinusOne? then
- exactInterval(-1, 1)
+ qinterval(-1, 1)
else
vals : List R := sort [sin inf u, sin sup u]
if hasOne? then
@@ -224,7 +266,7 @@ Interval(R : Join(FloatingPointSystem, TranscendentalFunctionCategory)) :
else if hasMinusOne? then
exactInfInterval(-1, last vals)
else
- qinterval(first vals, last vals)
+ rinterval(first vals, last vals)
@@ -234,7 +276,7 @@ Interval(R : Join(FloatingPointSystem, TranscendentalFunctionCategory)) :
hasMinusOne? : Boolean := hasTwoPiMultiple(ipi, ipi, u)
if hasOne? and hasMinusOne? then
- exactInterval(-1, 1)
+ qinterval(-1, 1)
else
vals : List R := sort [cos inf u, cos sup u]
if hasOne? then
@@ -242,7 +284,7 @@ Interval(R : Join(FloatingPointSystem, TranscendentalFunctionCategory)) :
else if hasMinusOne? then
exactInfInterval(-1, last vals)
else
- qinterval(first vals, last vals)
+ rinterval(first vals, last vals)
@@ -258,7 +300,7 @@ Interval(R : Join(FloatingPointSystem, TranscendentalFunctionCategory)) :
hi : R := tan sup u
lo > hi => error "Interval contains a singularity"
- qinterval(lo, hi)
+ rinterval(lo, hi)
@@ -276,7 +318,7 @@ Interval(R : Join(FloatingPointSystem, TranscendentalFunctionCategory)) :
else if hasTwoPiMultiple(3*ipi/(2::R), ipi, u) then
exactSupInterval(first vals, -1)
else
- qinterval(first vals, last vals)
+ rinterval(first vals, last vals)
@@ -295,7 +337,7 @@ Interval(R : Join(FloatingPointSystem, TranscendentalFunctionCategory)) :
else if hasTwoPiMultiple(ipi, ipi, u) then
exactSupInterval(first vals, -1)
else
- qinterval(first vals, last vals)
+ rinterval(first vals, last vals)
@@ -312,7 +354,7 @@ Interval(R : Join(FloatingPointSystem, TranscendentalFunctionCategory)) :
lo : R := cot sup u
lo > hi => error "Interval contains a singularity"
- qinterval(lo, hi)
+ rinterval(lo, hi)
@@ -322,26 +364,26 @@ Interval(R : Join(FloatingPointSystem, TranscendentalFunctionCategory)) :
lo : R := inf(u)
hi : R := sup(u)
if (lo < -1) or (hi > 1) then error "asin only defined on the region -1..1"
- qinterval(asin lo, asin hi)
+ rinterval(asin lo, asin hi)
acos(u : %) : % ==
lo : R := inf(u)
hi : R := sup(u)
if (lo < -1) or (hi > 1) then error "acos only defined on the region -1..1"
- qinterval(acos hi, acos lo)
+ rinterval(acos hi, acos lo)
- atan(u : %) : % == qinterval(atan inf u, atan sup u)
+ atan(u : %) : % == rinterval(atan inf u, atan sup u)
- acot(u : %) : % == qinterval(acot sup u, acot inf u)
+ acot(u : %) : % == rinterval(acot sup u, acot inf u)
acsc(u : %) : % ==
lo : R := inf(u)
hi : R := sup(u)
if ((lo <= -1) and (hi >= -1)) or ((lo <= 1) and (hi >= 1)) then
error "acsc not defined on the region -1..1"
- qinterval(acsc hi, acsc lo)
+ rinterval(acsc hi, acsc lo)
asec(u : %) : % ==
@@ -349,44 +391,44 @@ Interval(R : Join(FloatingPointSystem, TranscendentalFunctionCategory)) :
hi : R := sup(u)
if ((lo < -1) and (hi > -1)) or ((lo < 1) and (hi > 1)) then
error "asec not defined on the region -1..1"
- qinterval(asec lo, asec hi)
+ rinterval(asec lo, asec hi)
-- From HyperbolicFunctionCategory
- tanh(u : %) : % == qinterval(tanh inf u, tanh sup u)
+ tanh(u : %) : % == rinterval(tanh inf u, tanh sup u)
- sinh(u : %) : % == qinterval(sinh inf u, sinh sup u)
+ sinh(u : %) : % == rinterval(sinh inf u, sinh sup u)
sech(u : %) : % ==
- negative? u => qinterval(sech inf u, sech sup u)
- positive? u => qinterval(sech sup u, sech inf u)
+ negative? u => rinterval(sech inf u, sech sup u)
+ positive? u => rinterval(sech sup u, sech inf u)
vals : List R := sort [sech inf u, sech sup u]
exactSupInterval(first vals, 1)
cosh(u : %) : % ==
- negative? u => qinterval(cosh sup u, cosh inf u)
- positive? u => qinterval(cosh inf u, cosh sup u)
+ negative? u => rinterval(cosh sup u, cosh inf u)
+ positive? u => rinterval(cosh inf u, cosh sup u)
vals : List R := sort [cosh inf u, cosh sup u]
exactInfInterval(1, last vals)
csch(u : %) : % ==
contains?(u,0) => error "csch: singularity at zero"
- qinterval(csch sup u, csch inf u)
+ rinterval(csch sup u, csch inf u)
coth(u : %) : % ==
contains?(u,0) => error "coth: singularity at zero"
- qinterval(coth sup u, coth inf u)
+ rinterval(coth sup u, coth inf u)
-- From ArcHyperbolicFunctionCategory
acosh(u : %) : % ==
inf(u)<1 => error "invalid argument: acosh only defined on the region 1.."
- qinterval(acosh inf u, acosh sup u)
+ rinterval(acosh inf u, acosh sup u)
acoth(u : %) : % ==
@@ -394,12 +436,12 @@ Interval(R : Join(FloatingPointSystem, TranscendentalFunctionCategory)) :
hi : R := sup(u)
if ((lo <= -1) and (hi >= -1)) or ((lo <= 1) and (hi >= 1)) then
error "acoth not defined on the region -1..1"
- qinterval(acoth hi, acoth lo)
+ rinterval(acoth hi, acoth lo)
acsch(u : %) : % ==
contains?(u,0) => error "acsch: singularity at zero"
- qinterval(acsch sup u, acsch inf u)
+ rinterval(acsch sup u, acsch inf u)
asech(u : %) : % ==
@@ -407,17 +449,17 @@ Interval(R : Join(FloatingPointSystem, TranscendentalFunctionCategory)) :
hi : R := sup(u)
if (lo <= 0) or (hi > 1) then
error "asech only defined on the region 0 < x <= 1"
- qinterval(asech hi, asech lo)
+ rinterval(asech hi, asech lo)
- asinh(u : %) : % == qinterval(asinh inf u, asinh sup u)
+ asinh(u : %) : % == rinterval(asinh inf u, asinh sup u)
atanh(u : %) : % ==
lo : R := inf(u)
hi : R := sup(u)
if (lo <= -1) or (hi >= 1) then
error "atanh only defined on the region -1 < x < 1"
- qinterval(atanh lo, atanh hi)
+ rinterval(atanh lo, atanh hi)
-- From RadicalCategory
