Numeric error propagation: Difference between revisions

Specification of a generic type Approximation.Number, providing all the operations required to solve the task ... and some more operations, for completeness.

 

type Real is digits <>;

with function Sqrt(X: Real) return Real;

Sigma: Real;

end record;

 

The implementation:

 

 

package body Approximation is

end Put;

 

 

Instantiating the package with Float operations, to compute the distance:

 

 

procedure Test_Approximations is

((X1-X2)**2 + (Y1 - Y2)**2)**0.5,

Sigma_Fore => 1);

 

Output:

=={{header|ALGOL 68}}==

{{works with|ALGOL 68G|Any - tested with release 2.8.3.win32}}

MODE UNCERTAIN = STRUCT( REAL v, uncertainty );

 

UNCERTAIN d = ( ( ( x1 - x2 ) ^ 2 ) + ( y1 - y2 ) ^ 2 ) ^ 0.5;

 

{{out}}

<pre>

=={{header|C}}==

Rewrote code to make it more compact and added a nice formatting function for imprecise values so that they are printed out in a technically correct way i.e. with the symbol '±' . Output pasted after code.

#include <stdlib.h>

#include <string.h>

return 0;

}

 

<pre>

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

numeric_error.hpp

 

#include <cmath>

private:

double v, s;

numeric_error.cpp

#include <iostream>

#include "numeric_error.hpp"

 

return EXIT_SUCCESS;

{{out}}

<pre>111.803398874989 ±2.938366893361</pre>

 

=={{header|D}}==

 

const struct Imprecise {

writefln("Point p2: (%s, %s)", p2[0], p2[1]);

writeln("Distance(p1, p2): ", distance(p1, p2));

{{out}}

<pre>Point p1: (I(value=100, delta=1.1), I(value=50, delta=1.2))

 

=={{header|F_Sharp|F#}}==

let abs (x : float) = System.Math.Abs x

let pow = System.Math.Pow

 

printfn "Distance: %A" ((((x1 %- x2) %^ 2.) %+ ((y1 %- y2) %^ 2.)) %^ 0.5)

{{out}}

<pre>Distance: 111.80 ±2.49</pre>

 

=={{header|Fortran}}==

===Direct calculation===

REAL X1, Y1, X2, Y2 !The co-ordinates.

REAL X1E,Y1E,X2E,Y2E !Their standard deviation.

WRITE (6,2) D,C,E !Ahh, the relief.

2 FORMAT ("Distance",F6.1,A1,F4.2) !Sizes to fit the example.

This is old-style Fortran, except for the CHARACTER variable caused by problems with character codes and their screen glyphs. As can be seen, the formulae invite mistakes which is why there is no attempt to produce a single arithmetic expression for the result and its error estimate. Further, rather than attempt to emplace appropriate instances of the formula for a value raised to some power (squaring, and square root), risking all manner of misthinks, a function to do so is prepared, here using Fortran's "arithmetic statement function" protocol, expressly intended for such situations. And the results are...

<pre>

 

===More general===

INTEGER VSP,VMAX !Do so with an arithmetic stack.

PARAMETER (VMAX = 28) !Surely sufficient.

WRITE (VOUT,1) "Empty",VSP !My stack may be empty.

ELSE !But normally, it is not.

END IF !So much for protection.

END SUBROUTINE VSHOW !Could reveal all the stack...

WRITE (6,2) STACKV(1),PM,STACKE(1) !Ahh, the relief.

2 FORMAT ("Distance",F6.1,A1,F4.2) !Sizes to fit the example.

This is closer to the idea of extending the language to supply additional facilities. At the cost of hand-compiling the arithmetic expression into a sequence of pseudo-machine code subroutines, it is apparent that the mind-tangling associated error formulae need no longer be worried over. The various arithmetic subroutines have to be coded correctly with careful attention to the V or E statements in fact being for the V and E terms (cut&paste followed by inadequate adjustment the culprit here: such mistakes are less likely when using a card punch because copying is more troublesome), but this is a straightforward matter of checking. And indeed the VPOW(2) routine has the same effect as VSQUARE. Output:

<pre>

(100.0±1.1, 50.0±1.2)

(200.0±2.2,100.0±2.3)

 

Although the example source file uses the F90 MODULE protocol, this is merely for convenience. Otherwise, there would be the same collection of routines, all sharing a COMMON work area. It would be easy enough to prepare an interactive calculator using this scheme so that different calculations (and data) could be more easily experimented with. Inspection might suggest a return to the laboratory in order to measure some factor with greater precision because its error proves to be the largest contributor to the spread in the result.

 

===Fortran 90 ''et seq''.===

REAL VALUE

REAL SD

END TYPE POINT

TYPE(POINT) P1,P2

 

=={{header|Go}}==

Variance from task requirements is that the following does not "extend the language." It simply defines a type with associated functions and methods as required to solve the remainder of the task.

 

import (

fmt.Println("d: ", d.n)

fmt.Println("error:", d.errorTerm())

Output:

<pre>

 

=={{header|Haskell}}==

 

instance (Floating a) => Num (Error a) where

x2 = Error 200 2.2

y2 = Error 100 2.3

{{out}}

<pre>Error {value = 111.80339887498948, uncertainty = 2.4871670631463423}</pre>

The following solution works in both languages.

 

 

procedure main(a)

return (numeric(a)^numeric(b)) |

num(f := a.val^numeric(b), abs(f*b*(a.err/a.val)))

 

The output is:

First, we will need some utilities. <code>num</code> will extract the number part of a number. <code>unc</code> will extract the uncertainty part of a number, and will also be used to associate uncertainty with a number. <code>dist</code> will compute the distance between two numbers (which is needed for multiplicative uncertainty).

 

unc=: {:@}."1 : ,.

 

Note that if a number has no uncertainty assigned to it, we assume the uncertainty is zero.

Jumping into the example values, for illustration purposes:

 

y1=: 50 unc 1.2

 

x2=: 200 unc 2.2

 

Above, we see <code>unc</code> being used to associate a number with its uncertainty. Here's how to take them apart again:

 

100

unc x1

 

Note that these operations "do the right thing" for normal numbers:

 

100

unc 100

 

And, a quick illustration of the distance function:

 

Next, we need to define our arithmetic operations:

 

sub=: -&num unc dist&unc

mul=: *&num unc |@(*&num * dist&(unc%num))

div=: %&num unc |@(%&num * dist&(unc%num))

 

Finally, our required example:

 

 

=={{header|Java}}==

private double value;

private double error;

System.out.println(x1);

}

{{out}}

<pre>111.80339887498948±2.4871670631463423</pre>

 

=={{header|Kotlin}}==

{{trans|Java}}

 

data class Approx(val ν: Double, val σ: Double = 0.0) {

val y2 = Approx(100.0, 2.3)

println(((x1 - x2).pow(2.0) + (y1 - y2).pow(2.0)).pow(0.5))

{{out}}

<pre>111.80339887498948 ±2.4871670631463423</pre>

 

PlusMinus /: a_ ± σa_ + b_ ± σb_ := N[(a + b) ± Norm@{σa, σb}];

PlusMinus /: c_?NumericQ (a_ ± σa_) := N[c a ± Abs[c σa]];

PlusMinus /: (a_ ± σa_) (b_ ± σb_) := N[a b ± (a b Norm@{σa/a, σb/b})^2];

y1 = 50 ± 1.2;

x2 = 200 ± 2.2;

y2 = 100 ± 2.3;

{{Out}}

<pre>111.803 ± 2.48717</pre>

=={{header|PARI/GP}}==

''This is a work-in-progress.''

sub(a,b)=if(type(a)==type(b), [a[1]-b[1],a[2]+b[2]], if(type(a)=="t_VEC",a-[b,0],[a,0]-b));

mult(a,b)=if(type(a)=="t_VEC", if(type(b)=="t_VEC", [a[1]*b[1], abs(a[1]*b[1])*sqrt(norml2([a[2]/a[1],b[2]/b[1]]))], [b*a[1], abs(b)*a[2]]), [a*b[1], abs(a)*b[2]]);

pow(a,b)=[a[1]^b, abs(a[1]^b*b*a[2]/a[1])];

x1=[100,1.1];y1=[50,1.2];x2=[200,2.2];y2=[100,2.3];

 

=={{header|Perl}}==

package ErrVar;

 

}

 

use overload

;

 

# make a variable with mean value and a list of coefficient to

# variables providing independent errors

 

# mean value of the var, or just the input if it's not of this class

 

# return variance index array

 

}

 

}

 

 

# to determine if a var is probably zero. we use 1σ here

}

 

}

 

}

 

}

 

}

 

}

 

}

 

}

 

}

 

}

 

}

 

 

package main;

 

my $z1 = sqrt(($x1 - $x2) ** 2 + ($y1 - $y2) ** 2);

 

# this is not for task requirement

my $a = $x1 + $x2;

my $b = $y1 - 2 * $x2;

{{out}}

<pre>distance: 111.803±2.49

 

 

=={{header|PicoLisp}}==

 

The overloaded +, -, * and / operators look a bit complicated, because they must handle an arbitrary number of arguments to be compatible with the standard operators.

(load "@lib/math.l")

 

(round (car N) 10)

" ± "

Test:

(**

(+ (** (- X1 X2) 2.0) (** (- Y1 Y2) 2.0))

(unc 50. 1.2)

(unc 200. 2.2)

Output:

<pre>Distance: 111.8033988750 ± 2.48716706</pre>

 

=={{header|Python}}==

import math

p1, p2 = (x1, y1), (x2, y2)

print("Distance between points\n p1: %s\n and p2: %s\n = %r" % (

 

;Sample output:

=={{header|Racket}}==

{{trans|Mathematica}}

 

(struct ± (x dx) #:transparent

(define x2 (± 200 2.2))

(define y2 (± 100 2.3))

 

{{output}}

<pre>111.803 ± 2.487</pre>

 

=={{header|REXX}}==

{{trans|Fortran}}

parse arg a b . /*obtain arguments from the CL*/

if a=='' | a=="," then a= '100±1.1, 50±1.2' /*Not given? Then use default.*/

parse var ax ax '±' axe; parse var bx bx '±' bxE /* " err " Ax and Bx.*/

parse var ay ay '±' aye; parse var by by '±' byE /* " " " Ay " By.*/

say /*blank line for the eyeballs.*/

say 'distance=' D "±" #(D**2, .5, sqrt(xE**2 + yE**2)) /*display " " " */

/*──────────────────────────────────────────────────────────────────────────────────────*/

#: procedure; arg x,p,e; if p=.5 then z=1/sqrt(abs(x)); else z=abs(x)**(p-1); return p*e*z

m.=9; do j=0 while h>9; m.j=h; h=h%2+1; end /*j*/

do k=j+5 to 0 by -1; numeric digits m.k; g=(g+x/g)*.5; end /*k*/

<pre>

A point (x,y)= 100 ± 1.1, 50 ± 1.2

 

=={{header|Ruby}}==

def initialize(number, error)

@num = number

y2 = NumberWithUncertainty.new(100, 2.3)

 

<pre>111.803398874989 ± 2.48716706314634</pre>

 

=={{header|Scala}}==

{{trans|Kotlin}}

 

class Approx(val ν: Double, val σ: Double = 0.0) {

val y2 = Approx(100.0, 2.3)

println(√(((x1 - x2)^2.0) + ((y1 - y2)^2.0))) // => 111.80339887498948 ±2.938366893361004

{{out}}

<pre>111.80339887498948 ±2.938366893361004</pre>

 

=={{header|Tcl}}==

{{works with|Tcl|8.6}}

Firstly, some support code for doing RAII-like things, evolved from code in the [[Quaternion type#Tcl|quaternion]] solution:

oo::class create RAII-support {

constructor {} {

try $script on ok msg {$msg return}

] [uplevel 1 {namespace current}]] {*}$vals

The implementation of the number+error class itself:

variable N E

constructor {number {error 0.0}} {

 

export + - * / **

Demonstrating:

set x2 [Err 200 2.2]

set y1 [Err 50 1.2]

[[[$x1 - $x2] ** 2] + [[$y1 - $y2] ** 2]] ** 0.5

}]

Output:

<pre>

d = 111.80339887498948 ± 2.4871670631463423

</pre>