www.virtual-maxim.de

Wer sich auf meinem Blog umschaut, sieht dass viele meiner Artikel mathematische Formeln beinhalten. Bis jetzt verwende ich ein WordPress-Plugin, welches Formeln in LaTeX-Format verarbeiten kann. Damit werden aus dem LaTex-Code PNG-Grafiken erstellt und in den HTML-Code eingebunden. Ein zweiseitiger Artikel hat deswegen locker 30 eingebundene Grafiken.
Der Leser muss dadurch nicht nur lange warten bis die Grafiken geladen sind, nein, er bekommt auch noch ein inkonsistentes Erscheinungsbild geboten.

Bereits 1998 hat W3C einen Standard zur Darstellung von mathematischen Formeln in HTML veröffentlicht – MathML. Dabei werden Formeln in XML-Format einfach in den HTML-Code eingefügt. Der Code selbst kann zum Beispiel von einem speziellen Programm generiert werden.

Beispiel MathML-Code:

<math mode="display" xmlns="http://www.w3.org/1998/Math/MathML">
    <mover><mi>x</mi><mo>&rarr;</mo></mover>
    <mo>=</mo>
    <mfenced>
        <mtable><mtr>4<mtd/></mtr><mtr>2<mtd/></mtr><mtr>5<mtd/></mtr></mtable>
    </mfenced>
    <mo>+</mo><mo>&mu;</mo>
    <mfenced>
        <mtable><mtr>-3<mtd/></mtr><mtr>3<mtd/></mtr><mtr>4<mtd/></mtr></mtable>
    </mfenced>
</math>

Beispiel MathML-Code 2:

<math mode="display" xmlns="http://www.w3.org/1998/Math/MathML">
  <mrow>
    <mi>x</mi>
    <mo>=</mo>
    <mfrac>
      <mrow>
        <mo form="prefix">&#x2212;<!-- &minus; --></mo>
        <mi>b</mi>
        <mo>&#x00B1;<!-- &PlusMinus; --></mo>
        <msqrt>
          <msup>
            <mi>b</mi>
            <mn>2</mn>
          </msup>
          <mo>&#x2212;<!-- &minus; --></mo>
          <mn>4</mn>
          <mo>&#x2062;<!-- &InvisibleTimes; --></mo>
          <mi>a</mi>
          <mo>&#x2062;<!-- &InvisibleTimes; --></mo>
          <mi>c</mi>
        </msqrt>
      </mrow>
      <mrow>
        <mn>2</mn>
        <mo>&#x2062;<!-- &InvisibleTimes; --></mo>
        <mi>a</mi>
      </mrow>
    </mfrac>
  </mrow>
</math>

Und was sieht ihr? Eine schön gezeichnete Formeln, die sich der Schriftart und Schriftgröße automatisch anpassen (mit Strg und “+” ausprobieren) oder komisches Zahlenwirrwarr?
So sehen die beiden Beispiele in Firefox aus:

Benutzer von Google Chrome, Internet Explorer und Safari sehen wahrscheinlich nur sowas in der Art:

Praktisch alle Browserhersteller kämpfen um Marktanteile und versuchen immer weiter die Seitendarstellung oder die Ladezeit um paar zerquetschte Millisekunden schneller zu machen, aber keiner bemüht sich Standards, die bereits über 10 Jahre alt sind, zu implementieren. Firefox und Opera (nicht getestet) haben zumindest die Basis implementiert, die für viele Bereiche ausreichen. WebKit basierte Browser wie Chrome oder Safari spielen zwar in den Nightly Builds mit MathML rum, aber scheinen es auch nicht gerade eilig zu haben.

Etwas kurios ist die Sache mit Internet Explorer. Um MathML in IE darzustellen braucht man ein Plugin, dass es bereits seit IE 5.5 gibt (soweit ich weiß).
Das kuriose ist, das MS bereits super funktionierenden Parser und Engine zur Darstellung von MathML hat. Das aktuelle Office 2010 hat einen wunderbaren Formel-Editor drin, der die Formeln in MathML-Format einbinden. Das “Mathematik-Eingabebereich” Tool, dass bei Windows 7 standardmäßig dabei ist, mit dem man Formeln mit der Maus bzw. Stift schreiben kann um diese in Word einzufügen, erzeugt auch MathML-Code.
Die Darstellung ist auch super, wie man es zum Beispiel an diesem Artikel sehen kann.
Sie haben die Technologie und nutzen sie nicht. Wie es scheint, hat Microsoft ein großes Kommunikationsproblem zwischen den einzelnen Entwicklungsabteilungen.
Mit dem Internet Explorer 9 wird sich dabei auch nichts ändern. Zumindest erhält man mit IE9 Preview 3 den gleichen Mist wie mit IE8.

Meiner Meinung nach ist MathML Unterstützung eine Marktlücke, die nicht erkannt wird. Die Größe des wissenschaftlichen Sektors wird unterschätzt.
Viele Universitäten stellen Fachtexte online, bei denen Formeln als Grafiken eingebunden sind. Ersten sieht es total schlecht bis unleserlich aus, zweiten habe ich schon so viele Fachtexte gesehen, in denen die Grafiken einfach fehlen (kA warum), was dazu führt dass der ganze Text unbrauchbar wird. Wären die Formeln als MathML eingebunden, so könnte man sie auch heute noch betrachten.

Dies führt zum nächsten Punkt: Informationserhalt im Internet. Wollte man früher einer breiter Menschenmasse was sagen, so schrieb man ein Buch, heute betreibt man eine Internetpräsenz. Geht die Website offline, so verschwinden auch viele Informationen.
Es gibt Projekte, die versuchen diese Informationen zu erhalten. Auch Google sichert alte Inhalte. Das Problem ist, dass meistens nur der HTML-Code gesichert wird. Damit haben wir wieder das Problem, dass die Fachtexte ohne Formeln unbrauchbar sind.

Man schimpft ständig auf MS, dass IE Standards nicht überstürzt, aber sind denn andere Browserhersteller so viel besser?
Es bleibt nur zu hoffen, dass die Browserhersteller mehr an ihren Produkten arbeiten, als nur jeden Monat die Versionsnummer zu erhöhen, weil die Seitenladenzeit von 500ms auf 499ms verbessert werden konnte.

Speccy ist ein schlankes Programm für Windows 7, Vista, XP und 2000, welches die wichtigsten Informationen zur verbauter PC-Hardware anzeigt.
Dies kann sehr hilfreich sein, wenn man einen neuen Treiber sucht, aber die genaue Bezeichnung der Hardware nicht kennt.

Die aktuelle Version, die heute erschienen ist und eine Reihe von Fehlern korrigiert, hat auch eine native 64-Bit Unterstützung bekommen.

Irgendwie mag ich diese schlanken Programme aus dem Hause Piriform =)

Motivation

In der Spieleprogrammierung werden häufig trigonometrische Funktionen, insbesondere Kosinus und Sinus verwendet. Vor allem im 3D Bereich bauen viele Berechnungen auf trigonometrischen Funktionen auf. Demnach könnten diese Funktionen nicht schnell genug sein. Um mehr Leistung aus dem Programm herauszuholen greifen viele auf Lookup Tabellen zurück. Das ist keine schlechte Lösung, wenn das Ergebnis nicht sehr genau sein soll. Möchte man jedoch mindestens eine Genauigkeit von einer Nachkommastelle haben, so braucht man einen 3600 Werte großen float-Array und je größer das Array, desto kleiner wird der Geschwindigkeitsgewinn gegenüber der Standardimplementierung sein.

Es gibt einen anderen Weg: die Approximation(deutsch: Annäherung) einer Funktion. Durch diese mathematische Methode kann man eine Funktion gewinnen, die eine bestimmte Funktion in einem vordefinierten Intervall annähernd beschreibt.
In unserem Fall würden wir eine Funktion gewinnen, die den Sinus im Intervall von –PI bis + PI gut beschreibt. Ein größeres Intervall wird nicht gebraucht, da der Sinus 2Pi-periodisch ist.

Sinus

Wie genau man auf die Näherungsfunktion kommt, werde ich hier nicht behandeln.
Wer sich dafür interessiert, sollte nach Begriffen Modellierung, Approximation, Potenzreihenentwicklung, Taylor’sche Formel oder Maclaurin’sche Formel suchen. Das Thema Modellierung wird normalerweise auf einem Gymnasium in der zwölften Klasse behandelt.
Stattdessen nehmen wir gleich die fertige Potenzreihe und setzen sie in C++ um.

Diese Funktion modelliert Sinus um den Nullpunkt herum. Je mehr Potenzreihenglieder aufsummiert werden, desto besser beschreibt diese Funktion den Sinus. Dies kann man sehr schön an folgenden Funktionsplots sehen.

Legende:
Orange = sin(x)
Schwarz =

Wie man sieht, wird der Sinus im Koordinatenursprung gut beschrieben, nicht aber an den PI-Grenzen. Wir addieren einen weiteren Term dazu.

Legende:
Orange = sin(x)
Schwarz =
Schon besser, aber die Ungenauigkeit doch noch zu groß. Wir erweitern weiter.

Legende:
Orange = sin(x)
Schwarz =
Diese Näherungsfunktion beschreibt den Sinus im Intervall von –PI bis PI schon ziemlich genau, was man daran erkennt, dass in diesem Intervall der orangenfarbene Graph komplett überdeckt wird.
Hier noch mal ein größerer Ausschnitt.

Wie gut diese Modellierung wirklich ist, schauen wir uns später an konkreten Zahlen an.

Demnach lautet unsere Funktion . Bevor wir sie aber in C++ umsetzen, müssen wir bedenken, dass sie nur im Intervall von –PI bis PI gilt. Wird ein Argument übergeben, der nicht in diesem Intervall liegt, so liefert die Funktion falsche Werte, was man an dem letzten Graphen sehr gut zu erkennen ist.

Die Sinusfunktion ist genauso wie der Kosinus periodisch. Das heißt, wenn man zum Funktionsargument z-mal dazu addiert oder von dem Funktionsargument z-mal abzieht, wird sich der Wert der Funktion sin(x) nicht verändert (z ist eine Ganzzahl).

Man könnte auch schreiben:

oder

Für uns bedeutet das, dass ein Funktionsargument, der nicht in dem Intervall von –PI bis PI liegt, auf dieses Intervall zurückgerechnet werden kann ohne dass das Ergebnis der Funktion verfälscht wird.

Dazu einige Beispiele (Argument => Umrechnung => neuer Wert ):

85,222 => 85,222-14*2*PI => ~ -2,743
12,982 => 12,982-2*PI => ~ 0,416
-29, 314 => -29, 314+5*2*PI => ~ 2.102

Wie kommt man auf die z für die Umrechnung?
Es wird zwischen zwei Fällen unterschieden:

1. x ist größer als PI
• int T = x + PI
• T wird durch 2*PI geteilt, das Ergebnis ohne Nachkommastelle(!) ist z
2. x ist kleiner als PI
• int T = x – PI
• T wird durch 2*PI geteilt, das Ergebnis ohne Nachkommastelle(!) ist z

Als nächstes wird das Argument durch die Formel auf das Intervall –PI bis PI zurückgerechnet. Ob das z positiv oder negativ ist, spielt keine Rolle, den wenn z negativ ist wird der Term zu .

Nachdem das Argument im richtigen Intervall liegt, können wir unsere Taylor-Reihe benutzen. Doch zuerst schreiben wir die Taylor-Reihe so um, dass die Fakultät- und Division-Operationen verschwinden, da diese sehr langsamer sind.

Wir bekommen

und durch Umformung

An dieser Stelle kommt unsere erste Sinus-Funktion in C++.

const float PI = 3.14159265f;
const float PI2 = PI*2;

// schnelle Sinus Implementierung
inline float _fastcall fastSin(float x)
{
	if(PI < x)
	{
		x = x-static_cast<int>((x+PI)/(PI2))*PI2;
	}
	else if(x < -PI)
	{
		x = x-static_cast<int>((x-PI)/(PI2))*PI2;
	}

	return x*(1 - x*x*(0.16666667f - x*x*(0.00833333f - x*x*(0.0001984f - x*x*0.0000027f))));
}

Ich habe einen kleinen Geschwindigkeitstest durchgeführt, bei dem std::sin(x) und unsere Implementierung 10-Millionen Mal aufgerufen wurde. Das Ergebnis ist ein Mittelwert aus 5 Durchläufen: std::sin(): 1174ms, fastSin(): 656ms. Selbstverständlich können diese Werte nur qualitativ betrachtet werden.

Unsere Implementierung weist eine Genauigkeit von 2-6 Nachkommastellen auf, die für Spieleprogrammierung akzeptabel sein sollte. Wie man erkennt, ist die Ungenauigkeit an den -Grenzen am höchsten.

Um eine höhere Genauigkeit zu erzielen, muss man das Polynom um weitere Glieder erweitern: .

Als nächsten Schritt werden wir die Methode in Inline Assembler schreiben, wobei dabei Intels SSE Technik verwendet wird. SSE wird ab Pentium 3 und Athon XP unterstützt. Ich werde hier nicht erklären, wie SSE funktioniert, weil es den Rahmen dieses Artikels sprengen würde. Es handelt sich um 1:1 Konvertierung der C++ Version der Funktion.

// schnelle Sinus Implementierung mit Inline Assembler
inline float _fastcall fastSinSSE(float x)
{
	static const float _0_16666 = -1/6.0f;
	static const float _0_00833 = 1/120.0f;
	static const float _0_00019 = -1/5040.0f;
	static const float _0_0000027 = 1/362880.0f;

	_asm
	{
		movss	xmm0, x
		comiss	xmm0, PI	// if (x > PI)
		jbe	SHORT ln1

		addss	xmm0, PI
		jmp SHORT ln2

ln1:
		movss	xmm1, _PI
		comiss	xmm1, xmm0	// if ( x < -PI)
		jbe	SHORT ln3

		subss	xmm0, PI

ln2:
		divss	xmm0, PI2

		cvttss2si ecx, xmm0	// = z ohne rundung!
		cvtsi2ss xmm1, ecx	// z wieder in float umwandeln

		movss	xmm0, x
		mulss	xmm1, PI2
		subss	xmm0, xmm1
ln3:
		// xmm0 beinhaltet immer das aktuelle ergebnis
		// xmm1 beinhaltet zu jeder zeit x^2
		// xmm2 beinhaltet zu jeder zeit x
		// x + x*x*x*(_0_16666 + x*x*(_0_00833 + x*x*(_0_00019f + _0_0000027*x*x))));

		movss	xmm2, xmm0			// x in xmm0 legen
		movss	xmm1, xmm0		// x in xmm1 legen
		mulss	xmm1, xmm1		// xmm1 = x*x
		movss	xmm0, _0_0000027
		mulss	xmm0, xmm1		// xmm0 = _0_0000027*x*x
		addss	xmm0, _0_00019
		mulss	xmm0, xmm1		// xmm0 = x*x*(_0_00019f + _0_0000027*x*x)
		addss	xmm0, _0_00833
		mulss	xmm0, xmm1		// xmm0 = x*x*(_0_00833 + x*x*(_0_00019f + _0_0000027*x*x))
		addss	xmm0, _0_16666
		mulss	xmm0, xmm1		// xmm0 = x*x*(_0_16666 + x*x*(_0_00833 + x*x*(_0_00019 + _0_0000027*x*x)))
		mulss	xmm0, xmm2
		addss	xmm0, xmm2
		movss	x, xmm0			// das ergebnis in x rüberschieben
	}

	return x;
}

Durch Assembler wurde unsere Funktion noch ein wenig schneller: 10-Millionen Berechnungen in 493ms. Also rund zwei Mal schneller als die Standardimplementierung. Der Preis dafür ist die Ungenauigkeit.

Als nächstes schauen wir uns im Schnelldurchlauf den Kosinus an.

Kosinus

Das Schreiben der Kosinus-Funktion funktioniert analog zur Sinus-Funktion. Man kann aber auch die Gesetzmäßigkeit ausnutzen. Man schreibt einfach eine neue Funktion, die unsere Sinus-Funktion mit dem Argument aufruft:

float fastCos(float x)
{
	return fastSin(PI/2 – x);
}

Ja, so schnell kann es gehen. Der Nachteil von dieser Funktion ist, dass das Aufrufen der Funktion länger dauert als die Berechnung selbst, also machen wir es analog zum Sinus – Nährungsfunktion in Inline Assembler implementieren.

Die Taylor-Reihe für Kosinus:

Wir nehmen die ersten sechs Glieder.

Legende:
Orange = cos(x)
Schwarz =
Schon besser, aber die Ungenauigkeit doch noch zu groß. Wir erweitern weiter.

Wir schreiben das Polynom in leistungsstärkere Variante um:

// schnelle Kosinus Implementierung
inline float _fastcall fastCos(float x)
{
	if(PI < x)
	{
		x = x-static_cast<int>((x+PI)/(PI2))*PI2;
	}
	else if(x < -PI)
	{
		x = x-static_cast<int>((x-PI)/(PI2))*PI2;
	}

	return 1-x*x*(0.5f-x*x*(0.04166667f-x*x*(0.00138889f-x*x*(0.00002480f-x*x*0.000000275f))));
}

Die Genauigkeit von Kosinus ist besser als die von Sinus und liegt bei 3-6 Nachkommastellen. Das liegt daran, dass wir mehr Reihenglieder benutzt haben.

Beim Geschwindigkeitstest (10 Millionen Mal aufgerufen) braucht die Standardmplementierung von Kosinus 2305ms und fastCos 1337ms.

Die Inline Assembler Version schafft 10 Millionen Aufrufe in 996ms.

// schnelle Kosinus Implementierung mit Inline Assembler
inline float _fastcall fastCosSSE(float x)
{
	static const float _0_5 = -1/2.0f;
	static const float _0_0416 = 1/24.0f;
	static const float _0_001387 = -1/720.0f;
	static const float _0_0000248 = 1/40320.0f;
	static const float _0_000000275 = -1/3629000.0f;
	static const float _1 = 1.0f;

	// Wenn x groeßer oder kleiner als PI,
	// dann wird x auf das Intervall -PI bis PI zurückgerechnet
	_asm
	{
		movss	xmm0, x
		comiss	xmm0, PI	// if (x > PI)
		jbe	SHORT ln1

		addss	xmm0, PI
		jmp SHORT ln2

ln1:
		movss	xmm1, _PI
		comiss	xmm1, xmm0	// if ( x < -PI)
		jbe	SHORT ln3

		subss	xmm0, PI

ln2:
		divss	xmm0, PI2

		cvttss2si ecx, xmm0	// = z ohne rundung!
		cvtsi2ss xmm1, ecx	// z wieder in float umwandeln

		movss	xmm0, x
		mulss	xmm1, PI2
		subss	xmm0, xmm1
ln3:
		// xmm0 beinhaltet immer das aktuelle ergebnis
		// xmm1 beinhaltet zu jeder zeit x^2

		//  1-x*x*(0.5f-x*x*(0.04166667f-x*x*(0.00138889f-x*x*(0.00002480f-x*x*0.000000275f))));

		movss	xmm1, xmm0		// x in xmm1 legen
		mulss	xmm1, xmm0		// xmm1 = x*x
		movss	xmm0, _0_000000275
		mulss	xmm0, xmm1		// xmm0 = _0_000000275*x*x
		addss	xmm0, _0_0000248
		mulss	xmm0, xmm1		// xmm0 = x*x*(_0_0000248 + _0_000000275*x*x)
		addss	xmm0, _0_001387
		mulss	xmm0, xmm1		// xmm0 = x*x*(_0_001387 + x*x*(_0_0000248 + _0_000000275*x*x))
		addss	xmm0, _0_0416
		mulss	xmm0, xmm1		// xmm0 = x*x*(_0_0416 + x*x*(_0_001387 + x*x*(_0_0000248 + _0_000000275*x*x)))
		addss	xmm0, _0_5
		mulss	xmm0, xmm1		// xmm0 = x*x*(_0_5+x*x*(_0_0416 + x*x*(_0_001387 + x*x*(_0_0000248 + _0_000000275*x*x))))
		addss	xmm0, _1		// xmm0++
		movss	x, xmm0			// das ergebnis in x rüberschieben
	}

	return x;
}

Auch hier gilt, dass weitere Reihenglieder die Genauigkeit erhöhen. Aber auch schon jetzt setzt der Datentyp float eine Genauigkeitsgrenze fest, da dieser nur 7-8 Nachkommastellen und sogar nur ungenau speichern kann.

Als Abschluss kann man hier die Funktionen als C++ Code mit den dazugehörigen Tests downloaden.
Schnelle trigonometrische Funktionen 2.06 kB -