Programm läuft nicht korrekt

Ramanujan

Ja, man quasi sagen, dass ich das richtige Ergebnis kenne.

Das Linux-Ergebnis ist richtig.
Das Windows-Ergebnis ist zu klein (allerdings nur für große Matrizen).

Ramanujan

SeppJ schrieb:

Möglicherweise ein Programmierfehler.

Und was für ein Programmierfehler könnte dafür sorgen, dass es einmal funktioniert und einmal nicht?

Athar

Ramanujan schrieb:

Und was für ein Programmierfehler könnte dafür sorgen, dass es einmal funktioniert und einmal nicht?

Viele, darunter auch uninitialisierte Variablen und ungültige Speicherzugriffe, die valgrind in vielen Fällen aufdecken kann.

SeppJ

Athar schrieb:

Ramanujan schrieb:

Und was für ein Programmierfehler könnte dafür sorgen, dass es einmal funktioniert und einmal nicht?

Viele, darunter auch uninitialisierte Variablen und ungültige Speicherzugriffe, die valgrind in vielen Fällen aufdecken kann.

Genau. valgrind ist ziemlich gut. Wenn das Programm valgrind sauber ist (mit allen Prüfungen!), dann kann man sich das Problem mal näher ansehen. Wenn valgrind nichts findet, dann ist zwar ein Fehler dadurch nicht ausgeschlossen, aber es ist zumindest nichts dummes, was man auch sofort hätte entdecken können. Außerdem: Debugmodus der STL aktivieren:
http://gcc.gnu.org/onlinedocs/libstdc++/manual/bk01pt03ch17s03.html#debug_mode.using.mode
(falls du nicht GCC benutzt, gibt es für deine Implementierung sicherlich etwas ähnliches)

Bis du keinen solchen valgrind-check gemacht hast, muss man erst einmal von einem Programmierfehler ausgehen. Es mag auch andere Möglichkeiten geben, die Symptome zu erklären, aber dies ist einfach eine sehr wahrscheinliche Möglichkeit und sollte daher zuerst geprüft werden.

Captain Obvious

Wie groß ist der Fehler? Vielleicht wird unter Windows mit geringerer Floating Point-Genauigkeit kompiliert/gerechnet.

Ramanujan

So Faktor 25 bzw. 5 (je nachdem welche Werte man sich anguckt) für n = 10^6 Zeilen. Das ist eigentlich zu viel, als das es an der Floating Point-Genauigkeit liegen könnte. Aber das ist eine interessante Idee. Wenn ich den Code ein wenig umschreibe und vereinfache, könnte ich mal mit long long rechnen.

Das mit valgrind werde ich auch mal ausprobieren.

Gerade habe ich nur keine Zeit dafür. Erstmal danke für die bisherigen Tipps.

cooky451

Ich dachte der Typ ist int? Das hat nichts mit Floating-Point zu tun. Integer sind absolut präzise, auf jedem System.

Captain Obvious

cooky451 schrieb:

Ich dachte der Typ ist int? Das hat nichts mit Floating-Point zu tun. Integer sind absolut präzise, auf jedem System.

Was sollen denn cos und sin mit ints anfangen?

cooky451

Hab ich gar nicht gelesen.

dd++ 0

mit "gcc -Q -v (inputfile)" bekommt man alle default Settings des jeweiligen gcc angezeigt, darunter sind sicherlich auch einige, die Einfluß auf die Floating Point Arithmetik haben (CPU Architektur, MMX, SSE usw.)

Ramanujan

Ich hab mich unklar ausgedrückt:

Die Größe der Matrix ist vom Typ int. Die Daten der Matrix sind double. Da die Matrixklasse ein Template ist, kann ich auch mal mit long long durchrechnen. Die Stellen, in denen sin, cos und sqrt vorkommen kann ich weglassen. Division kann ich auch vermeiden.

Mit long long ist das Ergebnis in Ordnung.

Die vereinfachte Version funktioniert allerdings auch mit double korrekt...

Edit:
valgrind gibt es wohl nicht für Windows, deshalb hab ich auf Dr. Memory zurückgegriffen:

Dr. Memory version 1.5.0 build 5 built on Aug 31 2012 16:19:51
Application cmdline: ""C:\Program Files\cpp-eclipse\Projekte\Linear Algebra\Debug\Linear Algebra.exe""
Recorded 63 suppression(s) from default C:\Program Files (x86)\Dr. Memory/bin/suppress-default.txt

DUPLICATE ERROR COUNTS:

SUPPRESSIONS USED:

NO ERRORS FOUND:
      0 unique,     0 total unaddressable access(es)
      0 unique,     0 total uninitialized access(es)
      0 unique,     0 total invalid heap argument(s)
      0 unique,     0 total GDI usage error(s)
      0 unique,     0 total warning(s)
      0 unique,     0 total,      0 byte(s) of leak(s)
      0 unique,     0 total,      0 byte(s) of possible leak(s)
ERRORS IGNORED:
    347 still-reachable allocation(s)
         (re-run with "-show_reachable" for details)
Details: C:\Users\Admin\AppData\Roaming/Dr. Memory/DrMemory-Linear Algebra.exe.11204.000/results.txt

Skym0sh0

Ramanujan schrieb:

Edit:
valgrind gibt es wohl nicht für Windows, deshalb hab ich auf Dr. Memory zurückgegriffen:

Visual Studio 2012 bringt einige Analysetools mit. Damit könntest du vielleciht mal durchschauen.
Ich weiss allerdings nicht welche Versionen das haben, da ich hier nur die Ultimate Edition habe.

Ramanujan

Ich habe den Fehler (bzw. vielmehr eine Löusng) gefunden. Verstehen tue ich ihn nicht.

Wenn ich die Eingabematrix so erzeuge, läuft alles durch (es gilt d=2):

Matrix<double> a_pre(n, d);
for (int i = 0; i < n; i++) {
	if (i < n/2) {
		a_pre(i, 0) = useful::random_LL(-50, 49) + useful::rnd(0,1);
		a_pre(i, 1) = useful::rnd(-.001, .001);
	} else {
		a_pre(i, 1) = useful::random_LL(-50, 49) + useful::rnd(0,1);
		a_pre(i, 0) = useful::rnd(-.001, .001);
	}
}

Wenn ich sie so erzeuge, stimmen die Ergebnisse nicht mehr:

Matrix<double> a_pre(n, d);
for (int i = 0; i < n; i++) {
	if (i < n/2) {
		a_pre(i, 0) = useful::rnd(-50, 50);
		a_pre(i, 1) = useful::rnd(-.001, .001);
	} else {
		a_pre(i, 1) = useful::rnd(-50, 50);
		a_pre(i, 0) = useful::rnd(-.001, .001);
	}
}

Hinter den useful-random-Funktionen steckt das folgende:

double useful::rnd() {
	return ((double) random_LL() / RANDOM_LL_MAX);
}

double useful::rnd(double min, double max) {
	return rnd() * (max - min) + min;
}

const unsigned long long useful::RANDOM_LL_MAX = -1ULL;
unsigned long long useful::random_LL() {
	long long rnd = 0;
	// 64bit-Zufallszahl erzeugen
	const unsigned lng = 8 * sizeof(long long);
	for (unsigned j = 0; j < lng; j += 15) {
		long long t = rand() & 32767; // Mindestgröße von RAND_MAX = 32767
		rnd |= t << j;
	}
	return rnd;
}

long long useful::random_LL(long long min, long long max) {
	unsigned long long size = max - min;
	unsigned long long rnd;
	do {
		rnd = random_LL();
	} while (rnd > RANDOM_LL_MAX - RANDOM_LL_MAX % size);
	return min + rnd % size;
}

In beiden Fällen liegen alle Werte der Matrix zwischen -50 und 50 und es gibt kein NaN-Werte.

Eigentlich sehe ich nur noch eine Möglichkeit für das Verhalten:
Bei dem Programm geht es darum, eine sogenannte Einbettung zu berechnen. Dazu wird eine implizit gespeicherte (s x n)-Matrix S erstellt, die als Einträge -1 und 1 hat. Die Matrix S*A/sqrt(s) hat nun Ähnliche Eigenschaften wie A. Unter anderem gilt, dass (A^T)*A (hier steht A^T für die Transponierte von A) etwa gleich ((SA)^T)*SA ist. Wählt man z.B. s = 34 sind so Abweichungen von höchstens 20% zu erwarten. Für die Matrix S hab ich einen eigenen pseudo-Zufallsgenerator geschrieben (der vor allem das Element Sij in O(1) Zeit berechnen kann). Jetzt könnte es natürlich trotzdem sein, dass irgendwelche stochastischen Abhängigkeiten zwischen der rand()-basierten rnd()-Funktion und der Matrix S bestehen, die dazu führen, dass das Ergebnis nicht mehr korrekt ist. Das ist extremst unwahrscheinlich, aber ich kann mir das anders nicht erklären. Zumindest würde das auch erklären, warum es bei Linux geht: Andere rand()-Implementierung => keine Abhängigkeit mehr => alles ok.

Oder sieht jemand noch einen anderen möglichen Grund (der möglicherweise weniger an den Haaren herbei gezogen ist)?

volkard

Kannste mit wenig Aufwand nen anderen Zufallszahlengenerator unter Windows testen? Haste boost oder sowas rumfliegen?

SeppJ

Zwei Bemerkungen:
1. Dir ist aber schon klar, dass die Summe zweier Gleichverteilungen keine Gleichverteilung mehr ist, oder? Das heißt deine Codes 1 und 2 erzeugen unterschiedlich verteilte Zahlen.
2. Bezüglich deines Verdachts: Darum nimmt man rand() auch nicht für kritische Anwendungen. Nimm Generatoren von geprüfter Qualität. In Numerikbibliotheken, Boost und der C++11-Standardbibliothek wirst du fündig. Zur allergrößten Not implementier selber einen bekannten Generator - ist nicht schwierig. Dann entfällt auch der ganze Aufwand, den du derzeit treibst, der Nachbehandlung der rand()-Zahlen. Und ja: Typische rand()-Implementierungen haben Korrelationen drin.

Ramanujan

@1.) Die Funktion useful::random_LL(-50, 49) gibt mir ganze Zahlen von -50 bis 49. useful::rnd(0,1) addiert Werte zwischen 0 und 1 hinzu. Somit ergibt sich eine Gleichverteilung von -50 bis 50. Ok, wenn man es genau nimmt, sind die ganzen Werte leicht häufiger, da meine Implementierung von rnd() auch 0 und 1 als Werte zulässt. Aber das fällt nicht ins Gewicht.

@2.) Ja, ok, ich werd mal gucken, was C++11 so bietet.

Edit:
Wenn ich sowas hier benutze:

Matrix<double> a_pre(n, d); 
for (int i = 0; i < n; i++) { 
     if (i < n/2) { 
         a_pre(i, 0) = useful::rnd(-50, 50); 
         for (int k = 0; k < rand() % 100; k++)
             rand();
         a_pre(i, 1) = useful::rnd(-.001, .001); 
     } else { 
         a_pre(i, 1) = useful::rnd(-50, 50); 
         for (int k = 0; k < rand() % 100; k++)
             rand();
         a_pre(i, 0) = useful::rnd(-.001, .001); 
     } 
}

Also quasi eine zufällige Anzahl von rand()-Aufrufen dazwischen schiebe. Läuft wieder alles korrekt. Naja, da lag es wohl wirklich an der Abhängigkeit. Möglicherweise hatten die Vorzeichen von S und A ja vorher ein irgendwie gemeinsames Muster.

SeppJ

Ramanujan schrieb:

@1.) Die Funktion useful::random_LL(-50, 49) gibt mir ganze Zahlen von -50 bis 49. useful::rnd(0,1) addiert Werte zwischen 0 und 1 hinzu. Somit ergibt sich eine Gleichverteilung von -50 bis 50.

Nein. Wie schon gesagt, ist die Summe von zwei Gleichverteilungen keine Gleichverteilung. Das ist ganz elementar zu wissen, wenn du mit Verteilungen arbeitest, wie die Summe zweier Zufallsvariablen verteilt ist^*. Schnapp dir unbedingt mal ein Buch über Stochastik! Das musst du wissen!

Beweis durch dein eigenes Gegenbeispiel: Die Zahl -50 hat die Wahrscheinlichkeit von random_LL==-50, also 1/100, mal die Wahrscheinlichkeit von rnd==0, also 1/2, also insgesamt 1/200. Gleiches gilt für die Zahl 50. Es sollten aber 1/101 sein, wenn es eine Gleichverteilung wäre. Ebenso ist die Wahrscheinlichkeit für eine Zahl N zwischen -50 und 50 nicht mehr 1/101, sondern die Wahrscheinlichkeit
P(N) = P(random_LL==N-1) * P(rnd==1) + P(random_LL==N)*P(rnd==0) = 1/200 + 1/200 = 1/100 != 1/101

qed.

^*: Ganz interessant und wohl eine der wichtigsten Erkenntnisse der Mathematik (riesige Teile der Naturwissenschaft beruhen darauf) ist übrigens, was die Summe unendlich vieler Zufallszahlen ist, (fast) egal wie die konkrete Verteilung dieser Zahlen aussieht. Das ist nämlich im Limes eine Gaussverteilung.

Ramanujan

rnd() liefert reelle Werte zwischen 0 und 1 (jeweils einschließlich). Die Wahrscheinlichkeit für den Wert für 0 ist 1/2^64 (wenn sizeof(long long) == 8 gilt). Da allerdings rnd() auch 1 als Ergebnis liefern kann, sind ganze Werte um 1/2^64 wahrscheinlicher (also theoretisch keine Gleichverteilung). Wenn rnd() von 0 (einschließlich) bis 1 (ausschließlich) liefern würde, wäre aber alles ok. Da double ziemlich ungenau ist (im Vergleich zu long long), spielt das aber keine Rolle, man kann also praktisch die Verteilung von einer Gleichverteilung nicht unterscheiden.

Im Allgemeinen hast du natürlich Recht: Wenn X und Y gleichverteilt, ist X+Y nicht zwingend ebenfalls gleichverteilt.

SeppJ

Ramanujan schrieb:

rnd() liefert reelle Werte zwischen 0 und 1 (jeweils einschließlich). Die Wahrscheinlichkeit für den Wert für 0 ist 1/2^64 (wenn sizeof(long long) == 8 gilt). Da allerdings rnd() auch 1 als Ergebnis liefern kann, sind ganze Werte um 1/2^64 wahrscheinlicher. Wenn rnd() von 0 (einschließlich) bis 1 (ausschließlich) liefern würde, wäre aber alles ok. Da double ziemlich ungenau ist (im Vergleich zu long long), spielt das aber keine Rolle.

Deine Argumentation ist nicht überzeugend, der Datentyp spielt keine Rolle. Aber wenn du nicht zur Kenntnis nehmen möchtest, wenn man dich auf Fehler hinweist, dann will ich dich auch nicht zu deinem Glück zwingen.

hustbaer

[Unsinn]
~~Freund Ramanujan...
Wenn du mit zwei Würfeln würfelst, und das Ergebnis zusammenzählst, wie ist dann die Verteilung? Ist das dann gleichverteilt?

Nein, ist es nicht. Genau so wenig wie die Werte in deinem Beispiel gleichverteilt sind.~~[/Unsinn]

Davon abgesehen ist die Implementierung von random_LL furchtbar. Sie verlässt sich auf rand() , und da rand() meist ein ganz einfacher LCG ist, hagelt es nur so kurze Serien in den unteren Bits.
Die Werte die random_LL produziert sind also voll von diversen statistischen Auffälligkeiten. Das selbe gilt dann natürlich auch für rnd() .

Und wie diese "statistischen Auffälligkeiten" aussehen ändert sich natürlich auch, wenn du 1x mehr oder weniger random_LL aufrufst.