Małe, nieprzewidywalne wyniki w seriach modelu deterministycznego

10

Mam spory model (~ 5000 linii) napisany w C. Jest to program szeregowy, bez generowania liczb losowych w dowolnym miejscu. Wykorzystuje bibliotekę FFTW dla funkcji korzystających z FFT - nie znam szczegółów implementacji FFTW, ale zakładam, że funkcje w niej również są deterministyczne (poprawcie mnie, jeśli się mylę).

Problem, którego nie rozumiem, polega na tym, że otrzymuję małe różnice w wynikach dla identycznych uruchomień na tym samym komputerze (ten sam kompilator, te same biblioteki).

Używam zmiennych o podwójnej precyzji i na przykład w wyniku wypisuję zmienną value: fprintf(outFID, "%.15e\n", value);lub
fwrite(&value, 1, sizeof(double), outFID);

I ciągle otrzymywałbym różnice, takie jak:
2.07843469652206 4 e-16 vs. 2.07843469652206 3 e-16

Spędziłem dużo czasu próbując zrozumieć, dlaczego tak jest. Początkowo myślałem, że jeden z moich układów pamięci zepsuł się i zamówiłem je i wymieniłem, ale bezskutecznie. Następnie spróbowałem uruchomić mój kod na maszynie Linuksa kolegi i dostaję różnice o tym samym charakterze.

Co może być tego przyczyną? Jest to teraz niewielka kwestia, ale zastanawiam się, czy jest to „wierzchołek góry lodowej” (poważny problem).

Myślałem, że opublikuję tutaj zamiast StackOverflow na wypadek, gdyby ktoś pracujący z modelami numerycznymi napotkał ten problem. Gdyby ktoś mógł rzucić na to światło, byłbym bardzo zobowiązany.

Kontynuacja komentarzy:
Christian Clason i Vikram: po pierwsze, dziękuję za uwagę na moje pytanie. Artykuły, które podłączyłeś, sugerują, że: 1. błędy zaokrąglania ograniczają dokładność, oraz 2. inny kod (np. Wprowadzanie pozornie nieszkodliwych instrukcji drukowania) może wpływać na wyniki aż do epsilon maszyny. Powinienem wyjaśnić, że nie porównuję efektów fwritei fprintffunkcji. Używam jednego LUB drugiego. W szczególności ten sam plik wykonywalny jest używany dla obu przebiegów. Po prostu stwierdzam, że problem występuje bez względu na to, czy używam fprintfOR fwrite.

Zatem ścieżka do kodu (i wykonywalna) jest taka sama, a sprzęt jest taki sam. Przy tych wszystkich czynnikach zewnętrznych utrzymywanych na stałym poziomie, skąd zasadniczo bierze się przypadkowość? Podejrzewałem, że bit się odwrócił, ponieważ uszkodzona pamięć nie zachowuje się trochę poprawnie, dlatego wymieniłem układy pamięci, ale to nie wydaje się być problemem tutaj, zweryfikowałem i wskazałeś. Mój program generuje tysiące tych podwójnie precyzyjnych liczb w jednym przebiegu i zawsze istnieje losowa garść, która ma losowe odwrócenie bitów.

Kontynuacja pierwszego komentarza Christiana Clasona: Dlaczego tym samym co 0 w zakresie precyzji maszyny? Najmniejsza liczba dodatnia dla podwójnego to 2.22e-308, więc czy nie powinno to być równe 0? Mój program generuje tysiące wartości w zakresie 10 ^ -16 (od 1e-15 do 8e-17) i widzieliśmy znaczące różnice w naszym projekcie badawczym, więc mam nadzieję, że nie patrzyliśmy na nonsensowne liczby.2)10-16

Kontynuacja # 2 :
Jest to wykres szeregów czasowych generowanych przez model, aby pomóc w dyskusji na temat odejścia w komentarzach. wprowadź opis zdjęcia tutaj

boxofchalk1
źródło
Witamy w SciComp.SE! Czy zdajesz sobie sprawę, że liczby zmiennoprzecinkowe mają ograniczoną dokładność - w szczególności, że w podwójnej precyzji jest równy zero do precyzji maszyny? Tak więc zgłaszane różnice nie są tak naprawdę znaczące i prawdopodobnie wynikają z niewielkich różnic w implementacjach dwóch nazwanych funkcji, co prowadzi do nieco odmiennego kodu maszynowego. 2)10-16
Christian Clason,
Pytasz, dlaczego twoja maszyna nie jest dokładniejsza niż precyzja maszyny. en.wikipedia.org/wiki/Machine_epsilon
Vikram
1
Zobacz inf.ethz.ch/personal/gander/Heisenberg/paper.html, aby znaleźć pokrewny przykład subtelnego wpływu ścieżek kodu na arytmetykę zmiennoprzecinkową. I, oczywiście, ece.uwaterloo.ca/~dwharder/NumericalAnalysis/02Numerics/Double/…
Christian Clason
1
Może być tak, że skalowanie twojego problemu jest takie, że poprawne odpowiedzi są rzędu . W każdym razie powinieneś być zorientowany ze względną dokładnością swoich rozwiązań. 10-16
Brian Borchers,
2
@ boxofchalk1 Z pewnością nie wyglądają jak hałas; jak powiedział Brian, jeśli wszystkie twoje dane są tego rzędu wielkości, możesz być w porządku (znowu, chodzi o względną dokładność). Aby się upewnić, możesz albo przeskalować problem, aby był w porządku , lub ponownie uruchomić kod z większą precyzją (patrz fftw.org/doc/Precision.html ). 1
Christian Clason,

Odpowiedzi:

9

Istnieją aspekty współczesnych systemów obliczeniowych, które z natury są niedeterministyczne, które mogą powodować tego rodzaju różnice. Dopóki różnice są bardzo małe w porównaniu z wymaganą dokładnością twoich rozwiązań, prawdopodobnie nie ma powodu, aby się tym martwić.

Przykład tego, co może pójść nie tak, na podstawie własnego doświadczenia. Rozważ problem obliczania iloczynu kropkowego dwóch wektorów x i y.

re=ja=1nxjayja

Obliczenie tego produktu kropkowego wymaga obliczenia produktów a następnie zsumowania wyników. Mnożenie zmiennoprzecinkowe powinno za każdym razem dawać dokładnie takie same wyniki. Jeśli uzupełnienia zmiennoprzecinkowe są obliczane za każdym razem w tej samej kolejności, suma powinna być identyczna. Ponieważ jednak dodawanie zmiennoprzecinkowe nie jest asocjacyjne, możesz uzyskać różne wyniki, jeśli iloczyn dwóch wektorów jest obliczany w taki sposób, że dodawanie jest wykonywane w innej kolejności. xjayja

Na przykład możesz najpierw obliczyć iloczyn dwóch wektorów jako

re=((x1y1)+(x2)y2)))+(x3)y3))

a następnie jako

re=(x1y1)+((x2)y2))+(x3)y3))) .

Jak to się mogło stać? Oto dwie możliwości.

  1. Obliczenia wielowątkowe na rdzeniach równoległych. Nowoczesne komputery zwykle mają 2, 4, 8 lub nawet więcej rdzeni procesorów, które mogą pracować równolegle. Jeśli kod używa równoległych wątków do obliczania iloczynu na wielu procesorach, wówczas dowolne zaburzenie systemu (np. Użytkownik przesunął mysz i jeden z rdzeni procesora musi przetworzyć ten ruch myszy przed powrotem do iloczynu) spowodować zmianę kolejności dodawania.

  2. Wyrównanie danych i instrukcji wektorowych. Nowoczesne procesory Intel mają specjalny zestaw instrukcji, które mogą działać (na przykład) dla liczb zmiennoprzecinkowych jednocześnie. Te instrukcje wektorowe działają najlepiej, jeśli dane są wyrównane do 16 bajtów granic. Zwykle pętla iloczynu dzieli dane na sekcje po 16 bajtów (4 zmienne na raz.) Jeśli ponownie uruchomisz kod po raz drugi, dane mogą zostać wyrównane w inny sposób z 16-bajtowymi blokami pamięci, dzięki czemu dodawane są wykonywane w innej kolejności, co skutkuje inną odpowiedzią.

Możesz rozwiązać punkt 1, uruchamiając kod jako pojedynczy wątek i wyłączając wszystkie przetwarzanie równoległe. Możesz zająć się punktem 2, wymagając alokacji pamięci, aby wyrównać bloki pamięci (zwykle robisz to, kompilując kod za pomocą przełącznika, takiego jak -align.) Jeśli kod nadal daje różne wyniki, istnieją inne możliwości wyszukiwania w.

Ta dokumentacja firmy Intel omawia problemy, które mogą prowadzić do braku powtarzalności wyników w bibliotece Intel Math Kernel Library. Kolejny dokument firmy Intel omawiający przełączniki kompilatora do użycia z kompilatorami Intela.

Brian Borchers
źródło
Widzę, że uważasz, że twój kod działa w jednym wątku. Chociaż prawdopodobnie znasz dobrze swój kod, nie zdziwiłbym się, gdybyś wywoływał podprogramy (np. Procedury BLAS), które działają wielowątkowo. Powinieneś sprawdzić, aby zobaczyć dokładnie, jakich bibliotek używasz. Możesz także użyć narzędzi do monitorowania systemu, aby sprawdzić zużycie procesora.
Brian Borchers,
1
lub, jak stwierdzono, biblioteka FFTW ...
Christian Clason
@BrianBorchers, dziękuję. Przykład losowości wynikający z niesocjacyjnego charakteru dodawania zmiennoprzecinkowego jest pouczający. Christian Clason poruszył drugorzędną kwestię dotyczącą tego, czy wyniki mojego modelu są znaczące, biorąc pod uwagę wielkość liczb - może to być poważny problem, jeśli ma rację (i dobrze go rozumiem), więc teraz patrzę na to.
boxofchalk1,
2

Wspomniana biblioteka FFTW może działać w trybie niedeterministycznym.

Jeśli używasz trybu FFTW_MEASURE lub FFTW_PATIENT, programy sprawdzają w czasie wykonywania, które wartości parametrów działają najszybciej, a następnie będą używać tych parametrów w całym programie. Ponieważ czas pracy oczywiście nieco się waha, parametry będą różne, a wynik przekształceń Fouriera będzie niedeterministyczny. Jeśli chcesz deterministyczne FFTW, użyj trybu FFTW_ESTIMATE.

eimrek
źródło
1

Chociaż prawdą jest, że zmiany kolejności oceny wyrażeń mogą bardzo dobrze wystąpić ze względu na scenariusze przetwarzania wielordzeniowego / wielowątkowego, nie zapominaj, że może istnieć (nawet jeśli jest to długa szansa) jakaś usterka w projektowaniu sprzętu. Pamiętasz problem Pentium FDIV? (Zobacz https://en.wikipedia.org/wiki/Pentium_FDIV_bug ). Jakiś czas temu pracowałem nad oprogramowaniem komputerowym do symulacji obwodów analogowych. Część naszej metodologii obejmowała opracowanie pakietów testów regresji, które działałyby przeciwko nocnym wersjom oprogramowania. W przypadku wielu opracowanych modeli metody iteracyjne (np. Newton-Raphson ( https://en.wikipedia.org/wiki/Newton%27s_method) i Runge-Kutta) były szeroko stosowane w algorytmach symulacyjnych. W przypadku urządzeń analogowych często zdarza się, że artefakty wewnętrzne, takie jak napięcia, prądy itp., Mają wyjątkowo małe wartości liczbowe. Wartości te, jako część procesu symulacji, zmieniają się stopniowo w czasie (symulowanym). Skala tych zmian może być bardzo mała i często obserwowaliśmy, że kolejne operacje FPU na takich wartościach delta graniczyły z progiem „szumu” precyzji FPU (pływanie 64-bitowe ma 53-bitową mantysę, IIRC). To, w połączeniu z faktem, że często musieliśmy wprowadzać do modeli kod rejestrujący „PrintF”, aby umożliwić debugowanie (ach, dobre stare dni!), Praktycznie gwarantowało sporadyczne wyniki, na co dzień! Więc co' czy to wszystko znaczy? W takich okolicznościach należy spodziewać się różnic, a najlepszą rzeczą jest zdefiniowanie i wdrożenie sposobu decydowania (wielkość, częstotliwość, trend itp.), Kiedy / jak je zignorować.

Jim
źródło
Dziękuję Jim za wgląd. Jakieś pomysły na to, jakie fundamentalne zjawiska spowodowałyby takie „artefakty wewnętrzne”? Myślałem, że mogą to być zakłócenia elektromagnetyczne, ale wtedy znaczące bity również by na to wpłynęły, prawda?
boxofchalk1,
1

Podczas gdy zaokrąglanie zmiennoprzecinkowe z operacji asynchronicznych może być problemem, podejrzewam, że jest to coś bardziej banalnego. Zastosowanie niezainicjowanej zmiennej, która dodaje losowość do twojego kodu deterministycznego. Jest to częsty problem, który jest często pomijany przez programistów, ponieważ podczas uruchamiania w trybie debugowania wszystkie zmienne są inicjowane na 0 w momencie deklaracji. Gdy pamięć nie jest uruchomiona w trybie debugowania, pamięć przypisana do zmiennej ma dowolną wartość, jaką pamięć miała przed przypisaniem. Pamięć nie jest zerowana przy przypisaniu jako optymalizacja. Jeśli dzieje się tak w twoim kodzie, łatwo będzie to naprawić, w mniejszym stopniu w kodzie bibliotek.

brent.payne
źródło