Dlaczego wynik tej wyraźnej obsady różni się od niejawnej?

Question

Dlaczego wynik tej wyraźnej obsady różni się od niejawnej?

#include <stdio.h> 

double a; 
double b; 
double c; 

long d; 

double e; 

int main() { 
    a = 1.0; 
    b = 2.0; 
    c = .1; 

    d = (b - a + c)/c; 
    printf("%li\n", d);  // 10 

    e = (b - a + c)/c; 
    d = (long) e; 
    printf("%li\n", d);  // 11 
    } 

Jeśli wykonam d = (długa) ((b - a + c)/c); Otrzymuję również 10. Dlaczego przypisanie podwójnemu stanowi różnicę?

Answer 1

podejrzewać, że różnica jest konwersja z 80-bitowej wartości zmiennoprzecinkowej na długi vs konwersji z 80-bitowej wartości zmiennoprzecinkowej na 64-bitowy jednego i następnie konwersję długi.

(Powodem dla 80 bitów wymyślanie w ogóle jest to, że jest to typowy precyzja wykorzystywane do rzeczywistej arytmetyki, a szerokość przestawne rejestrów punktowych.)

Załóżmy, że wynik 80-bit jest coś 10.999999999999999 - w konwersja z tego na długą wydajność 10. Jednak najbliższa 64-bitowa wartość zmiennoprzecinkowa do 80-bitowej wartości wynosi faktycznie 11,0, więc konwersja dwustopniowa kończy się dając 11.

EDYCJA: Aby nadać temu trochę więcej wagi ...

Oto program w języku Java, który wykorzystuje arytmetykę arytmetyki arbitralnej do wykonywania sam e obliczenia. Zwróć uwagę, że konwertuje podwójną wartość najbliższą 0.1 w BigDecimal - ta wartość to 0,1000000000000000055511151231257827021181583404541015625. (Innymi słowy, dokładny wynik obliczeń jest nie 11 anyway.)

import java.math.*; 

public class Test 
{ 
    public static void main(String[] args) 
    { 
     BigDecimal c = new BigDecimal(0.1d);   
     BigDecimal a = new BigDecimal(1d); 
     BigDecimal b = new BigDecimal(2d); 

     BigDecimal result = b.subtract(a) 
          .add(c) 
          .divide(c, 40, RoundingMode.FLOOR); 
     System.out.println(result); 
    } 
} 

Oto wynik:

10.9999999999999994448884876874217606030632 

Innymi słowy, jest to poprawne do około 40 cyfr po przecinku (droga więcej niż 64 lub 80-bitowy zmiennoprzecinkowy może obsłużyć).

Rozważmy teraz, jak wygląda ten numer w systemie binarnym. Nie mam żadnych narzędzi, które umożliwiałyby łatwą konwersję, ale znowu możemy pomóc w Javie. Zakładając znormalizowaną liczbę, część "10" kończy się na trzech bitach (jeden mniej niż jedenastu = 1011). Pozostawia to 60 bitów mantysy dla rozszerzonej precyzji (80 bitów) i 48 bitów dla podwójnej precyzji (64 bity).

Jaka jest największa liczba do 11 w każdej precyzji? Znowu użyjmy Java:

import java.math.*; 

public class Test 
{ 
    public static void main(String[] args) 
    { 
     BigDecimal half = new BigDecimal("0.5");   
     BigDecimal eleven = new BigDecimal(11); 

     System.out.println(eleven.subtract(half.pow(60))); 
     System.out.println(eleven.subtract(half.pow(48)));   
    } 
} 

Wyniki:

10.999999999999999999132638262011596452794037759304046630859375 
10.999999999999996447286321199499070644378662109375 

Tak, trzy numery mamy są:

Correct value: 10.999999999999999444888487687421760603063... 
11-2^(-60): 10.999999999999999999132638262011596452794037759304046630859375 
11-2^(-48): 10.999999999999996447286321199499070644378662109375 

Teraz wyszło najbliższą wartość do prawidłowego dla każdej precyzji - dla większej precyzji jest mniej niż 11. Zaokrąglij każdą z tych wartości do wartości długiej, a otrzymasz odpowiednio 10 i 11.

Mam nadzieję, że jest wystarczająco dużo dowodów, aby przekonać niedowiarków;)

Answer 2

Prosto kopiuj/wklej i skompilować Linux daje mi 11 dla obu stron. Dodanie d = (long) ((b - a + c)/c); daje również 11. To samo dotyczy OpenBSD.

Answer 3

codepad.org (gcc 4.1.2) odwraca wyniki twojego przykładu, podczas gdy w moim lokalnym systemie (gcc 4.3.2) dostaję 11 w obu przypadkach. Sugeruje to mi, że jest to kwestia zmiennoprzecinkowa. Alternatywnie może teoretycznie być obcięciem (b - a + c), które w kontekście całkowitym oceniałoby (2 - 1 + 0)/.1, co oznaczałoby 10, podczas gdy w kontekście zmiennoprzecinkowym (2.0 - 1.0 + 0.1)/.1 = 1.1/.1 = 11. To byłoby dziwne.

Answer 4

Here is a bunch of detail on floating point issues and a really good article. Ale nie wszystkie wartości zmiennoprzecinkowe mogą być reprezentowane przez określoną liczbę bitów (32-bitowych lub 64-bitowych lub cokolwiek innego). Jest to temat głęboki, ale który lubię, ponieważ przypomina mi Prof. Kahan. :)

Answer 5

Otrzymuję 10 & 11 na moim 32-bitowym systemie x86 Linux działającym z gcc 4.3.2.

Odpowiedni C/ASM jest tutaj:

26:foo.c   ****  d = (b - a + c)/c;            
    42       .loc 1 26 0 
    43 0031 DD050000    fldl b 
    43  0000 
    44 0037 DD050000    fldl a 
    44  0000 
    45 003d DEE9     fsubrp %st, %st(1) 
    46 003f DD050000    fldl c 
    46  0000 
    47 0045 DEC1     faddp %st, %st(1) 
    48 0047 DD050000    fldl c 
    48  0000 
    49 004d DEF9     fdivrp %st, %st(1) 
    50 004f D97DFA    fnstcw -6(%ebp) 
    51 0052 0FB745FA    movzwl -6(%ebp), %eax 
    52 0056 B40C     movb $12, %ah 
    53 0058 668945F8    movw %ax, -8(%ebp) 
    54 005c D96DF8    fldcw -8(%ebp) 
    55 005f DB5DF4    fistpl -12(%ebp) 
    56 0062 D96DFA    fldcw -6(%ebp) 
    57 0065 8B45F4    movl -12(%ebp), %eax 
    58 0068 A3000000    movl %eax, d 
    58  00 
    27:foo.c   **** 
    28:foo.c   ****  printf("%li\n", d);             
    59       .loc 1 28 0 
    60 006d A1000000    movl d, %eax 
    60  00 
    61 0072 89442404    movl %eax, 4(%esp) 
    62 0076 C7042400    movl $.LC3, (%esp) 
    62  000000 
    63 007d E8FCFFFF    call printf 
    63  FF 
    29:foo.c   ****  // 10               
    30:foo.c   **** 
    31:foo.c   ****  e = (b - a + c)/c;            
    64       .loc 1 31 0 
    65 0082 DD050000    fldl b 
    65  0000 
    66 0088 DD050000    fldl a 
    66  0000 
    67 008e DEE9     fsubrp %st, %st(1) 
    68 0090 DD050000    fldl c 
    68  0000 
    69 0096 DEC1     faddp %st, %st(1) 
    70 0098 DD050000    fldl c 
    70  0000 
    71 009e DEF9     fdivrp %st, %st(1) 
    72 00a0 DD1D0000    fstpl e 
    72  0000 
    32:foo.c   **** 
    33:foo.c   ****  d = (long) e;              
    73       .loc 1 33 0 
    74 00a6 DD050000    fldl e 
    74  0000 
    75 00ac D97DFA    fnstcw -6(%ebp) 
    76 00af 0FB745FA    movzwl -6(%ebp), %eax 
    77 00b3 B40C     movb $12, %ah 
    78 00b5 668945F8    movw %ax, -8(%ebp) 
    79 00b9 D96DF8    fldcw -8(%ebp) 
    80 00bc DB5DF4    fistpl -12(%ebp) 
    81 00bf D96DFA    fldcw -6(%ebp) 
    82 00c2 8B45F4    movl -12(%ebp), %eax 
    83 00c5 A3000000    movl %eax, d 
    83  00 

Odpowiedź pozostawiamy jako ćwiczenie dla czytelnika zainteresowanego.