help-octave
[Top][All Lists]

## Re: what's the ULP in Octave

 From: Kai Torben Ohlhus Subject: Re: what's the ULP in Octave Date: Wed, 2 Dec 2020 11:50:14 +0900 User-agent: Mozilla/5.0 (X11; Linux x86_64; rv:78.0) Gecko/20100101 Thunderbird/78.5.0

```>> On Dec 1, 2020, at 7:55 AM, 康 珊 <kangshan0910@hotmail.com> wrote:
>>
>> Hi,
```
>> I wonder what's the ULP(unit of least precision) for those double operations such as pow in octave.
```

On 12/2/20 6:01 AM, Victor Norton via Help-octave wrote:
> Try this
>    octave> x = 1; eps = 1;
>    octave> while(x + eps > x)
>    > eps /= 10;
>    > end
>    octave> eps
>    eps =    1.0000e-16
>

```
There is one important change necessary. Octave's double is binary not decimal [1,2], thus it should be "eps /= 2;". Otherwise the program will give the correct answer.
```
```
However, there is an Octave function called "eps" [3] which gives the result as well. No need to write own code, but for clarification it can be beneficial to run some code.
```
```````````````````````````````````````````````
x = 1;
e = 1;
i = 0;
while (x + e > x)
printf ("1 + 2^(%3d) = %s\n", ...
i, num2hex (x + e))
e /= 2;
i--;
endwhile
printf ("1 + 2^(%3d) = %s  <-- 'e' too small\n", ...
i, num2hex (x + e))
```````````````````````````````````````````````

```
The output gives a hexadecimal representation of double numbers "1.0 + e" for shrinking values of "e":
```
1 + 2^(  0) = 4000000000000000
1 + 2^( -1) = 3ff8000000000000
1 + 2^( -2) = 3ff4000000000000
...
1 + 2^(-51) = 3ff0000000000002
1 + 2^(-52) = 3ff0000000000001
1 + 2^(-53) = 3ff0000000000000  <-- 'e' too small

Thus:

```
- 2^(-52) == eps() (approx. "2.2204e-16") is the smallest number that can be added to "1.0" to "change" this number (the ULP unit in the last place).
```
```
- 2^(-53) == eps()/2 (approx. "1.1102e-16") is the number that added to 1.0 is too small to "change" this number. The precision of double is not enough to store the exact result of this addition.
```
```
A final note: eps() is *relative* to 1.0. To obtain eps relative to another double number, say 5.0, call "eps (5.0) == eps() * 4 == 2^(-50)"
```
HTH,
Kai

[1] https://en.wikipedia.org/wiki/Unit_in_the_last_place
[2] https://en.wikipedia.org/wiki/Double-precision_floating-point_format
[3] https://www.octave.org/doc/v6.1.0/XREFeps.html

```