Kleine, unvorhersehbare Ergebnisse in Läufen eines deterministischen Modells

Ich habe ein großes Modell (~ 5000 Zeilen) in C geschrieben. Es ist ein serielles Programm, bei dem nirgendwo Zufallszahlen generiert werden. Es nutzt die FFTW-Bibliothek für Funktionen, die FFT verwenden. Ich kenne die Details der FFTW-Implementierung nicht, gehe jedoch davon aus, dass die darin enthaltenen Funktionen auch deterministisch sind (korrigieren Sie mich, wenn ich mich irre).

Das Problem, das ich nicht verstehen kann, ist, dass ich kleine Unterschiede in den Ergebnissen für identische Läufe auf demselben Computer (demselben Compiler, denselben Bibliotheken) erhalte.

Ich verwende Variablen mit doppelter Genauigkeit, und um das Ergebnis beispielsweise in Variablen auszugeben value, gebe ich Folgendes aus: fprintf(outFID, "%.15e\n", value);oder
fwrite(&value, 1, sizeof(double), outFID);

Und ich würde ständig Unterschiede bekommen wie:
2.07843469652206 4 e-16 vs. 2.07843469652206 3 e-16

Ich habe viel Zeit damit verbracht herauszufinden, warum das so ist. Anfangs dachte ich, einer meiner Speicherchips sei defekt, und ich habe sie ohne Erfolg bestellt und ersetzt. Anschließend habe ich auch versucht, meinen Code auf dem Linux-Computer eines Kollegen auszuführen, und ich erhalte Unterschiede derselben Art.

Was könnte das verursachen? Es ist jetzt ein kleines Problem, aber ich frage mich, ob es die "Spitze des Eisbergs" ist (ein ernstes Problem).

Ich dachte, ich würde hier anstelle von StackOverflow posten, falls jemand, der mit numerischen Modellen arbeitet, auf dieses Problem stoßen könnte. Wenn jemand Licht ins Dunkel bringen kann, wäre ich sehr dankbar.

Follow-up zu Kommentaren:
Christian Clason und Vikram: Zunächst einmal vielen Dank für Ihre Aufmerksamkeit auf meine Frage. Die Artikel, auf die Sie verlinkt haben, schlagen Folgendes vor: 1. Rundungsfehler schränken die Genauigkeit ein und 2. Unterschiedlicher Code (z. B. das Einführen scheinbar harmloser Druckanweisungen) kann die Ergebnisse bis zum Maschinen-Epsilon beeinflussen. Ich sollte klarstellen, dass ich die Effekte fwriteund fprintfFunktionen nicht vergleiche . Ich benutze das eine oder das andere. Insbesondere wird für beide Läufe dieselbe ausführbare Datei verwendet. Ich sage nur, dass das Problem auftritt, ob ich fprintfOR verwende fwrite.

Der Codepfad (und die ausführbare Datei) sind also identisch, und die Hardware ist identisch. Woher kommt die Zufälligkeit bei all diesen konstanten externen Faktoren im Grunde genommen? Ich vermutete, dass der Bit-Flip passiert ist, weil der fehlerhafte Speicher ein bisschen nicht richtig beibehalten wurde, weshalb ich die Speicherchips ausgetauscht habe, aber das scheint hier nicht das Problem zu sein, habe ich überprüft und Sie haben angegeben. Mein Programm gibt Tausende dieser Zahlen mit doppelter Genauigkeit in einem einzigen Lauf aus, und es gibt immer eine zufällige Handvoll mit zufälligen Bitflips.

$2\cdot 10^{-16}$

Follow-up Nr. 2 :
Dies ist eine grafische Darstellung der vom Modell ausgegebenen Zeitreihen, um die Ablegerdiskussion in den Kommentaren zu unterstützen.

Es gibt Aspekte moderner Computersysteme, die von Natur aus nicht deterministisch sind und solche Unterschiede verursachen können. Solange die Unterschiede im Vergleich zur erforderlichen Genauigkeit Ihrer Lösungen sehr gering sind, gibt es wahrscheinlich keinen Grund, sich darüber Sorgen zu machen.

Ein Beispiel dafür, was aufgrund meiner eigenen Erfahrung schief gehen kann. Betrachten Sie das Problem der Berechnung des Punktprodukts zweier Vektoren x und y.

$d=\sum_{i=1}^{n} x_{i}y_{i}$

$x_{i}y_{i}$

Beispielsweise könnten Sie zuerst das Produkt zweier Vektoren als berechnen

$d=((x_{1}y_{1})+(x_{2}y_{2}))+(x_{3}y_{3})$

und dann als

$d=(x_{1}y_{1})+((x_{2}y_{2})+(x_{3}y_{3}))$

Wie könnte das passieren? Hier sind zwei Möglichkeiten.

1. Multithread-Berechnungen auf parallelen Kernen. Moderne Computer verfügen normalerweise über 2, 4, 8 oder sogar mehr Prozessorkerne, die parallel arbeiten können. Wenn Ihr Code parallele Threads verwendet, um ein Punktprodukt auf mehreren Prozessoren zu berechnen, kann dies zu einer zufälligen Störung des Systems führen (z. B. wenn der Benutzer seine Maus bewegt hat und einer der Prozessorkerne diese Mausbewegung verarbeiten muss, bevor er zum Punktprodukt zurückkehrt) führen zu einer Änderung in der Reihenfolge der Ergänzungen.

2. Ausrichtung von Daten und Vektoranweisungen. Moderne Intel-Prozessoren verfügen über einen speziellen Befehlssatz, der beispielsweise für Gleitkommazahlen gleichzeitig ausgeführt werden kann. Diese Vektoranweisungen funktionieren am besten, wenn die Daten an 16-Byte-Grenzen ausgerichtet sind. In der Regel teilt eine Punktproduktschleife die Daten in Abschnitte von 16 Byte auf (jeweils 4 Gleitkommazahlen). Wenn Sie den Code ein zweites Mal ausführen, werden die Daten möglicherweise anders an den 16-Byte-Speicherblöcken ausgerichtet, sodass die Hinzufügungen erfolgen in einer anderen Reihenfolge durchgeführt, was zu einer anderen Antwort führt.

Sie können Punkt 1 ansprechen, indem Sie Ihren Code als einzelnen Thread ausführen und die gesamte parallele Verarbeitung deaktivieren. Sie können Punkt 2 adressieren, indem Sie eine Speicherzuweisung zum Ausrichten von Speicherblöcken benötigen (normalerweise würden Sie dies tun, indem Sie den Code mit einem Schalter wie -align kompilieren.) Wenn Ihr Code immer noch Ergebnisse liefert, die variieren, gibt es andere Möglichkeiten, nachzuschauen beim.

In dieser Dokumentation von Intel werden Probleme erläutert, die dazu führen können, dass die Ergebnisse mit der Intel Math Kernel Library nicht reproduzierbar sind. Ein weiteres Dokument von Intel, in dem Compiler-Switches für die Verwendung mit Intel-Compilern erläutert werden.

Die erwähnte FFTW-Bibliothek wird möglicherweise in einem nicht deterministischen Modus ausgeführt.

Wenn Sie den Modus FFTW_MEASURE oder FFTW_PATIENT verwenden, prüfen die Programme zur Laufzeit, welche Parameterwerte am schnellsten funktionieren, und verwenden diese Parameter dann im gesamten Programm. Da die Laufzeit offensichtlich etwas schwankt, sind die Parameter unterschiedlich und das Ergebnis der Fourier-Transformationen ist nicht deterministisch. Wenn Sie deterministisches FFTW möchten, verwenden Sie den Modus FFTW_ESTIMATE.

Zwar kann es aufgrund von Multi-Core- / Multi-Thread-Verarbeitungsszenarien durchaus zu Änderungen der Auswertungsreihenfolge von Ausdrucksterms kommen, aber vergessen Sie nicht, dass bei der Arbeit möglicherweise (auch wenn es sich um einen langen Schuss handelt) ein Hardware-Designfehler vorliegt. Erinnern Sie sich an das Pentium FDIV-Problem? (Siehe https://en.wikipedia.org/wiki/Pentium_FDIV_bug ). Vor einiger Zeit habe ich an einer PC-basierten Software zur Simulation analoger Schaltkreise gearbeitet. Ein Teil unserer Methodik bestand in der Entwicklung von Regressionstestsuiten, die wir gegen nächtliche Builds der Software ausführen würden. Bei vielen der von uns entwickelten Modelle wurden iterative Methoden angewendet (z. B. Newton-Raphson ( https://en.wikipedia.org/wiki/Newton%27s_method)) und Runge-Kutta) wurden ausgiebig in den Simulationsalgorithmen verwendet. Bei analogen Geräten treten häufig interne Artefakte wie Spannungen, Ströme usw. mit extrem kleinen numerischen Werten auf. Diese Werte werden im Rahmen des Simulationsprozesses schrittweise über die (simulierte) Zeit variiert. Das Ausmaß dieser Änderungen kann extrem gering sein, und wir haben oft beobachtet, dass nachfolgende FPU-Operationen mit solchen Delta-Werten an die "Rausch" -Schwelle der FPU-Genauigkeit grenzten (64-Bit-Floating hat eine 53-Bit-Mantisse, IIRC). Dies, zusammen mit der Tatsache, dass wir häufig "PrintF" -Protokollierungscode in Modelle einführen mussten, um das Debuggen (ah, die guten alten Tage!) Zu ermöglichen, garantierte praktisch täglich sporadische Ergebnisse! Na und' Ist das alles gemein? Unter solchen Umständen müssen Sie mit Unterschieden rechnen. Am besten definieren und implementieren Sie einen Weg, um zu entscheiden (Größe, Häufigkeit, Trend usw.), wann und wie Sie diese ignorieren.

Während Gleitkomma-Rundungen aus asynchronen Operationen das Problem sein können, vermute ich, dass es etwas Banaleres ist. Die Verwendung einer nicht initialisierten Variablen, die Ihrem ansonsten deterministischen Code Zufälligkeit hinzufügt. Dies ist ein häufiges Problem, das von Entwicklern häufig übersehen wird, da beim Ausführen im Debug-Modus alle Variablen bei der Deklaration auf 0 initialisiert werden. Wenn der Speicher nicht im Debug-Modus ausgeführt wird, hat der einer Variablen zugewiesene Speicher den Wert, den der Speicher vor der Zuweisung hatte. Der Speicher wird bei Zuweisung als Optimierung nicht auf Null gesetzt. Wenn dies in Ihrem Code geschieht, kann es leicht behoben werden, weniger im Bibliothekscode.