Safe addition: Difference between revisions

=={{header|Ada}}==

An interval type based on Float:

Lower : Float;

Upper : Float;

Implementation of safe addition:

Result : constant Float := A + B;

begin

return (Float'Adjacent (0.0, Float'First), Float'Adjacent (0.0, Float'Last));

end if;

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

begin

Put (1.14 + 2000.0);

Sample output:

<pre> 2.00113989257813E+03, 2.00114013671875E+03</pre>

=={{header|AutoHotkey}}==

SetFormat, FloatFast, 0.20

err:=0.1**(SubStr(A_FormatFloat,3) > 15 ? 15 : SubStr(A_FormatFloat,3))

Return "[" a+b-err ","a+b+err "]"

=={{header|C}}==

=== C99 with fesetround() ===

#include <stdio.h> /* printf() */

}

return 0;

Output:

<pre>

{{works with|MinGW}}

#include <stdio.h> /* printf() */

}

return 0;

=== fpsetround() ===

{{works with|OpenBSD|4.8/amd64}}

#include <stdio.h> /* printf() */

}

return 0;

Output from OpenBSD: <pre>1 + 2 =

=={{header|C sharp|C#}}==

{{trans|Java}}

namespace SafeAddition {

}

}

{{out}}

<pre>(1.2 + 0.03) is in the range (1.23, 1.23)</pre>

=={{header|C++}}==

{{trans|C#}}

#include <tuple>

return 0;

{{out}}

<pre>(1.200000 + 0.030000) is in the range (1.2299998998641968, 1.2300001382827759)</pre>

=={{header|D}}==

{{trans|Kotlin}}

auto safeAdd(T)(T a, T b)

if (isFloatingPoint!T) {

auto r = safeAdd(a, b);

writefln("(%s + %s) is in the range %0.16f .. %0.16f", a, b, r.d, r.u);

{{out}}

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 {

}

return interval

E provides only 64-bit "double precision" floats, and always prints them with sufficient decimal places to reproduce the original floating point number exactly.

=={{header|Forth}}==

{{trans|Tcl}}

s" m" add-lib

\c #include <math.h>

: savef+ ( F: r1 r2 -- r3 r4 ) \ r4 <= r1+r2 <= r3

{{output}}

<pre>

=={{header|Go}}==

{{trans|Tcl}}

import (

a, b := 1.2, .03

fmt.Println(a, b, safeAdd(a, b))

Output:

<pre>

=={{header|Hare}}==

{{trans|C}}

use math;

math::setround(orig);

return interval{a = r0, b = r1};

{{out}}

<pre>

=={{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

=={{header|Java}}==

{{trans|Kotlin}}

private static double stepDown(double d) {

return Math.nextAfter(d, Double.NEGATIVE_INFINITY);

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>

=={{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:

julia> using IntervalArithmetic

julia> a + b

[0.399999, 0.900001]

=={{header|Kotlin}}==

{{trans|Tcl}}

fun stepDown(d: Double) = Math.nextAfter(d, Double.NEGATIVE_INFINITY)

val b = 0.03

println("($a + $b) is in the range ${safeAdd(a, b)}")

{{out}}

=={{header|Nim}}==

{{trans|C}}

proc `++`(a, b: float): tuple[lower, upper: float] =

echo x.ff, " + ", y.ff, " ="

echo " [", d.ff, ", ", u.ff, "]"

Output:

<pre>1.0000000000000000 + 2.0000000000000000 =

=={{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>;

}

{{out}}

<pre>0.25396825396825395 (0.25396825396825390, 0.25396825396825401)</pre>

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>

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

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])

=={{header|Racket}}==

#lang racket

;; literals

(bf+ (bf "1.14") (bf "2000.0"))

=={{header|Raku}}==

say "Error: " ~ (2**-53).Num;

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)

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.*/

=={{header|Ruby}}==

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

["0.1", "0.00002"],

["0.1", "-0.00002"],

Output: <pre>1 + 2 = 0.3E1..0.3E1

=={{header|Scala}}==

val (a, b) = (1.2, 0.03)

val result = safeAdd(a, b)

println(f"($a%.2f + $b%.2f) is in the range ${result.head}%.16f .. ${result.last}%.16f")

=={{header|Swift}}==

let b = 0.03

{{out}}

<br>

{{libheader|critcl}}

package provide stepaway 1.0

critcl::ccode {

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 : {

(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]

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) {

var a = 1.2

var b = 0.03

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>

free(script);

return 0;

{{out}}