Safe addition: Difference between revisions
Content added Content deleted
(Added Hare) |
Thundergnat (talk | contribs) m (syntax highlighting fixup automation) |
||
Line 28:
=={{header|Ada}}==
An interval type based on Float:
<
Lower : Float;
Upper : Float;
end record;</
Implementation of safe addition:
<
Result : constant Float := A + B;
begin
Line 51:
return (Float'Adjacent (0.0, Float'First), Float'Adjacent (0.0, Float'Last));
end if;
end "+";</
The implementation uses the attribute T'Machine_Rounds to determine if rounding is performed on inexact results. If the machine rounds a symmetric interval around the result is used. When the machine does not round, it truncates. Truncating is rounding towards zero. In this case the implementation is twice better (in terms of the result interval width), because depending on its sign, the outcome of addition can be taken for either of the bounds. Unfortunately most of modern processors are rounding.
Test program:
<
procedure Test_Interval_Addition is
-- Definitions from above
Line 64:
begin
Put (1.14 + 2000.0);
end Test_Interval_Addition;</
Sample output:
<pre> 2.00113989257813E+03, 2.00114013671875E+03</pre>
=={{header|AutoHotkey}}==
<
SetFormat, FloatFast, 0.20
Line 82:
err:=0.1**(SubStr(A_FormatFloat,3) > 15 ? 15 : SubStr(A_FormatFloat,3))
Return "[" a+b-err ","a+b+err "]"
}</
=={{header|C}}==
Line 97:
=== C99 with fesetround() ===
<
#include <stdio.h> /* printf() */
Line 146:
}
return 0;
}</
Output:
<pre>
Line 169:
{{works with|MinGW}}
<
#include <stdio.h> /* printf() */
Line 217:
}
return 0;
}</
=== fpsetround() ===
{{works with|OpenBSD|4.8/amd64}}
<
#include <stdio.h> /* printf() */
Line 270:
}
return 0;
}</
Output from OpenBSD: <pre>1 + 2 =
Line 291:
=={{header|C sharp|C#}}==
{{trans|Java}}
<
namespace SafeAddition {
Line 330:
}
}
}</
{{out}}
<pre>(1.2 + 0.03) is in the range (1.23, 1.23)</pre>
Line 336:
=={{header|C++}}==
{{trans|C#}}
<
#include <tuple>
Line 378:
return 0;
}</
{{out}}
<pre>(1.200000 + 0.030000) is in the range (1.2299998998641968, 1.2300001382827759)</pre>
Line 384:
=={{header|D}}==
{{trans|Kotlin}}
<
auto safeAdd(T)(T a, T b)
if (isFloatingPoint!T) {
Line 399:
auto r = safeAdd(a, b);
writefln("(%s + %s) is in the range %0.16f .. %0.16f", a, b, r.d, r.u);
}</
{{out}}
Line 416:
E currently inherits [[Java]]'s choices of [[IEEE]] floating point behavior, including round to nearest. Therefore, as in the Ada example, given ignorance of the actual direction of rounding, we must take one step away in ''both'' directions to get a correct interval.
<
require(a <= b)
def interval {
Line 432:
}
return interval
}</
<
# value: [2001.1399999999999, 2001.1400000000003]</
E provides only 64-bit "double precision" floats, and always prints them with sufficient decimal places to reproduce the original floating point number exactly.
Line 441:
=={{header|Forth}}==
{{trans|Tcl}}
<
s" m" add-lib
\c #include <math.h>
Line 455:
: savef+ ( F: r1 r2 -- r3 r4 ) \ r4 <= r1+r2 <= r3
f+ fdup fstepup fswap fstepdown ;</
{{output}}
<pre>
Line 465:
=={{header|Go}}==
{{trans|Tcl}}
<
import (
Line 494:
a, b := 1.2, .03
fmt.Println(a, b, safeAdd(a, b))
}</
Output:
<pre>
Line 502:
=={{header|Hare}}==
{{trans|C}}
<
use math;
Line 528:
math::setround(orig);
return interval{a = r0, b = r1};
};</
{{out}}
<pre>
Line 539:
=={{header|J}}==
J uses 64 bit IEEE floating points, providing 53 binary digits of accuracy.
<
safeadd =. + (-,+) +&err
0j15": 1.14 safeadd 2000.0 NB. print with 15 digits after the decimal
2001.139999999999873 2001.140000000000327</
=={{header|Java}}==
{{trans|Kotlin}}
<
private static double stepDown(double d) {
return Math.nextAfter(d, Double.NEGATIVE_INFINITY);
Line 565:
System.out.printf("(%.2f + %.2f) is in the range %.16f..%.16f", a, b, result[0], result[1]);
}
}</
{{out}}
<pre>(1.20 + 0.03) is in the range 1.2299999999999998..1.2300000000000002</pre>
Line 571:
=={{header|Julia}}==
Julia has the IntervalArithmetic module, which provides arithmetic with defined precision along with the option of simply computing with actual two-number intervals:
<syntaxhighlight lang="julia">
julia> using IntervalArithmetic
Line 591:
julia> a + b
[0.399999, 0.900001]
</syntaxhighlight>
=={{header|Kotlin}}==
{{trans|Tcl}}
<
fun stepDown(d: Double) = Math.nextAfter(d, Double.NEGATIVE_INFINITY)
Line 607:
val b = 0.03
println("($a + $b) is in the range ${safeAdd(a, b)}")
}</
{{out}}
Line 619:
=={{header|Nim}}==
{{trans|C}}
<
proc `++`(a, b: float): tuple[lower, upper: float] =
Line 638:
echo x.ff, " + ", y.ff, " ="
echo " [", d.ff, ", ", u.ff, "]"
echo " size ", (u - d).ff, "\n"</
Output:
<pre>1.0000000000000000 + 2.0000000000000000 =
Line 658:
=={{header|Perl}}==
There are ways avoid this whole problem (e.g. <code>Math::Decimal</code>), but another module suffices for this task.
<
use warnings;
use Data::IEEE754::Tools <nextUp nextDown>;
Line 668:
}
printf "%.17f (%.17f, %.17f)\n", safe_add (1/9,1/7);</
{{out}}
<pre>0.25396825396825395 (0.25396825396825390, 0.25396825396825401)</pre>
Line 680:
Not surprisingly you have to get a bit down and dirty to manage this sort of stuff in Phix.
<!--<
<span style="color: #008080;">include</span> <span style="color: #000000;">builtins</span><span style="color: #0000FF;">\</span><span style="color: #000000;">VM</span><span style="color: #0000FF;">\</span><span style="color: #000000;">pFPU</span><span style="color: #0000FF;">.</span><span style="color: #000000;">e</span> <span style="color: #000080;font-style:italic;">-- :%down53 etc</span>
Line 734:
<span style="color: #7060A8;">printf</span><span style="color: #0000FF;">(</span><span style="color: #000000;">1</span><span style="color: #0000FF;">,</span><span style="color: #008000;">" size %.16g\n\n"</span><span style="color: #0000FF;">,</span> <span style="color: #000000;">high</span> <span style="color: #0000FF;">-</span> <span style="color: #000000;">low</span><span style="color: #0000FF;">);</span>
<span style="color: #008080;">end</span> <span style="color: #008080;">for</span>
<!--</
{{out}}
<pre>
Line 763:
Python doesn't include a module that returns an interval for safe addition, however, it does [http://docs.python.org/library/math.html#math.fsum include] a routine for performing additions of floating point numbers whilst preserving precision:
<
0.9999999999999999
>>> from math import fsum
>>> fsum([.1, .1, .1, .1, .1, .1, .1, .1, .1, .1])
1.0</
=={{header|Racket}}==
<syntaxhighlight lang="racket">
#lang racket
Line 799:
;; literals
(bf+ (bf "1.14") (bf "2000.0"))
</syntaxhighlight>
=={{header|Raku}}==
Line 831:
<syntaxhighlight lang="raku"
say "Error: " ~ (2**-53).Num;
Line 858:
say "\nRaku default stringification for 1.5**63:\n" ~ $rat; # standard stringification
say "\nRat::Precise stringification for 1.5**63:\n" ~$rat.precise; # full precision
}</
{{out}}
<pre>Floating points: (Nums)
Line 887:
constrained by how much virtual memory is (or can be) allocated.
<br>Eight million digits seems about a practical high end, however.
<
y=digits() /*sets Y to existing number of decimal digits.*/
numeric digits y + y%10 /*increase the (numeric) decimal digits by 10%.*/</
=={{header|Ruby}}==
Line 898:
When adding <tt>BigDecimal</tt> values, <tt>a + b</tt> is always safe. This example uses <tt>a.add(b, prec)</tt>, which is not safe because it rounds to <tt>prec</tt> digits. This example computes a safe interval by rounding to both floor and ceiling.
<
require 'bigdecimal/util' # String#to_d
Line 920:
["0.1", "0.00002"],
["0.1", "-0.00002"],
].each { |a, b| puts "#{a} + #{b} = #{safe_add(a, b, 3)}" }</
Output: <pre>1 + 2 = 0.3E1..0.3E1
Line 928:
=={{header|Scala}}==
<
val (a, b) = (1.2, 0.03)
val result = safeAdd(a, b)
Line 940:
println(f"($a%.2f + $b%.2f) is in the range ${result.head}%.16f .. ${result.last}%.16f")
}</
=={{header|Swift}}==
<
let b = 0.03
print("\(a) + \(b) is in the range \((a + b).nextDown)...\((a + b).nextUp)")</
{{out}}
Line 955:
<br>
{{libheader|critcl}}
<
package provide stepaway 1.0
critcl::ccode {
Line 966:
critcl::cproc stepdown {double value} double {
return nextafter(value, -DBL_MAX);
}</
With that package it's then trivial to define a "safe addition" that returns an interval as a list (lower bound, upper bound).
<
proc safe+ {a b} {
set val [expr {double($a) + $b}]
return [list [stepdown $val] [stepup $val]]
}</
=={{header|Transd}}==
<
MainModule : {
Line 988:
(lout "(+ " a " " b ") is in the range: " prec: 20 (safeAdd a b))
)
}</
<pre>
(+ 1.2 0.03) is in the range: [1.2299999999999997602, 1.2300000000000002043]
Line 997:
However, if Wren is embedded in a suitable C program, then we can ask the latter to call it for us and pass back the result.
<
class Interval {
construct new(lower, upper) {
Line 1,021:
var a = 1.2
var b = 0.03
System.print("(%(a) + %(b)) is in the range %(Interval.safeAdd(a,b))")</
which we embed in the following C program and run it (using GCC 9.3.0).
Note that Wren's built-in print statement never displays numbers with more than 14 digit accuracy. However, if the interval bounds had been displayed directly from C with 17 digit accuracy (the maximum for the 'double' type), there would have been a marginal difference between them, namely: [1.2299999999999998, 1.2300000000000002].
<
#include <stdio.h>
#include <math.h>
Line 1,108:
free(script);
return 0;
}</
{{out}}
|