Sinus-Eingabe skalieren
Markus Wichmann
ich versuche hier gerade, den Sinus anzunähern (wird mal ein Algorithmus
für die dietlibc). OK, bisher hab ich schon interne Abhängigkeiten zu
floor(), ceil() und fabs(), was ja eigentlich zu vermeiden ist, aber
naja.
So, der Algorithmus, der den Sinus annähert funktioniert nur mit einem
Eingabewert im Bereich von -Pi bis Pi. Alle anderen Eingaben müssen
entsprechend skaliert werden. Nun hab ich mir folgendes gedacht:
Sei x > 0, dann ist
sin x = sin (x - 2*n*Pi)
mit n = max {n \in N | 2*n*Pi < x}
= max {n \in N| n < x/(2*Pi)}
Das skaliert ein positives x auf den Bereich zwischen 0 und 2*pi. Nun
gut, sollte x immer noch größer als Pi sein, ziehe ich wieder 2*pi ab
und das Problem ist gelöst. So hab ich das ganze auch implementiert:
if (x > 2*M_PI)
x -= 2*M_PI*floor(x / 2*M_PI);
if (x > M_PI)
x -= 2*M_PI;
Die Definition von M_PI hab ich aus der math.h kopiert. Außerdem linke
ich das ganze (momentan) gegen die (g)libm.
Ach ja, für negative x dachte ich mir:
sin x = sin (x + 2*n*Pi)
mit n = min {n \in N | 2*n*Pi > -x}
= min {n \in N | n > -x/(2*Pi)}
Nur sagt mir mein Testprogramm irgendwie, dass das nicht klappt:
sin(6.00) = -0.278463
sin(6.25) = -0.032793
sin(6.50) = 332740.314793
Hier mal der komplette Algorithmus:
#define M_PI 3.14159265358979323846 /* pi */
double ceil(double);
double floor(double);
double fabs(double);
static const double A = -4. / (M_PI * M_PI);
static const double B = 4. / M_PI;
static const double P = 0.775;
static const double Q = 0.225;
double sin(double x)
{
if (x < -2*M_PI)
x += 2*M_PI*ceil(-x / 2*M_PI);
if (x > 2*M_PI)
x -= 2*M_PI*floor(x / 2*M_PI);
if (x > M_PI)
x -= 2*M_PI;
double y = A * x * fabs(x) + B * x;
return P * y + Q * y * fabs(y);
}
#ifdef WITH_MAIN
#include <stdio.h>
int main()
{
for (double d = -2; d < 3; d += 0.0625)
printf("sin(%f * Pi) = %f\n", d, sin(d*M_PI));
return 0;
}
#endif
Tschö,
Markus
P.S.: Mein Compiler scheint buggy zu sein: Wenn ich die ersten zwei
Zeilen des Sinus durch
if (x < -M_PI)
return -sin(-x);
ersetze, kommt eine Segmentation fault. Ich habe das nachgeforscht: Es
wird an dieser Stelle immer wieder sin(x) aufgerufen, nicht sin(-x).
WTF?
Ach ja: Angaben zu meinem System:
$ uname -a
Linux debbox 2.6.30-2-amd64 #1 SMP Fri Sep 25 22:16:56 UTC 2009 x86_64
GNU/Linux
$ gcc --version |& head -1
gcc (Debian 4.3.4-6) 4.3.4
--
Nur weil ein Genie nix reißt, muß ja nun nicht gleich jeder Idiot
pausieren... Bully hats ja auch geschafft.
-- gUnter nanonüm in de.alt.anime
Thomas Richter
> Hi all,
> ich versuche hier gerade, den Sinus anzunähern (wird mal ein Algorithmus
> für die dietlibc). OK, bisher hab ich schon interne Abhängigkeiten zu
> floor(), ceil() und fabs(), was ja eigentlich zu vermeiden ist, aber
> naja.
Erstmal einige Punkte vorweg: Eine numerisch stabile Implementierung von
transzendenten Funktionen gehört nicht zu den Dingen, die man im
"Hau-Ruck" Verfahren implementiert. Falls das Resultat wirklich generell
nutzbar sein soll, also für allgemeine Anwendungen taugen soll, so
überlasse man das besser den Experten der Numerik.
Einige wünscheswerte Eigenschaften wie Symmetrie oder Stetigkeit sind
mit Sicherheit zu berücksichtigen, aber darüberhinaus auch noch das
richtiges Grenzwertverhalten etc...
Zum Problem:
> double sin(double x)
> {
> if (x < -2*M_PI)
> x += 2*M_PI*ceil(-x / 2*M_PI);
Hiermit bleibt das x in (-2Pi,0] immer noch negativ. Soll das so sein?
Ich würde stattdessen vorschlagen, zunächst negative Eingaben abzufangen
und -sin(-x) aufzurufen, damit die Symmetrie exakt stimmt.
Ferner vermute ich, dass die obige Berechnung des Modulus nicht recht
stabil ist, da Du einmal hin, und dann wieder zurückskalierst, so dass
Rundungsfehler akkumuliert werden. Ich würde stattdessen zunächst das
Argument so skalieren, dass M_PI auf 1 fällt, und würde dann den Rest
abseparieren. Dazu reicht C-Integerarithmetik, und floor() und ceil()
entfallen. Das Approximationspolynom sollte dann natürlich entsprechend
angepasst werden, besser noch mit einer Fallunterscheidung für sehr
kleine Werte (dort genauer appoximieren) und für Werte in der Nähe von 1
(ebenso genauer), der Rest durch ein Minimax-Polynom.
> if (x > 2*M_PI)
> x -= 2*M_PI*floor(x / 2*M_PI);
Ditto. Statt floor zu benutzen, kannst Du auch (für positive x) nach int
casten, dann ist die Abhängigkeit weg.
> if (x > M_PI)
> x -= 2*M_PI;
Hiermit reduzierst Du nur das Intervall (Pi,2 Pi] auf (-Pi,0], aber
inkorrekt, der Sinus ist dort ja spiegelsymmetrisch. D.h. für x > Pi
müsstest Du stattdessen sin(2Pi - x) ausrechnen, und nicht das obige.
Der Fall für negatives x, also x < -Pi fehlt.
> double y = A * x * fabs(x) + B * x;
> return P * y + Q * y * fabs(y);
Mit der vorgeschlagenen Reduktion sin(x) = -sin(-x) würde fabs hier
entfallen. Wie gesagt, ich würde hier drei Fälle unterscheiden und die
Polynome sehr sorgfältig wählen.
> P.S.: Mein Compiler scheint buggy zu sein: Wenn ich die ersten zwei
> Zeilen des Sinus durch
>
> if (x < -M_PI)
> return -sin(-x);
>
> ersetze, kommt eine Segmentation fault. Ich habe das nachgeforscht: Es
> wird an dieser Stelle immer wieder sin(x) aufgerufen, nicht sin(-x).
> WTF?
Kann ich nicht nachvollziehen. Abgesehen davon sollte es wohl besser
heißen: if (x < 0) (!)
Grüße,
Thomas
Markus Wichmann
> Markus Wichmann wrote:
>> Hi all,
>> ich versuche hier gerade, den Sinus anzunähern (wird mal ein Algorithmus
>> für die dietlibc). OK, bisher hab ich schon interne Abhängigkeiten zu
>> floor(), ceil() und fabs(), was ja eigentlich zu vermeiden ist, aber
>> naja.
>
> Erstmal einige Punkte vorweg: Eine numerisch stabile Implementierung von
> transzendenten Funktionen gehört nicht zu den Dingen, die man im
> "Hau-Ruck" Verfahren implementiert. Falls das Resultat wirklich generell
> nutzbar sein soll, also für allgemeine Anwendungen taugen soll, so
> überlasse man das besser den Experten der Numerik.
>
> Einige wünscheswerte Eigenschaften wie Symmetrie oder Stetigkeit sind
> mit Sicherheit zu berücksichtigen, aber darüberhinaus auch noch das
> richtiges Grenzwertverhalten etc...
>
Ich hab den Algorithmus von einer Seite im Internet, die damit warb,
dass dieser Algorithmus (immerhin lediglich ein biquadratisches Polynom)
einem Taylor-Polynom siebten Grades gerade in puncto Stetigkeit
überlegen ist. Es ist mehr oder weniger eine Parabel, die auf bestimmte
Punkte festgelegt ist. Wie es der Zufall will, ist die Parabel selbst
fast immer zu groß, während ihr Quadrat fast immer zu niedrig ist. Also
addiere ich beides gewichtet.
> Zum Problem:
>
>> double sin(double x)
>> {
>> if (x < -2*M_PI)
>> x += 2*M_PI*ceil(-x / 2*M_PI);
>
> Hiermit bleibt das x in (-2Pi,0] immer noch negativ. Soll das so sein?
> Ich würde stattdessen vorschlagen, zunächst negative Eingaben abzufangen
> und -sin(-x) aufzurufen, damit die Symmetrie exakt stimmt.
>
Das hab ich ja versucht, aber bei diesem Statement ist mein Compiler
bockig.
Aber du hast recht, eigentlich will ich das x ja nach [-Pi; Pi]
skalieren. Also
if (x < -2*M_PI)
x += 2*M_PI*(ceil(-x / 2*M_PI) + 1)
?
> Ferner vermute ich, dass die obige Berechnung des Modulus nicht recht
> stabil ist, da Du einmal hin, und dann wieder zurückskalierst, so dass
> Rundungsfehler akkumuliert werden. Ich würde stattdessen zunächst das
> Argument so skalieren, dass M_PI auf 1 fällt, und würde dann den Rest
> abseparieren. Dazu reicht C-Integerarithmetik, und floor() und ceil()
> entfallen. Das Approximationspolynom sollte dann natürlich entsprechend
> angepasst werden, besser noch mit einer Fallunterscheidung für sehr
> kleine Werte (dort genauer appoximieren) und für Werte in der Nähe von 1
> (ebenso genauer), der Rest durch ein Minimax-Polynom.
>
Fallunterscheidungen sind fast immer Performance-Killer. Was meinst du,
warum ich so absonderliche Konstruktionen wie "x * fabs(x)" sonst
verwende? Ich könnte natürlich auch
if (x < 0)
y = -A * x * x + B * x;
else
y = A * x * x + B * x;
verwenden, aber das wollte ich eben vermeiden.
Aber deine Idee ist nicht schlecht. Vielleicht veröffentliche ich sogar
die reskalierte Sinusfunktion. Sowas hier:
/**
* sine function
* x - angle in radians
*/
double sin(double x)
{
return tsin(x / (2*M_PI));
}
/**
* sine function
* x - angle in degrees
*/
double dsin(double x)
{
return tsin(x / 360);
}
/**
* sine function
* x - angle in turns
*/
double tsin(double x)
{
x -= trunc(x); /* kann man die Nachkommastellen direkt bekommen? */
if (x > .5)
x -= 1.;
else if (x < -.5)
x += 1.;
//...
}
Und schon ist der Wert zwischen -0,5 und 0.5.
>> if (x > 2*M_PI)
>> x -= 2*M_PI*floor(x / 2*M_PI);
>
> Ditto. Statt floor zu benutzen, kannst Du auch (für positive x) nach int
> casten, dann ist die Abhängigkeit weg.
>
Stimmt, allerdings könnte es sein, dass ich mir damit in den Fuß
schieße: Wenn der Compiler das ganze nämlich strikt ausführt, könnte es
Datenverlust geben, weil in ein double (IEEE 754 double precision)
wesentlich größere Zahlen passen als in ein int (in meinem Fall 8 Byte
Signed Integer mit Zweierkomplement).
>> P.S.: Mein Compiler scheint buggy zu sein: Wenn ich die ersten zwei
>> Zeilen des Sinus durch
>>
>> if (x < -M_PI)
>> return -sin(-x);
>>
>> ersetze, kommt eine Segmentation fault. Ich habe das nachgeforscht: Es
>> wird an dieser Stelle immer wieder sin(x) aufgerufen, nicht sin(-x).
>> WTF?
>
> Kann ich nicht nachvollziehen. Abgesehen davon sollte es wohl besser
> heißen: if (x < 0) (!)
>
> Grüße,
> Thomas
Ich habe meinen Compiler nebst Version genannt, und zumindest hier ist
das Ergebnis deterministisch. Der ganze Kram funktioniert übrigens
schlagartig, wenn ich den Compiler das ganze in Assembler übersetzen
lasse, anschließend das Assembler-Listing entsprechend verändere (obwohl
ich eigentlich AT&T-Syntax-Verächter bin) und anschließend vom
Assembler-Listing an selber weiter übersetzen lasse.
Ich habe übrigens die Compiler-Magie entdeckt: Der gcc war der Ansicht,
dass sin(x) = -sin(-x), weswegen neben dem togglen des vordersten Bits
des Argumentes auch das togglen des vordersten Bits der Rückgabe
vonnötenwar.
Das Problem ist mit der Option -fno-builtin-sin zu beheben.
Tschö,
Markus
Thomas Richter
> Thomas Richter (th...@math.tu-berlin.de) schrieb:
>>> ich versuche hier gerade, den Sinus anzunähern (wird mal ein Algorithmus
>>> für die dietlibc). OK, bisher hab ich schon interne Abhängigkeiten zu
>>> floor(), ceil() und fabs(), was ja eigentlich zu vermeiden ist, aber
>>> naja.
>> Erstmal einige Punkte vorweg: Eine numerisch stabile Implementierung von
>> transzendenten Funktionen gehört nicht zu den Dingen, die man im
>> "Hau-Ruck" Verfahren implementiert. Falls das Resultat wirklich generell
>> nutzbar sein soll, also für allgemeine Anwendungen taugen soll, so
>> überlasse man das besser den Experten der Numerik.
>>
>> Einige wünscheswerte Eigenschaften wie Symmetrie oder Stetigkeit sind
>> mit Sicherheit zu berücksichtigen, aber darüberhinaus auch noch das
>> richtiges Grenzwertverhalten etc...
>>
>
> Ich hab den Algorithmus von einer Seite im Internet, die damit warb,
> dass dieser Algorithmus (immerhin lediglich ein biquadratisches Polynom)
> einem Taylor-Polynom siebten Grades gerade in puncto Stetigkeit
> überlegen ist.
*Stetig* sind beide, denn Polynome sind nunmal stetig. Nicht nur das,
sondern sogar analytisch. (-: Wichtig ist, wie die Polynome am Ende
ihres Definitionsbereiches zusammengeklebt werden, so dass es passt.
> Es ist mehr oder weniger eine Parabel, die auf bestimmte
> Punkte festgelegt ist. Wie es der Zufall will, ist die Parabel selbst
> fast immer zu groß, während ihr Quadrat fast immer zu niedrig ist. Also
> addiere ich beides gewichtet.
Und das ist warum gut? Wie ich schon schrieb, numerische Mathematik
gehört nicht zu den Dingen, die man in Heimwerkermanier mit Schraubstock
und Hammer löst. Dazu gehört schon ein gerüttelt Maß an Mathematik.
Davon abgesehen, außer als akademische Übung erschließt sich mir der
Sinn des ganzen nicht so recht. Der numerische Koprozessor, der in jeder
halbwegs modernen Desktop-CPU eingebaut ist, hat eine
Sinus-Approximation eingebaut. Gebaut von Experten, also hoffentlich
besser, als was man als Laie so zusammenschustern kann.
>> Zum Problem:
>>
>>> double sin(double x)
>>> {
>>> if (x < -2*M_PI)
>>> x += 2*M_PI*ceil(-x / 2*M_PI);
>> Hiermit bleibt das x in (-2Pi,0] immer noch negativ. Soll das so sein?
>> Ich würde stattdessen vorschlagen, zunächst negative Eingaben abzufangen
>> und -sin(-x) aufzurufen, damit die Symmetrie exakt stimmt.
>>
>
> Das hab ich ja versucht, aber bei diesem Statement ist mein Compiler
> bockig.
Siehe unten. Dein Compiler ist eigentlich schlau und weiß, was eine
Funktion namens sin() tun soll.
> Aber du hast recht, eigentlich will ich das x ja nach [-Pi; Pi]
> skalieren.
Eigentlich sogar nur auf [0,Pi/2).
> if (x < -2*M_PI)
> x += 2*M_PI*(ceil(-x / 2*M_PI) + 1)
>
> ?
Damit wird x aber auch nur auf das Interval (-2Pi,infinity) gemappt.
>> Ferner vermute ich, dass die obige Berechnung des Modulus nicht recht
>> stabil ist, da Du einmal hin, und dann wieder zurückskalierst, so dass
>> Rundungsfehler akkumuliert werden. Ich würde stattdessen zunächst das
>> Argument so skalieren, dass M_PI auf 1 fällt, und würde dann den Rest
>> abseparieren. Dazu reicht C-Integerarithmetik, und floor() und ceil()
>> entfallen. Das Approximationspolynom sollte dann natürlich entsprechend
>> angepasst werden, besser noch mit einer Fallunterscheidung für sehr
>> kleine Werte (dort genauer appoximieren) und für Werte in der Nähe von 1
>> (ebenso genauer), der Rest durch ein Minimax-Polynom.
>>
>
> Fallunterscheidungen sind fast immer Performance-Killer. Was meinst du,
> warum ich so absonderliche Konstruktionen wie "x * fabs(x)" sonst
> verwende?
Und? Die Frage ist, worauf es bei dem Spaß ankommt. Wenn Geschwindigkeit
eine Rolle spielt, würde ich zu einer Lookup-Tabelle greifen und, wenn's
hoch kommt, zwischen den Einträgen interpolieren.
> Ich könnte natürlich auch
>
> if (x < 0)
> y = -A * x * x + B * x;
> else
> y = A * x * x + B * x;
>
> verwenden, aber das wollte ich eben vermeiden.
Die Abfrage für negatives x gehört aber an den Anfang, nicht ans Ende. (-:
> Aber deine Idee ist nicht schlecht. Vielleicht veröffentliche ich sogar
> die reskalierte Sinusfunktion.
> double tsin(double x)
> {
> x -= trunc(x); /* kann man die Nachkommastellen direkt bekommen? */
Latürnich. Schrieb ich ja bereits, die C-Umwandlung double->int ist für
postive int ein Abschneiden.
> if (x > .5)
> x -= 1.;
> else if (x < -.5)
> x += 1.;
> //...
> }
>
> Und schon ist der Wert zwischen -0,5 und 0.5.
>
>>> if (x > 2*M_PI)
>>> x -= 2*M_PI*floor(x / 2*M_PI);
>> Ditto. Statt floor zu benutzen, kannst Du auch (für positive x) nach int
>> casten, dann ist die Abhängigkeit weg.
Was bei der ganzen Sache noch nicht stimmt: Da M_PI nicht exakt Pi ist,
sondern eben nur fast, wird durch die Skalierung eigentlich die Periode
des sin() verändert; sie ist dann M_PI und nicht Pi, was bedeutet, dass
Du für große Argumente hierdurch zusätzlich Genauigkeit verlierst. Ich
wüsste nicht, was man dagegen tun könnte.
> Stimmt, allerdings könnte es sein, dass ich mir damit in den Fuß
> schieße: Wenn der Compiler das ganze nämlich strikt ausführt, könnte es
> Datenverlust geben, weil in ein double (IEEE 754 double precision)
> wesentlich größere Zahlen passen als in ein int (in meinem Fall 8 Byte
> Signed Integer mit Zweierkomplement).
Alles richtig. Nur überlege, *wann* es Datenverlust gibt. Ein IEEE
double hat IIRC 52 Bit Mantisse, d.h. maximal kann es Ganzzahlen bis
2^52 korrekt darstellen. Bei Zahlen größer als dieser geht jegliche
Genauigkeit verloren, d.h. man könnte auch vorher eine Abfrage starten,
ob x > 2^52, falls ja, kann man eigentlich alles im Bereich [-1,1]
zurückliefern, denn das Resultat ist von Stellen abhängig, die die
Eingabe gar nicht darstellen kann. Falls nicht, so reicht immer noch ein
64-bit Integer aus.
> Ich habe meinen Compiler nebst Version genannt, und zumindest hier ist
> das Ergebnis deterministisch. Der ganze Kram funktioniert übrigens
> schlagartig, wenn ich den Compiler das ganze in Assembler übersetzen
> lasse, anschließend das Assembler-Listing entsprechend verändere (obwohl
> ich eigentlich AT&T-Syntax-Verächter bin) und anschließend vom
> Assembler-Listing an selber weiter übersetzen lasse.
>
> Ich habe übrigens die Compiler-Magie entdeckt: Der gcc war der Ansicht,
> dass sin(x) = -sin(-x), weswegen neben dem togglen des vordersten Bits
> des Argumentes auch das togglen des vordersten Bits der Rückgabe
> vonnötenwar.
Eben, der Compiler kennt die Symmetrie von sin(). Vorsorglich hatte ich
meine Funktion natürlich nicht sin() genannt.
Also, so als Hausaufgabe zum Nachgrübeln: Was macht man eigentlich für
sehr große Argumente, etwa wenn man vieleicht auch den cos() berechnen
will? Es bietet sich an, cos(x) = sin(x + Pi/2) zu verwenden. Nur für
sehr große x ist das eigentlich dasselbe, d.h. der Genauigkeitsverlust
ist so groß, dass x etwa x + Pi/2. Insbesondere hat man dann nicht mehr
sin^2 + cos^2 = 1, was vielleicht eine wünschenswerte Eigenschaft wäre.
Ich weiß nicht, was man da tun könnte.
Wie gesagt, das ganze ist ein höchst nichttriviales Problem. Ich würde
jetzt erstmal vorschlagen, in die Literatur zu gehen, etwa "Numerical
Recipes in C"(*). Ich habe derweil man in einige antike Quellen(**)
geguckt, und nein, es ist nicht nichttrivial wie das richtig
funktioniert. Der Trick mit dem Skalieren auf Integer, den ich
vorschlug, scheint aber auf jeden Fall der richtige Ansatz - wird dort
auch verwendet. Alles andere ist leider komplizierter.
Grüße,
Thomas
(*): Etwas Vorsicht mit dem Buch. Die Algorithmen aus dem Buch darfst Du
nicht in OpenSource Software einsetzen, siehe die Lizenzbestimmungen der
Autoren. Insofern gefährlich.
(**): Leider nicht zitationsfähig. Das sind Quellen, die ein (c) tragen,
und die ich nicht offiziell habe. Obwohl sich vermutlich nach den Jahren
niemand mehr darum schert.
Rainer Weikusat
[...]
> Davon abgesehen, au�er als akademische �bung erschlie�t sich mir der
> Sinn des ganzen nicht so recht. Der numerische Koprozessor, der in
> jeder halbwegs modernen Desktop-CPU eingebaut ist, hat eine
> Sinus-Approximation eingebaut.
Es gibt Leute, die der Ansicht sind, dass es sinnvoll waere, 'high
performance software floating point calculations' mit ARM-CPUs ohne
FPU zu betreiben :->. Ob deswegen, weil man dann den ganzen Kram, den
es ohnehin schon gibt, nochmal implementieren kann, weil man sonst
keine Probleme hat, die geloest werden muessten, sei
dahingestellt. Jedenfalls gibt es brownie points :->.
Markus Wichmann
> Markus Wichmann schrieb:
>> Thomas Richter (th...@math.tu-berlin.de) schrieb:
>
>>>> ich versuche hier gerade, den Sinus anzunähern (wird mal ein Algorithmus
>>>> für die dietlibc). OK, bisher hab ich schon interne Abhängigkeiten zu
>>>> floor(), ceil() und fabs(), was ja eigentlich zu vermeiden ist, aber
>>>> naja.
>>> Erstmal einige Punkte vorweg: Eine numerisch stabile Implementierung von
>>> transzendenten Funktionen gehört nicht zu den Dingen, die man im
>>> "Hau-Ruck" Verfahren implementiert. Falls das Resultat wirklich generell
>>> nutzbar sein soll, also für allgemeine Anwendungen taugen soll, so
>>> überlasse man das besser den Experten der Numerik.
>>>
>>> Einige wünscheswerte Eigenschaften wie Symmetrie oder Stetigkeit sind
>>> mit Sicherheit zu berücksichtigen, aber darüberhinaus auch noch das
>>> richtiges Grenzwertverhalten etc...
>>>
>>
>> Ich hab den Algorithmus von einer Seite im Internet, die damit warb,
>> dass dieser Algorithmus (immerhin lediglich ein biquadratisches Polynom)
>> einem Taylor-Polynom siebten Grades gerade in puncto Stetigkeit
>> überlegen ist.
>
> *Stetig* sind beide, denn Polynome sind nunmal stetig. Nicht nur das,
> sondern sogar analytisch. (-: Wichtig ist, wie die Polynome am Ende
> ihres Definitionsbereiches zusammengeklebt werden, so dass es passt.
>
Ja, und genau da war ja das Problem. Mit dem Taylor-Polynom gilt:
sin(Pi) = -0.075. Knapp daneben. Und wenn man es zusammenklebt, ist das
Taylorpolynom an der Stelle eben _nicht_ mehr stetig.
>> Es ist mehr oder weniger eine Parabel, die auf bestimmte
>> Punkte festgelegt ist. Wie es der Zufall will, ist die Parabel selbst
>> fast immer zu groß, während ihr Quadrat fast immer zu niedrig ist. Also
>> addiere ich beides gewichtet.
>
> Und das ist warum gut?
Weil ich auf die weise nur die Grundrechenarten brauche, um den Sinus zu
berechnen. Und wie du unten sehen wirst, ist das durchaus nötig.
> Wie ich schon schrieb, numerische Mathematik gehört nicht zu den
> Dingen, die man in Heimwerkermanier mit Schraubstock und Hammer löst.
> Dazu gehört schon ein gerüttelt Maß an Mathematik.
>
So, ich hab jetzt dein "gerüttelt Maß Mathematik" genutzt:
Erstmal hab ich eine Parabel mit den Fixpunkten 0|0, 0,25|1 und 0,5|0
gebastelt. Da kam heraus:
y = -16 * x^2 + 8 * x
Anschließend hab ich
integral(P*y^2 + (1-P)*y - sin(2*pi*x)) = 0
gesetzt. Es kam heraus: P = 5/37.
So, diese Konstante benutze ich jetzt auch. Dass die Fläche zwischen den
Funktionen 0 ist, bedeutet ja, dass sich die Fehler gegenseitig
aufheben.
> Davon abgesehen, außer als akademische Übung erschließt sich mir der
> Sinn des ganzen nicht so recht. Der numerische Koprozessor, der in jeder
> halbwegs modernen Desktop-CPU eingebaut ist, hat eine
> Sinus-Approximation eingebaut. Gebaut von Experten, also hoffentlich
> besser, als was man als Laie so zusammenschustern kann.
>
Muhaha.
Weißt du, Intel unterstelle ich, seit ich Assembler für AMD64 schreibe,
nicht mehr unbedingt, zu wissen, was sie tun. Du weißt ja vielleicht,
das x86-kompatible CPUs aus Kompatibilitätsgründen aus Unmengen an
Bausteinen bestehen. Da gibt es zum Bleistift die SSE-Einheit und die
ALU und die FPU. Super.
Das für meine Plattform gültige ABI definiert jetzt, dass double und
float in der SSE-Einheit übergeben werden. Dummerweise kann diese nur
die Grundrechenarten und Bitlogik ohne immediate data[*] und die
Quadratwurzel. Das wars. Wollte ich nun die FPU nutzen (da ist eine
Sinus-Approximation drin), müsste ich zunächst das Argument in den
Speicher verfrachten, von dort aus in die FPU laden, dort verarbeiten
(mit fsin und evtl. fprem1), anschließend zurück mit dem Ergebnis in den
Speicher und wieder rein in die SSE-Einheit. Mich schüttelts.
Außerdem ist der C-Quelltext im Gegensatz zum Assembler-Listing
portabel.
[*] Das bedeutet: Ich kann nicht sagen: Kipp mal das 63. Bit um, sondern
ich muss sagen:
- lade in ein GP-Register eine Zahl bei der nur das 63. Bit gesetzt ist,
- lade dieses GP-Register in den Speicher und
- führe ein xor von diesem Speicher auf das SSE-Register aus.
Unschön bis zum dorthinaus. Ich kann auch nicht sagen: Halbiere diese
Zahl, sondern ich muss irgendwo eine 2 oder 0.5 speichern.
>>> Zum Problem:
>>>
>>>> double sin(double x)
>>>> {
>>>> if (x < -2*M_PI)
>>>> x += 2*M_PI*ceil(-x / 2*M_PI);
>>> Hiermit bleibt das x in (-2Pi,0] immer noch negativ. Soll das so sein?
>>> Ich würde stattdessen vorschlagen, zunächst negative Eingaben abzufangen
>>> und -sin(-x) aufzurufen, damit die Symmetrie exakt stimmt.
>>>
>>
>> Das hab ich ja versucht, aber bei diesem Statement ist mein Compiler
>> bockig.
>
> Siehe unten. Dein Compiler ist eigentlich schlau und weiß, was eine
> Funktion namens sin() tun soll.
>
Ja, und er wollte schlauer sein, als ich.
>> Aber du hast recht, eigentlich will ich das x ja nach [-Pi; Pi]
>> skalieren.
>
> Eigentlich sogar nur auf [0,Pi/2).
>
Nein, auf [-Pi; Pi], denn auf dem Intervall gilt meine Approximation.
>> if (x < -2*M_PI)
>> x += 2*M_PI*(ceil(-x / 2*M_PI) + 1)
>>
>> ?
>
> Damit wird x aber auch nur auf das Interval (-2Pi,infinity) gemappt.
>
>>> Ferner vermute ich, dass die obige Berechnung des Modulus nicht recht
>>> stabil ist, da Du einmal hin, und dann wieder zurückskalierst, so dass
>>> Rundungsfehler akkumuliert werden. Ich würde stattdessen zunächst das
>>> Argument so skalieren, dass M_PI auf 1 fällt, und würde dann den Rest
>>> abseparieren. Dazu reicht C-Integerarithmetik, und floor() und ceil()
>>> entfallen. Das Approximationspolynom sollte dann natürlich entsprechend
>>> angepasst werden, besser noch mit einer Fallunterscheidung für sehr
>>> kleine Werte (dort genauer appoximieren) und für Werte in der Nähe von 1
>>> (ebenso genauer), der Rest durch ein Minimax-Polynom.
>>>
>>
>> Fallunterscheidungen sind fast immer Performance-Killer. Was meinst du,
>> warum ich so absonderliche Konstruktionen wie "x * fabs(x)" sonst
>> verwende?
>
> Und? Die Frage ist, worauf es bei dem Spaß ankommt. Wenn Geschwindigkeit
> eine Rolle spielt, würde ich zu einer Lookup-Tabelle greifen und, wenn's
> hoch kommt, zwischen den Einträgen interpolieren.
>
Es kommt in eine Lib rein. Die ist eigentlich für General
Purpose-Programmierung gedacht, d.h. sie sollte nicht ganz riesig und
auch nicht ganz lahm werden. Der beste Kompromiss, sozusagen.
>
>> Ich könnte natürlich auch
>>
>> if (x < 0)
>> y = -A * x * x + B * x;
>> else
>> y = A * x * x + B * x;
>>
>> verwenden, aber das wollte ich eben vermeiden.
>
> Die Abfrage für negatives x gehört aber an den Anfang, nicht ans Ende. (-:
>
Wieso? 3*Pi/2 < 0, das weiß doch jeder :-) (zumindest nach dem
skalieren, denn sin(3*PI/2) = sin(-Pi/2).)
>> if (x > .5)
>> x -= 1.;
>> else if (x < -.5)
>> x += 1.;
>> //...
>> }
>>
>> Und schon ist der Wert zwischen -0,5 und 0.5.
>>
>>>> if (x > 2*M_PI)
>>>> x -= 2*M_PI*floor(x / 2*M_PI);
>>> Ditto. Statt floor zu benutzen, kannst Du auch (für positive x) nach int
>>> casten, dann ist die Abhängigkeit weg.
>
> Was bei der ganzen Sache noch nicht stimmt: Da M_PI nicht exakt Pi ist,
> sondern eben nur fast, wird durch die Skalierung eigentlich die Periode
> des sin() verändert; sie ist dann M_PI und nicht Pi, was bedeutet, dass
> Du für große Argumente hierdurch zusätzlich Genauigkeit verlierst. Ich
> wüsste nicht, was man dagegen tun könnte.
>
Da hast du recht. Und da ist wohl wirklich kein Kraut dagegen gewachsen.
Ich meine, selbst wenn M_PI genauer definiert ist, als in ein double
hineinpasst, ist 10*M_PI wieder nicht die double-Approximation von
10*Pi.
Nun gut, für den Gutteil der Argumente funktioniert mein Sinus jetzt.
Ach ja: Ich habe jetzt eine Tabelle mit Funktionswerten. Hast du eine
Ahnung wo ich die gegenprüfen (lassen) kann, ob sie wirklich "nahe
genug" am Sinus dran sind?
>> Stimmt, allerdings könnte es sein, dass ich mir damit in den Fuß
>> schieße: Wenn der Compiler das ganze nämlich strikt ausführt, könnte es
>> Datenverlust geben, weil in ein double (IEEE 754 double precision)
>> wesentlich größere Zahlen passen als in ein int (in meinem Fall 8 Byte
>> Signed Integer mit Zweierkomplement).
>
> Alles richtig. Nur überlege, *wann* es Datenverlust gibt. Ein IEEE
> double hat IIRC 52 Bit Mantisse, d.h. maximal kann es Ganzzahlen bis
> 2^52 korrekt darstellen. Bei Zahlen größer als dieser geht jegliche
> Genauigkeit verloren, d.h. man könnte auch vorher eine Abfrage starten,
> ob x > 2^52, falls ja, kann man eigentlich alles im Bereich [-1,1]
> zurückliefern, denn das Resultat ist von Stellen abhängig, die die
> Eingabe gar nicht darstellen kann. Falls nicht, so reicht immer noch ein
> 64-bit Integer aus.
>
Wo du recht hast... verflixt, man sollte halt doch nachdenken, bevor man
etwas schreibt. Ich sollte an der Stelle wohl eine Precision Exception
auslösen und 0 zurückgeben, oder so.
>
>> Ich habe meinen Compiler nebst Version genannt, und zumindest hier ist
>> das Ergebnis deterministisch. Der ganze Kram funktioniert übrigens
>> schlagartig, wenn ich den Compiler das ganze in Assembler übersetzen
>> lasse, anschließend das Assembler-Listing entsprechend verändere (obwohl
>> ich eigentlich AT&T-Syntax-Verächter bin) und anschließend vom
>> Assembler-Listing an selber weiter übersetzen lasse.
>>
>> Ich habe übrigens die Compiler-Magie entdeckt: Der gcc war der Ansicht,
>> dass sin(x) = -sin(-x), weswegen neben dem togglen des vordersten Bits
>> des Argumentes auch das togglen des vordersten Bits der Rückgabe
>> vonnötenwar.
>
> Eben, der Compiler kennt die Symmetrie von sin(). Vorsorglich hatte ich
> meine Funktion natürlich nicht sin() genannt.
>
Wenn ich einen Sinus implementiere, dann nenne ich die Funktion nunmal
sin() :-).
> Also, so als Hausaufgabe zum Nachgrübeln: Was macht man eigentlich für
> sehr große Argumente, etwa wenn man vieleicht auch den cos() berechnen
> will? Es bietet sich an, cos(x) = sin(x + Pi/2) zu verwenden. Nur für
> sehr große x ist das eigentlich dasselbe, d.h. der Genauigkeitsverlust
> ist so groß, dass x etwa x + Pi/2. Insbesondere hat man dann nicht mehr
> sin^2 + cos^2 = 1, was vielleicht eine wünschenswerte Eigenschaft wäre.
> Ich weiß nicht, was man da tun könnte.
Erst das x skalieren und dann das halbe Pi addieren?
Wobei, 2*Pi hat 3 binäre Vorkommastellen. Das löst nur dann beim
Argument keine Änderung mehr aus, wenn x > 2^55. Da ich aber bei 2^52
eine Precision Exception auslösen muss, weil die Genauigkeit nicht mehr
ausreicht, kommt der Fall wohl nicht vor.
> Wie gesagt, das ganze ist ein höchst nichttriviales Problem. Ich würde
> jetzt erstmal vorschlagen, in die Literatur zu gehen, etwa "Numerical
> Recipes in C"(*). Ich habe derweil man in einige antike Quellen(**)
> geguckt, und nein, es ist nicht nichttrivial wie das richtig
> funktioniert. Der Trick mit dem Skalieren auf Integer, den ich
> vorschlug, scheint aber auf jeden Fall der richtige Ansatz - wird dort
> auch verwendet. Alles andere ist leider komplizierter.
>
Hmm... ja, stimmt. Eine Sache, die mir noch aufgefallen ist: Die meisten
Funktionswerte des Sinus sind transzendente Zahlen, d.h. es gibt kein
Polynom mit rationalen Koeffizienten, bei dem eine Nullstelle ein
solcher Funktionswert wäre. Die Approximierung über ein Polynom hingegen
sorgt dafür, dass alle Funktionswerte algebraisch sind.
> Grüße,
> Thomas
>
> (*): Etwas Vorsicht mit dem Buch. Die Algorithmen aus dem Buch darfst Du
> nicht in OpenSource Software einsetzen, siehe die Lizenzbestimmungen der
> Autoren. Insofern gefährlich.
>
Hrmpf... nicht gefährlich, sondern unbrauchbar: Die dietlibc steht unter
der LGPL, und da ich nicht fefe bin, werde ich das nicht ändern können.
> (**): Leider nicht zitationsfähig. Das sind Quellen, die ein (c) tragen,
> und die ich nicht offiziell habe. Obwohl sich vermutlich nach den Jahren
> niemand mehr darum schert.
Wie das halt immer so ist. In 75 Jahren kannst du es dann posten, oder
so...
Thomas Richter
> Ja, und genau da war ja das Problem. Mit dem Taylor-Polynom gilt:
> sin(Pi) = -0.075. Knapp daneben. Und wenn man es zusammenklebt, ist das
> Taylorpolynom an der Stelle eben _nicht_ mehr stetig.
Ehem:
\begin{Mathematiker}
Polynome sind immer stetig. Eine stückweise aus Polynomen bestehende
Funktion hingegen ist nicht unbedingt stetig.
\end{Mathematiker}
>>> Es ist mehr oder weniger eine Parabel, die auf bestimmte
>>> Punkte festgelegt ist. Wie es der Zufall will, ist die Parabel selbst
>>> fast immer zu groß, während ihr Quadrat fast immer zu niedrig ist. Also
>>> addiere ich beides gewichtet.
>> Und das ist warum gut?
>
> Weil ich auf die weise nur die Grundrechenarten brauche, um den Sinus zu
> berechnen. Und wie du unten sehen wirst, ist das durchaus nötig.
Achso? Also, da darf ich auch noch einige andere Algorithmen anbieten,
wie wäre es etwa mit dem CORDIC-Algorithmus?
>> Wie ich schon schrieb, numerische Mathematik gehört nicht zu den
>> Dingen, die man in Heimwerkermanier mit Schraubstock und Hammer löst.
>> Dazu gehört schon ein gerüttelt Maß an Mathematik.
>>
>
> So, ich hab jetzt dein "gerüttelt Maß Mathematik" genutzt:
>
> Erstmal hab ich eine Parabel mit den Fixpunkten 0|0, 0,25|1 und 0,5|0
> gebastelt. Da kam heraus:
>
> y = -16 * x^2 + 8 * x
>
> Anschließend hab ich
>
> integral(P*y^2 + (1-P)*y - sin(2*pi*x)) = 0
>
> gesetzt. Es kam heraus: P = 5/37.
>
> So, diese Konstante benutze ich jetzt auch. Dass die Fläche zwischen den
> Funktionen 0 ist, bedeutet ja, dass sich die Fehler gegenseitig
> aufheben.
Nö. Schau mal, das bedeutet nur, dass der mittlere Fehler null ist, was
aber nicht viel heißt. Da kann ich noch jede beliebige Funktion mit
Mittelwert 0 'draufaddieren (etwa 10^6 * sin(20*x)!), und die
Approximation bleibt "genauso gut". Was sie natürlich nicht ist.
Mal einige Antworten aus Sicht eines Mathematikers:
Erstmal gilt es zu überlegen, in welchem Sinne Du den Sinus überhaupt
approximieren willst. Da gibt es diverse Möglichkeiten, und je nach
Fragestellung gibt es diverse Antworten.
Eine Möglichkeit wäre, den mittleren *quadratischen* Fehler zu
minimieren, also |sin(x)-f(x)|^2 minimal. Hier wäre die richtige
Antwort, den Sin auf dem Interval durch Legendre-Polynome zu
approximieren. Der *mittlere* Fehler sin(x)-f(x) (ohne Betrag) ist
offenbar *kein* Fehlermaß, denn wenn er null ist, ist das Resultat noch
lange nicht gut.
Üblicherweise möchte man aber den l^\infty-Fehler minimieren, will
sagen, der maximale Fehler soll minimal werden, also sup |sin(x)-f(x)| =
min. Hierzu hilft ein Integral leider wenig. Das hierzu nötige Polynom
heißt "Minimax-Polynom" und es gibt fertige Algorithmen, die Dir das
ausrechnen. Diese Algorithmen sind in mathematischer Software wie Maple
oder Matlab integriert und können benutzt werden.
Google-> Minimax-Polynom
Google-> CORDIC
Darüberhinaus hängt das Minimax-Polynom natürlich von der Intervalgröße
ab, auf der man interpolieren will, und vermutlich ist es günstiger,
dieses Intervall auch zusammenzustückeln. Etwa will man in der Nähe von
0 das Verhalten sin(x) approx x reproduzieren, ebenso wie cos(x) = 1 -
x^2/2. Beides wird das Minimax-Polynom "naiv" nicht leisten, so dass man
hier "stückelt". Vermutlich braucht man drei Abschnitte, um das recht zu
machen.
>> Davon abgesehen, außer als akademische Übung erschließt sich mir der
>> Sinn des ganzen nicht so recht. Der numerische Koprozessor, der in jeder
>> halbwegs modernen Desktop-CPU eingebaut ist, hat eine
>> Sinus-Approximation eingebaut. Gebaut von Experten, also hoffentlich
>> besser, als was man als Laie so zusammenschustern kann.
>>
>
> Muhaha.
>
> Weißt du, Intel unterstelle ich, seit ich Assembler für AMD64 schreibe,
> nicht mehr unbedingt, zu wissen, was sie tun. Du weißt ja vielleicht,
> das x86-kompatible CPUs aus Kompatibilitätsgründen aus Unmengen an
> Bausteinen bestehen. Da gibt es zum Bleistift die SSE-Einheit und die
> ALU und die FPU. Super.
Und? Was sagt das über die Qualität der Mathematiker bei intel? Schau,
der Ingenieur, der die Transistoren plaziert (oder eher das VHDL
schreibt), wird auch nur die Gatter nehmen, die er vom Hausmathematiker
bekommt. Die Algorithmen sind über Jahre erprobt worden, bzw. sind aus
entsprechender Fachliteratur entnommen und vermutlich von Hand optimiert
worden.
> Das für meine Plattform gültige ABI definiert jetzt, dass double und
> float in der SSE-Einheit übergeben werden. Dummerweise kann diese nur
> die Grundrechenarten und Bitlogik ohne immediate data[*] und die
> Quadratwurzel. Das wars. Wollte ich nun die FPU nutzen (da ist eine
> Sinus-Approximation drin), müsste ich zunächst das Argument in den
> Speicher verfrachten, von dort aus in die FPU laden, dort verarbeiten
> (mit fsin und evtl. fprem1), anschließend zurück mit dem Ergebnis in den
> Speicher und wieder rein in die SSE-Einheit. Mich schüttelts.
Und? warum? Wo ist da das Problem?
> Außerdem ist der C-Quelltext im Gegensatz zum Assembler-Listing
> portabel.
Portabel auf was? Will sagen, jetzt, wo Du die Probleme der
Approximation kennst, sollte klar sein, dass die Polynome, die Du
brauchst, auch von der Rechengenauigkeit abhängen. Die garantiert Dir C
aber überhaupt nicht.
>> Und? Die Frage ist, worauf es bei dem Spaß ankommt. Wenn Geschwindigkeit
>> eine Rolle spielt, würde ich zu einer Lookup-Tabelle greifen und, wenn's
>> hoch kommt, zwischen den Einträgen interpolieren.
>>
>
> Es kommt in eine Lib rein. Die ist eigentlich für General
> Purpose-Programmierung gedacht, d.h. sie sollte nicht ganz riesig und
> auch nicht ganz lahm werden. Der beste Kompromiss, sozusagen.
Hmm. Für wen? Für mich ist der beste Kompromiss, ein autoconf-Skript zu
benutzen welches mir die Verfügbarkeit einer built-in "sin" ermittelt,
von mir aus mit gcc-Syntax. Wenn es das nicht gibt, hängt der Kompromis
stark davon ab, was mein Problem ist.
> Nun gut, für den Gutteil der Argumente funktioniert mein Sinus jetzt.
> Ach ja: Ich habe jetzt eine Tabelle mit Funktionswerten. Hast du eine
> Ahnung wo ich die gegenprüfen (lassen) kann, ob sie wirklich "nahe
> genug" am Sinus dran sind?
Logo. Maple oder Matlab. Berechnet Dir den sin(x) bei Bedarf auf 100
Stellen genau aus.
>> Alles richtig. Nur überlege, *wann* es Datenverlust gibt. Ein IEEE
>> double hat IIRC 52 Bit Mantisse, d.h. maximal kann es Ganzzahlen bis
>> 2^52 korrekt darstellen. Bei Zahlen größer als dieser geht jegliche
>> Genauigkeit verloren, d.h. man könnte auch vorher eine Abfrage starten,
>> ob x > 2^52, falls ja, kann man eigentlich alles im Bereich [-1,1]
>> zurückliefern, denn das Resultat ist von Stellen abhängig, die die
>> Eingabe gar nicht darstellen kann. Falls nicht, so reicht immer noch ein
>> 64-bit Integer aus.
>>
>
> Wo du recht hast... verflixt, man sollte halt doch nachdenken, bevor man
> etwas schreibt. Ich sollte an der Stelle wohl eine Precision Exception
> auslösen und 0 zurückgeben, oder so.
sin() wirft keine Exceptions nach C-Konvention. Vermutlich darf man
nicht einmal ein NAN zurückliefern, was sich eigentlich anböte. Mein
altes "Motorola"-Manual sagt dazu nur, dass für große Argumente jegliche
Genauigkeit verloren geht, was auch immer das heißen mag.
>> Also, so als Hausaufgabe zum Nachgrübeln: Was macht man eigentlich für
>> sehr große Argumente, etwa wenn man vieleicht auch den cos() berechnen
>> will? Es bietet sich an, cos(x) = sin(x + Pi/2) zu verwenden. Nur für
>> sehr große x ist das eigentlich dasselbe, d.h. der Genauigkeitsverlust
>> ist so groß, dass x etwa x + Pi/2. Insbesondere hat man dann nicht mehr
>> sin^2 + cos^2 = 1, was vielleicht eine wünschenswerte Eigenschaft wäre.
>> Ich weiß nicht, was man da tun könnte.
>
>
> Erst das x skalieren und dann das halbe Pi addieren?
Hmm. Möglich. Zumindest konsistent.
>> Wie gesagt, das ganze ist ein höchst nichttriviales Problem. Ich würde
>> jetzt erstmal vorschlagen, in die Literatur zu gehen, etwa "Numerical
>> Recipes in C"(*). Ich habe derweil man in einige antike Quellen(**)
>> geguckt, und nein, es ist nicht nichttrivial wie das richtig
>> funktioniert. Der Trick mit dem Skalieren auf Integer, den ich
>> vorschlug, scheint aber auf jeden Fall der richtige Ansatz - wird dort
>> auch verwendet. Alles andere ist leider komplizierter.
>>
>
> Hmm... ja, stimmt. Eine Sache, die mir noch aufgefallen ist: Die meisten
> Funktionswerte des Sinus sind transzendente Zahlen, d.h. es gibt kein
> Polynom mit rationalen Koeffizienten, bei dem eine Nullstelle ein
> solcher Funktionswert wäre. Die Approximierung über ein Polynom hingegen
> sorgt dafür, dass alle Funktionswerte algebraisch sind.
Der Rechner kann keine transzendenten Zahlen. Nicht mal reelle. Nicht
mal alle rationalen. Nur einen Bruchteil der rationalen. Insofern ist
das noch das geringste Problem.
Achja, mir fällt noch eine Quelle ein, falls es Dich interessiert: Es
gibt von Motorla die FPSP (Floating Point Software Package) für den
68060 und 68040. Beide Prozessoren machten den Sinus nicht mehr in
Hardware, sondern in Software über emulierte Befehle. Die Quellen habe
ich, und sie wurden von Mot damals herausgegeben. Falls Du also MC68K
Assembler lesen kannst, wäre das auch eine Möglichkeit, und *diese*
Quelle ist benutzbar. Vielleicht auch nur als Abschreckung. (-: Meine
Mail ist gültig, wenn Du den Algorithmus dort haben willst...
Aber, wie ich schon sagte, das ganze ist kein Spaziergang, und der
Algorithmus ist ziemlich "convoluted".
Grüße,
Thomas
Markus Wichmann
> Markus Wichmann schrieb:
>
>> Ja, und genau da war ja das Problem. Mit dem Taylor-Polynom gilt:
>> sin(Pi) = -0.075. Knapp daneben. Und wenn man es zusammenklebt, ist das
>> Taylorpolynom an der Stelle eben _nicht_ mehr stetig.
>
> Ehem:
> \begin{Mathematiker}
> Polynome sind immer stetig. Eine stückweise aus Polynomen bestehende
> Funktion hingegen ist nicht unbedingt stetig.
> \end{Mathematiker}
>
Eben dieses meinte ich. Wenn du den Krümelkackermodus haben willst:
Sei f_n: R -> R ein Taylorpolynom n-ten Grades. Es sei
\forall x \in R: \lim_{n->\infty} f_n(x) = \sin x
Sei weiterhin g: R -> R definiert mit
g(x) = f_7(x) für x \in [-\pi; \pi] und
g(x) = g(x + 2\pi) sonst
Dann gilt: g ist an der Stelle x = Pi unstetig. Beweis:
\lim_{x->\Pi+0} g(x) = 0, aber
\lim_{x->\Pi-0} g(x) = -0.075
Krümelkackermodus aus.
>>>> Es ist mehr oder weniger eine Parabel, die auf bestimmte
>>>> Punkte festgelegt ist. Wie es der Zufall will, ist die Parabel selbst
>>>> fast immer zu groß, während ihr Quadrat fast immer zu niedrig ist. Also
>>>> addiere ich beides gewichtet.
>>> Und das ist warum gut?
>>
>> Weil ich auf die weise nur die Grundrechenarten brauche, um den Sinus zu
>> berechnen. Und wie du unten sehen wirst, ist das durchaus nötig.
>
> Achso? Also, da darf ich auch noch einige andere Algorithmen anbieten,
> wie wäre es etwa mit dem CORDIC-Algorithmus?
>
Was ist das? Ach, moment, Wikipedia...
>>> Wie ich schon schrieb, numerische Mathematik gehört nicht zu den
>>> Dingen, die man in Heimwerkermanier mit Schraubstock und Hammer löst.
>>> Dazu gehört schon ein gerüttelt Maß an Mathematik.
>>>
>>
>> So, ich hab jetzt dein "gerüttelt Maß Mathematik" genutzt:
>>
>> Erstmal hab ich eine Parabel mit den Fixpunkten 0|0, 0,25|1 und 0,5|0
>> gebastelt. Da kam heraus:
>>
>> y = -16 * x^2 + 8 * x
>>
>> Anschließend hab ich
>>
>> integral(P*y^2 + (1-P)*y - sin(2*pi*x)) = 0
>>
>> gesetzt. Es kam heraus: P = 5/37.
>>
>> So, diese Konstante benutze ich jetzt auch. Dass die Fläche zwischen den
>> Funktionen 0 ist, bedeutet ja, dass sich die Fehler gegenseitig
>> aufheben.
>
> Nö. Schau mal, das bedeutet nur, dass der mittlere Fehler null ist, was
> aber nicht viel heißt. Da kann ich noch jede beliebige Funktion mit
> Mittelwert 0 'draufaddieren (etwa 10^6 * sin(20*x)!), und die
> Approximation bleibt "genauso gut". Was sie natürlich nicht ist.
>
Wo du recht hast...
> Mal einige Antworten aus Sicht eines Mathematikers:
>
> Erstmal gilt es zu überlegen, in welchem Sinne Du den Sinus überhaupt
> approximieren willst. Da gibt es diverse Möglichkeiten, und je nach
> Fragestellung gibt es diverse Antworten.
>
Und lass mich raten: "General Purpose" steht nicht zur Debatte?
> Eine Möglichkeit wäre, den mittleren *quadratischen* Fehler zu
> minimieren, also |sin(x)-f(x)|^2 minimal. Hier wäre die richtige
> Antwort, den Sin auf dem Interval durch Legendre-Polynome zu
> approximieren. Der *mittlere* Fehler sin(x)-f(x) (ohne Betrag) ist
> offenbar *kein* Fehlermaß, denn wenn er null ist, ist das Resultat noch
> lange nicht gut.
>
> Üblicherweise möchte man aber den l^\infty-Fehler minimieren, will
> sagen, der maximale Fehler soll minimal werden, also sup |sin(x)-f(x)| =
> min. Hierzu hilft ein Integral leider wenig. Das hierzu nötige Polynom
> heißt "Minimax-Polynom" und es gibt fertige Algorithmen, die Dir das
> ausrechnen. Diese Algorithmen sind in mathematischer Software wie Maple
> oder Matlab integriert und können benutzt werden.
>
> Google-> Minimax-Polynom
> Google-> CORDIC
>
> Darüberhinaus hängt das Minimax-Polynom natürlich von der Intervalgröße
> ab, auf der man interpolieren will, und vermutlich ist es günstiger,
> dieses Intervall auch zusammenzustückeln. Etwa will man in der Nähe von
> 0 das Verhalten sin(x) approx x reproduzieren, ebenso wie cos(x) = 1 -
> x^2/2. Beides wird das Minimax-Polynom "naiv" nicht leisten, so dass man
> hier "stückelt". Vermutlich braucht man drei Abschnitte, um das recht zu
> machen.
>
Uff. Und das nur wegen Faulheit, den Speicher zu belasten.
>>> Davon abgesehen, außer als akademische Übung erschließt sich mir der
>>> Sinn des ganzen nicht so recht. Der numerische Koprozessor, der in jeder
>>> halbwegs modernen Desktop-CPU eingebaut ist, hat eine
>>> Sinus-Approximation eingebaut. Gebaut von Experten, also hoffentlich
>>> besser, als was man als Laie so zusammenschustern kann.
>>>
>>
>> Muhaha.
>>
>> Weißt du, Intel unterstelle ich, seit ich Assembler für AMD64 schreibe,
>> nicht mehr unbedingt, zu wissen, was sie tun. Du weißt ja vielleicht,
>> das x86-kompatible CPUs aus Kompatibilitätsgründen aus Unmengen an
>> Bausteinen bestehen. Da gibt es zum Bleistift die SSE-Einheit und die
>> ALU und die FPU. Super.
>
> Und? Was sagt das über die Qualität der Mathematiker bei intel? Schau,
> der Ingenieur, der die Transistoren plaziert (oder eher das VHDL
> schreibt), wird auch nur die Gatter nehmen, die er vom Hausmathematiker
> bekommt. Die Algorithmen sind über Jahre erprobt worden, bzw. sind aus
> entsprechender Fachliteratur entnommen und vermutlich von Hand optimiert
> worden.
>
Das mag ja sein, aber wenn der Ingenieur blöde vorgaben bekommt...
>> Das für meine Plattform gültige ABI definiert jetzt, dass double und
>> float in der SSE-Einheit übergeben werden. Dummerweise kann diese nur
>> die Grundrechenarten und Bitlogik ohne immediate data[*] und die
>> Quadratwurzel. Das wars. Wollte ich nun die FPU nutzen (da ist eine
>> Sinus-Approximation drin), müsste ich zunächst das Argument in den
>> Speicher verfrachten, von dort aus in die FPU laden, dort verarbeiten
>> (mit fsin und evtl. fprem1), anschließend zurück mit dem Ergebnis in den
>> Speicher und wieder rein in die SSE-Einheit. Mich schüttelts.
>
> Und? warum? Wo ist da das Problem?
>
Die Transistoren, die den Sinus berechnen, sitzen an der falschen
Stelle, und die einzige Verbindungsleitung ist ausreichend dünn.
>> Außerdem ist der C-Quelltext im Gegensatz zum Assembler-Listing
>> portabel.
>
> Portabel auf was? Will sagen, jetzt, wo Du die Probleme der
> Approximation kennst, sollte klar sein, dass die Polynome, die Du
> brauchst, auch von der Rechengenauigkeit abhängen. Die garantiert Dir C
> aber überhaupt nicht.
>
Nein, C garantiert mir die nicht, aber so weit wollte ich auch gar nicht
gehen. Ich habe für meinen Dunstkreis weitere Garantien: Die ABI meines
Systems und aller Systeme, zu denen ich kompatibel zu sein wünsche. Und
diese ABIs definieren, dass es eine x87-FPU gibt.
>>> Und? Die Frage ist, worauf es bei dem Spaß ankommt. Wenn Geschwindigkeit
>>> eine Rolle spielt, würde ich zu einer Lookup-Tabelle greifen und, wenn's
>>> hoch kommt, zwischen den Einträgen interpolieren.
>>>
>>
>> Es kommt in eine Lib rein. Die ist eigentlich für General
>> Purpose-Programmierung gedacht, d.h. sie sollte nicht ganz riesig und
>> auch nicht ganz lahm werden. Der beste Kompromiss, sozusagen.
>
> Hmm. Für wen? Für mich ist der beste Kompromiss, ein autoconf-Skript zu
> benutzen welches mir die Verfügbarkeit einer built-in "sin" ermittelt,
> von mir aus mit gcc-Syntax.
gcc-Syntax?
#include <stdio.h>
int main()
{
printf("%d\n", __builtin_sin(0.75));
return 0;
}
ist gcc-Syntax?
> Wenn es das nicht gibt, hängt der Kompromis
> stark davon ab, was mein Problem ist.
>
>> Nun gut, für den Gutteil der Argumente funktioniert mein Sinus jetzt.
>> Ach ja: Ich habe jetzt eine Tabelle mit Funktionswerten. Hast du eine
>> Ahnung wo ich die gegenprüfen (lassen) kann, ob sie wirklich "nahe
>> genug" am Sinus dran sind?
>
> Logo. Maple oder Matlab. Berechnet Dir den sin(x) bei Bedarf auf 100
> Stellen genau aus.
>
>>> Alles richtig. Nur überlege, *wann* es Datenverlust gibt. Ein IEEE
>>> double hat IIRC 52 Bit Mantisse, d.h. maximal kann es Ganzzahlen bis
>>> 2^52 korrekt darstellen. Bei Zahlen größer als dieser geht jegliche
>>> Genauigkeit verloren, d.h. man könnte auch vorher eine Abfrage starten,
>>> ob x > 2^52, falls ja, kann man eigentlich alles im Bereich [-1,1]
>>> zurückliefern, denn das Resultat ist von Stellen abhängig, die die
>>> Eingabe gar nicht darstellen kann. Falls nicht, so reicht immer noch ein
>>> 64-bit Integer aus.
>>>
>>
>> Wo du recht hast... verflixt, man sollte halt doch nachdenken, bevor man
>> etwas schreibt. Ich sollte an der Stelle wohl eine Precision Exception
>> auslösen und 0 zurückgeben, oder so.
>
> sin() wirft keine Exceptions nach C-Konvention. Vermutlich darf man
> nicht einmal ein NAN zurückliefern, was sich eigentlich anböte. Mein
> altes "Motorola"-Manual sagt dazu nur, dass für große Argumente jegliche
> Genauigkeit verloren geht, was auch immer das heißen mag.
>
Exception im Sinne von Floating Point Exception kann sin() durchaus
werfen, so wie jede andere Floating Point-Anweisung. Ich meine, was soll
den sin(NAN) sein? sin(NaNq) = NaNq und keine Invalid Exception,
sin(NaNs) = NaNs und Invalid Exception sind doch gültige
Implementierungen, oder? Da liegt es doch nahe, eine Precision Exception
auszulösen, wenn das Ergebnis unexakt ist, oder?
Dummerweise bezieht sich 7.12.1 größtenteils auf das Ergebnis einer
Operation. Allerdings ist ja 7.12.1p5 anwendbar, wenn man die Sache so
interpretiert, dass mit dem kleinen Ergebnis die Größe in Relation zum
Argument gemeint ist. Von einer Precision Exception steht dort nix, also
ist nicht definiert, ob eine kommt oder nicht. Es steht allerdings da,
dass, wenn math_errhandling & MATH_ERREXCEPT ungleich Null ist, es
implementationsdefiniert ist, ob eine Underflow Exception ausgelöst
wird.
>>> Also, so als Hausaufgabe zum Nachgrübeln: Was macht man eigentlich für
>>> sehr große Argumente, etwa wenn man vieleicht auch den cos() berechnen
>>> will? Es bietet sich an, cos(x) = sin(x + Pi/2) zu verwenden. Nur für
>>> sehr große x ist das eigentlich dasselbe, d.h. der Genauigkeitsverlust
>>> ist so groß, dass x etwa x + Pi/2. Insbesondere hat man dann nicht mehr
>>> sin^2 + cos^2 = 1, was vielleicht eine wünschenswerte Eigenschaft wäre.
>>> Ich weiß nicht, was man da tun könnte.
>>
>>
>> Erst das x skalieren und dann das halbe Pi addieren?
>
> Hmm. Möglich. Zumindest konsistent.
>
>>> Wie gesagt, das ganze ist ein höchst nichttriviales Problem. Ich würde
>>> jetzt erstmal vorschlagen, in die Literatur zu gehen, etwa "Numerical
>>> Recipes in C"(*). Ich habe derweil man in einige antike Quellen(**)
>>> geguckt, und nein, es ist nicht nichttrivial wie das richtig
>>> funktioniert. Der Trick mit dem Skalieren auf Integer, den ich
>>> vorschlug, scheint aber auf jeden Fall der richtige Ansatz - wird dort
>>> auch verwendet. Alles andere ist leider komplizierter.
>>>
>>
>> Hmm... ja, stimmt. Eine Sache, die mir noch aufgefallen ist: Die meisten
>> Funktionswerte des Sinus sind transzendente Zahlen, d.h. es gibt kein
>> Polynom mit rationalen Koeffizienten, bei dem eine Nullstelle ein
>> solcher Funktionswert wäre. Die Approximierung über ein Polynom hingegen
>> sorgt dafür, dass alle Funktionswerte algebraisch sind.
>
> Der Rechner kann keine transzendenten Zahlen. Nicht mal reelle.
Ich hoffe, du weist, dass die reellen Zahlen eine Obermenge der
transzendenten Zahlen sind. Genauer sind die reellen Zahlen die
Vereinigung der algebraischen und transzendenten Zahlen.
> Nicht
> mal alle rationalen. Nur einen Bruchteil der rationalen. Insofern ist
> das noch das geringste Problem.
>
Stimmt auch wieder. Mit ist nur eingefallen, dass das ein Beweis dafür
ist, dass der Algorithmus außer an den ausgewählten Stellen nie exakt
sein kann. Prinzipbedingt.
> Achja, mir fällt noch eine Quelle ein, falls es Dich interessiert: Es
> gibt von Motorla die FPSP (Floating Point Software Package) für den
> 68060 und 68040. Beide Prozessoren machten den Sinus nicht mehr in
> Hardware, sondern in Software über emulierte Befehle. Die Quellen habe
> ich, und sie wurden von Mot damals herausgegeben. Falls Du also MC68K
> Assembler lesen kannst, wäre das auch eine Möglichkeit, und *diese*
> Quelle ist benutzbar. Vielleicht auch nur als Abschreckung. (-: Meine
> Mail ist gültig, wenn Du den Algorithmus dort haben willst...
>
> Aber, wie ich schon sagte, das ganze ist kein Spaziergang, und der
> Algorithmus ist ziemlich "convoluted".
>
Bei den Aussichten lehne ich das Angebot dankend ab. Ich hab nämlich
nicht den Hauch einer Ahnung vom Innenleben eines MC68K-Prozessors,
geschweige denn seiner Assembler-Sprache.
> Grüße,
> Thomas
Stefan Reuther
> Thomas Richter (th...@math.tu-berlin.de) schrieb:
>>Davon abgesehen, außer als akademische Übung erschließt sich mir der
>>Sinn des ganzen nicht so recht. Der numerische Koprozessor, der in jeder
>>halbwegs modernen Desktop-CPU eingebaut ist, hat eine
>>Sinus-Approximation eingebaut. Gebaut von Experten, also hoffentlich
>>besser, als was man als Laie so zusammenschustern kann.
>
> Muhaha.
>
> Weißt du, Intel unterstelle ich, seit ich Assembler für AMD64 schreibe,
> nicht mehr unbedingt, zu wissen, was sie tun. Du weißt ja vielleicht,
> das x86-kompatible CPUs aus Kompatibilitätsgründen aus Unmengen an
> Bausteinen bestehen. Da gibt es zum Bleistift die SSE-Einheit und die
> ALU und die FPU. Super.
>
> Das für meine Plattform gültige ABI definiert jetzt, dass double und
> float in der SSE-Einheit übergeben werden. Dummerweise kann diese nur
> die Grundrechenarten und Bitlogik ohne immediate data[*] und die
> Quadratwurzel. Das wars. Wollte ich nun die FPU nutzen (da ist eine
> Sinus-Approximation drin), müsste ich zunächst das Argument in den
> Speicher verfrachten, von dort aus in die FPU laden, dort verarbeiten
> (mit fsin und evtl. fprem1), anschließend zurück mit dem Ergebnis in den
> Speicher und wieder rein in die SSE-Einheit. Mich schüttelts.
Ganz ehrlich: wenn du *deswegen* den Sinus so implementierst, wie in
deinem ersten Posting, würde ich dir auch nicht unterstellen, zu wissen,
was du tust. ceil, floor - mich schüttelt's. Wir haben doch auch fmod.
> Außerdem ist der C-Quelltext im Gegensatz zum Assembler-Listing
> portabel.
Das Assembler-fsin ist aber, was die Performance angeht (Geschwindigkeit
wie Genauigkeit) das Maß der Dinge und durch eine handgestrickte
Polynom-Approximation nicht zu übertreffen.
> [*] Das bedeutet: Ich kann nicht sagen: Kipp mal das 63. Bit um, sondern
> ich muss sagen:
>
> - lade in ein GP-Register eine Zahl bei der nur das 63. Bit gesetzt ist,
> - lade dieses GP-Register in den Speicher und
> - führe ein xor von diesem Speicher auf das SSE-Register aus.
Die "Zahl, in der un das 63. Bit gesetzt ist" kannst du auch gleich aus
dem Speicher holen. Genau wie die Floating-Point-Literale beim x87.
>>>>Hiermit bleibt das x in (-2Pi,0] immer noch negativ. Soll das so sein?
>>>>Ich würde stattdessen vorschlagen, zunächst negative Eingaben abzufangen
>>>>und -sin(-x) aufzurufen, damit die Symmetrie exakt stimmt.
>>>
>>>Das hab ich ja versucht, aber bei diesem Statement ist mein Compiler
>>>bockig.
>>
>>Siehe unten. Dein Compiler ist eigentlich schlau und weiß, was eine
>>Funktion namens sin() tun soll.
>
> Ja, und er wollte schlauer sein, als ich.
Du hast ihm vermutlich nicht gesagt, dass er die Finger davon lassen
solle. Ein "normales" C-Programm (hosted environment) darf halt keine
Funktion 'sin' umdefinieren.
>>>Fallunterscheidungen sind fast immer Performance-Killer. Was meinst du,
>>>warum ich so absonderliche Konstruktionen wie "x * fabs(x)" sonst
>>>verwende?
>>
>>Und? Die Frage ist, worauf es bei dem Spaß ankommt. Wenn Geschwindigkeit
>>eine Rolle spielt, würde ich zu einer Lookup-Tabelle greifen und, wenn's
>>hoch kommt, zwischen den Einträgen interpolieren.
>
> Es kommt in eine Lib rein. Die ist eigentlich für General
> Purpose-Programmierung gedacht, d.h. sie sollte nicht ganz riesig und
> auch nicht ganz lahm werden. Der beste Kompromiss, sozusagen.
Hast du denn wenigstens nachgemessen, dass das durch 'fabs' ersetzte
'if' auch tatsächlich eine kostspielige Fallunterscheidung ergibt? Das
ist heutzutage schon lange nicht mehr selbstverständlich.
>>Was bei der ganzen Sache noch nicht stimmt: Da M_PI nicht exakt Pi ist,
>>sondern eben nur fast, wird durch die Skalierung eigentlich die Periode
>>des sin() verändert; sie ist dann M_PI und nicht Pi, was bedeutet, dass
>>Du für große Argumente hierdurch zusätzlich Genauigkeit verlierst. Ich
>>wüsste nicht, was man dagegen tun könnte.
>
> Da hast du recht. Und da ist wohl wirklich kein Kraut dagegen gewachsen.
> Ich meine, selbst wenn M_PI genauer definiert ist, als in ein double
> hineinpasst, ist 10*M_PI wieder nicht die double-Approximation von
> 10*Pi.
Wenn du das in Assembler machst, schon. Weil 'fldpi' mehr Genauigkeit
generiert als double hat.
> Nun gut, für den Gutteil der Argumente funktioniert mein Sinus jetzt.
> Ach ja: Ich habe jetzt eine Tabelle mit Funktionswerten. Hast du eine
> Ahnung wo ich die gegenprüfen (lassen) kann, ob sie wirklich "nahe
> genug" am Sinus dran sind?
Ich würde mit der Ausgabe des Assemblerbefehls oder von MFPR
vergleichen. Allerdings würde ich mich für eine allgemeine Library nicht
mit Abweichungen zufriedengeben.
>>Wie gesagt, das ganze ist ein höchst nichttriviales Problem. Ich würde
>>jetzt erstmal vorschlagen, in die Literatur zu gehen, etwa "Numerical
>>Recipes in C"(*). [...]
>>(*): Etwas Vorsicht mit dem Buch. Die Algorithmen aus dem Buch darfst Du
>>nicht in OpenSource Software einsetzen, siehe die Lizenzbestimmungen der
>>Autoren. Insofern gefährlich.
>
> Hrmpf... nicht gefährlich, sondern unbrauchbar: Die dietlibc steht unter
> der LGPL, und da ich nicht fefe bin, werde ich das nicht ändern können.
Man darf die Algorithmen nicht rüberkopieren, aber ich wüsste nichts,
was dagegenspräche, sie nachzuimplementieren.
Stefan
Georg Bauhaus
>> (*): Etwas Vorsicht mit dem Buch. Die Algorithmen aus dem
>> Buch darfst Du nicht in OpenSource Software einsetzen, siehe
>> die Lizenzbestimmungen der Autoren. Insofern gef�hrlich.
>
> Vermutlich sind die Quelltexte gemeint, also spezielle
> /Formulierungen/ von Algorithmen in einer
> Programmiersprache; die beschrieben (abstrakten) Algorithmen
> k�nnen ja selber nicht in dieser Weise gesch�tzt werden.
>
(OT) Look & Feel, 1-Click.
Georg Bauhaus
> gesch�tzt, sondern patentrechtlich (insofern sie �berhaupt einer
> Anfechtung standhielten). Die Autoren des erw�hnten Buchs nahmen
> aber gar keinen patentrechtlichen Schutz f�r ihre Algorithmen in
> Anspruch, soweit der Vorposter dies berichtet hat.
>
Stimmt schon, nur braucht man zur wie immer gesch�tzten
Nutzung i.A. einen Vertrag...
(1-Click ist m.W. in EU nicht mehr standfest.)
Georg Bauhaus
> (Algorithmen) eben gerade nicht.
Dann ist es kein gesch�tztes Wissen.
Im Geltungsbereich des deutschen Urheberrechts
Software zu entwickeln ist eine Sache f�r sich,
nicht f�r die Welt... die anders damit umgeht.
Sinnhaft ist f�r Patentanw�lte, Rechtsabteilungen,
oder Betriebswirte je etwas anderes, als f�r Didakten,
weil jeweis andere Rationalit�t entscheidet.
Wir wissen ja nicht mal, was Wissen ist, sonst k�nnten
wir das allgemein kodifizieren.
Nur was ein Werk ist, wird gelegentlich zu bestimmen
versucht und daf�r braucht man vor Gericht eine Rechtsvertretung.
Thomas Richter
>>> Nutzung i.A. einen Vertrag...
>> Solch eines Vertrages bedarf es zur Nutzung von Wissen
>> (Algorithmen) eben gerade nicht.
>
> Im Geltungsbereich des deutschen Urheberrechts
> nicht f�r die Welt... die anders damit umgeht.
>
> Sinnhaft ist f�r Patentanw�lte, Rechtsabteilungen,
> oder Betriebswirte je etwas anderes, als f�r Didakten,
> weil jeweis andere Rationalit�t entscheidet.
>
> Wir wissen ja nicht mal, was Wissen ist, sonst k�nnten
> wir das allgemein kodifizieren.
> Nur was ein Werk ist, wird gelegentlich zu bestimmen
Wer sich die Details antun m�chte, hier die Kurzzusammenfassung:
---
When Do I Need a License?
The Numerical Recipes� book, and its electronic subscription version
Numerical Recipes Electronic, are sold as text and reference works, for
reading and educational purposes. If you intend to use Numerical Recipes
code on a computer, you need a separate license.
What Kinds of Licenses Are There?
There are two kinds of licenses.
1. With the purchase of a Numerical Recipes Code product, either a code
download or a Numerical Recipes Code CD-ROM, you automatically get a
non-expiring Numerical Recipes Personal, Single-User License that allows
one individual to use the code on any number of computers.
2. With the purchase of a Numerical Recipes institutional subscription,
you get a Numerical Recipes Institutional Subscriber License that allows
use of the code by any number of individuals on computers within a
subscriber's fixed IP address range, during the term of the subscription.
See complete license terms for both license types here.
Can I Distribute Numerical Recipes Routines with My Application
Bound "invisibly" into your executable .exe file? Generally yes. As
source code, linkable object libraries, or DLLs? Generally no. See more
details here.
--
Die Langversion ist hier:
http://www.nr.com/aboutNR3license.html
Grob, das Buch ist gut, um einen steilen Einstieg in die Materie zu
erhalten. Es ist nicht gut, um die Hintergr�nde zu beleuchten - es ist
eben ein Buch von Kochrezepten - und bei den Algorithmen haben sich die
Verfasser hier und da auch schon einen Schnitzer erlaubt. M.a.W,
gef�hrlich in den falschen H�nden, aber gut geeignet f�r einen
kritischen Leser mit einem Minimum an Hintergrund. Wenn man ein "ja wie
ging dass denn gleich..." braucht, genau das passende Buch.
Gr��e,
Thomas
Georg Bauhaus
> http://www.nr.com/aboutNR3license.html
Na, das kl�rt es doch. "All rights reserved", Urheberrecht
und allf�lliges Nutzungsrecht ist um diese Bestimmungen
erg�nzt. Es wird �berdeutlich, dass man sich weiters
schriftlich zu verst�ndigen hat.
Thomas Richter
/* snip */
> Krümelkackermodus aus.
Genau so. (-: Mathematiker eben. (Ich darf das sagen, ich gehöre zu der
Gesellschaft...)
>>> Weil ich auf die weise nur die Grundrechenarten brauche, um den Sinus zu
>>> berechnen. Und wie du unten sehen wirst, ist das durchaus nötig.
>> Achso? Also, da darf ich auch noch einige andere Algorithmen anbieten,
>> wie wäre es etwa mit dem CORDIC-Algorithmus?
>>
>
> Was ist das? Ach, moment, Wikipedia...
CORDIC ist nicht schnell, hat aber den Vorteil, eine garantierte
Genauigkeit liefern zu können. Ein entsprechendes Polynom zu designen
ist schwieriger. Darum fragte ich ja, was Deine Anwendung sein soll.
>> Erstmal gilt es zu überlegen, in welchem Sinne Du den Sinus überhaupt
>> approximieren willst. Da gibt es diverse Möglichkeiten, und je nach
>> Fragestellung gibt es diverse Antworten.
>>
>
> Und lass mich raten: "General Purpose" steht nicht zur Debatte?
Eigentlich nicht. Es ist immer eine Frage "wie genau oder wie schnell
hätten's denn gern." Die Frage ist also, was ist der Zweck der dietlibc?
Wenn der Zweck ist "so kurzen Code wie möglich", so wäre die Antwort,
die "fsin" Instruktion zu nutzen und nur dann eigenen Code zu verwenden,
wenn der Compiler kein "builtin" für sin() hat.
Wenn es garantierte Genauigkeit wäre, vermute ich ist CORDIC noch am
einfachsten zu handhaben.
Wenn es maximale Geschwindigkeit wäre, würde ich zu einer Lookup-Tabelle
greifen.
M.a.W, treffe eine wohlbegründete Entscheidung. Die üblichen
Polynomimplementierungen versuchen etwa einen Kompromiss zwischen
Genauigkeit und Codelänge.
> Die Transistoren, die den Sinus berechnen, sitzen an der falschen
> Stelle, und die einzige Verbindungsleitung ist ausreichend dünn.
(-:
>> Portabel auf was? Will sagen, jetzt, wo Du die Probleme der
>> Approximation kennst, sollte klar sein, dass die Polynome, die Du
>> brauchst, auch von der Rechengenauigkeit abhängen. Die garantiert Dir C
>> aber überhaupt nicht.
>>
>
> Nein, C garantiert mir die nicht, aber so weit wollte ich auch gar nicht
> gehen. Ich habe für meinen Dunstkreis weitere Garantien: Die ABI meines
> Systems und aller Systeme, zu denen ich kompatibel zu sein wünsche. Und
> diese ABIs definieren, dass es eine x87-FPU gibt.
Jetzt verstehe ich Dich immer weniger. x87-FPU *hat* ein FSIN. Also, wo
ist eigentlich das Problem? Die kürzeste Implementierung für eine
"diet"libc ist eine einzige Instruktion.
>>> Es kommt in eine Lib rein. Die ist eigentlich für General
>>> Purpose-Programmierung gedacht, d.h. sie sollte nicht ganz riesig und
>>> auch nicht ganz lahm werden. Der beste Kompromiss, sozusagen.
>> Hmm. Für wen? Für mich ist der beste Kompromiss, ein autoconf-Skript zu
>> benutzen welches mir die Verfügbarkeit einer built-in "sin" ermittelt,
>> von mir aus mit gcc-Syntax.
>
> gcc-Syntax?
>
> #include <stdio.h>
> int main()
> {
> printf("%d\n", __builtin_sin(0.75));
> return 0;
> }
>
> ist gcc-Syntax?
Naja, zum Beispiel. Man könnte auch direkt mittels __asm() die
Instruktion in den Code einbauen.
>> sin() wirft keine Exceptions nach C-Konvention. Vermutlich darf man
>> nicht einmal ein NAN zurückliefern, was sich eigentlich anböte. Mein
>> altes "Motorola"-Manual sagt dazu nur, dass für große Argumente jegliche
>> Genauigkeit verloren geht, was auch immer das heißen mag.
>>
>
> Exception im Sinne von Floating Point Exception kann sin() durchaus
> werfen, so wie jede andere Floating Point-Anweisung. Ich meine, was soll
> den sin(NAN) sein?
NAN.
> sin(NaNq) = NaNq und keine Invalid Exception,
> sin(NaNs) = NaNs und Invalid Exception sind doch gültige
> Implementierungen, oder?
Schon wahr, allerdings ist eben alles zwischen "groß" und "riesig" ein
valides Argument vom Sinus und sollte eben kein NAN liefert.
> Da liegt es doch nahe, eine Precision Exception
> auszulösen, wenn das Ergebnis unexakt ist, oder?
"Unexakt" ist das Ergebnis sowieso immer, außer das Argument ist null,
oder Du misst das Argument in "turns".
> Dummerweise bezieht sich 7.12.1 größtenteils auf das Ergebnis einer
> Operation. Allerdings ist ja 7.12.1p5 anwendbar, wenn man die Sache so
> interpretiert, dass mit dem kleinen Ergebnis die Größe in Relation zum
> Argument gemeint ist. Von einer Precision Exception steht dort nix, also
> ist nicht definiert, ob eine kommt oder nicht. Es steht allerdings da,
> dass, wenn math_errhandling & MATH_ERREXCEPT ungleich Null ist, es
> implementationsdefiniert ist, ob eine Underflow Exception ausgelöst
> wird.
Wir haben kein Overflow oder Underflow...
>>> Hmm... ja, stimmt. Eine Sache, die mir noch aufgefallen ist: Die meisten
>>> Funktionswerte des Sinus sind transzendente Zahlen, d.h. es gibt kein
>>> Polynom mit rationalen Koeffizienten, bei dem eine Nullstelle ein
>>> solcher Funktionswert wäre. Die Approximierung über ein Polynom hingegen
>>> sorgt dafür, dass alle Funktionswerte algebraisch sind.
>> Der Rechner kann keine transzendenten Zahlen. Nicht mal reelle.
>
> Ich hoffe, du weist, dass die reellen Zahlen eine Obermenge der
> transzendenten Zahlen sind. Genauer sind die reellen Zahlen die
> Vereinigung der algebraischen und transzendenten Zahlen.
Sowas passiert abends. Ja.
> Stimmt auch wieder. Mit ist nur eingefallen, dass das ein Beweis dafür
> ist, dass der Algorithmus außer an den ausgewählten Stellen nie exakt
> sein kann. Prinzipbedingt.
Was ich oben sagen wollte, richtig.
>> Achja, mir fällt noch eine Quelle ein, falls es Dich interessiert: Es
>> gibt von Motorla die FPSP (Floating Point Software Package) für den
>> 68060 und 68040. Beide Prozessoren machten den Sinus nicht mehr in
>> Hardware, sondern in Software über emulierte Befehle. Die Quellen habe
>> ich, und sie wurden von Mot damals herausgegeben. Falls Du also MC68K
>> Assembler lesen kannst, wäre das auch eine Möglichkeit, und *diese*
>> Quelle ist benutzbar. Vielleicht auch nur als Abschreckung. (-: Meine
>> Mail ist gültig, wenn Du den Algorithmus dort haben willst...
>>
>> Aber, wie ich schon sagte, das ganze ist kein Spaziergang, und der
>> Algorithmus ist ziemlich "convoluted".
>
> Bei den Aussichten lehne ich das Angebot dankend ab. Ich hab nämlich
> nicht den Hauch einer Ahnung vom Innenleben eines MC68K-Prozessors,
> geschweige denn seiner Assembler-Sprache.
Das Innenleben ist rasch erklärt und beschränkt sich auf acht
FPU-Register, acht Daten und acht Addressregister mit einem
Programmiermodel, welches um einiges kanonischer und einfacher als das
vom intel ist. Es wäre nur die einzige Quelle, die ich ohne das
Lizenzgeraffel anbieten kann.
Achja, nochwas zum Lesen und zum Nachschlagen von präzisen Werten, da Du
ja danach fragtest:
Abramowitz, Stegun: Handbook of Mathematical Functions(aus dem Gedächtnis)
Dort findest Du auch einige Approximationspolynome mit garantierten
Genauigkeiten (vermutlich nicht genug für double, aber sei's drum) sowie
Tabellen mit den auf N Stellen präzisen Werten.
Ich sehe gerade, dass sich das auch im Web findet, wobei ich mir das
allerdings nicht näher angesehen habe:
http://www.math.ucla.edu/~cbm/aands//
Grüße,
Thomas
Rainer Weikusat
> Thomas Richter schrieb:
>
>> http://www.nr.com/aboutNR3license.html
>
> und allf�lliges Nutzungsrecht ist um diese Bestimmungen
> erg�nzt.
Da steht etwas ueber die Lizenzsierung von 'Numerical Recipes Code'
und der Text macht deutlich, dass darunter etwas zu verstehen ist,
das man (in Form einer Bibliothek von Unterroutinen) zusammen mit
'anderem Code' zu einer Anwendung compilieren koenne:
You may, under this license, transfer precompiled, executable
applications incorporating the code to other, unlicensed,
persons,
Interessanter ist aber das hier:
If you are a licensed NR user, you can make any personal use you want
on a licensed screen, including translating NR to another computer
language. However, you can't redistribute the translations in any
manner, since they are "derivative works" and still subject to our
copyright.
Was jetzt genau eine 'Uebersetzung' und was eine unabhaengige
Reimplementierung des zugrundeliegenden Algorithmus' ist, darueber
waere mutmasslich vor Gericht zu streiten. Ich habe keine 'numerischen
Probleme', die ich irgendwie loesen muesste (jedenfalls nicht mehr,
seitdem ich keine Vorlesungen ueber 'Loesungen fuer numerische
Problem' habe und das ist schon ziemlich lange her), aber haette ich
welche, dann wuerde ich mir dieses Buch wegen oa Einschraenkung
jedenfalls sich nicht kaufen.
Georg Bauhaus
> Was jetzt genau eine 'Uebersetzung' und was eine unabhaengige
> Reimplementierung des zugrundeliegenden Algorithmus' ist, darueber
> waere mutmasslich vor Gericht zu streiten. Ich habe keine 'numerischen
> Probleme', die ich irgendwie loesen muesste (jedenfalls nicht mehr,
> seitdem ich keine Vorlesungen ueber 'Loesungen fuer numerische
> Problem' habe und das ist schon ziemlich lange her), aber haette ich
> welche, dann wuerde ich mir dieses Buch wegen oa Einschraenkung
> jedenfalls sich nicht kaufen.
Oder es k�nnte sein, dass kundige Patentanwaltschaften
sich �hnlichkeitsklagen ins Arsenal legen m�chten,
insofern "Unwissenheit vor Strafe nicht sch�tzt" und sich
auf diese Unwissenheit aufbauen l�sst. Es k�nnte also fatal
sein, aus dem Verzicht auf das NR-Buch zu schlie�en, dass
die selbstgemachten L�sungen vor den Anspr�chen fremder
Juristen und Betriebswirte sch�tzen k�nnen. Seien die Ideen
auch offensichtlich. Jedenfalls in den USA.
Ich glaube, dass diese Art der Besch�ftigung von Anw�lten
bei 1-Click durchgezogen worden ist.
Stefan Reuther
> Thomas Richter schrieb:
>>http://www.nr.com/aboutNR3license.html
>
> Na, das kl�rt es doch. "All rights reserved", Urheberrecht
> und allf�lliges Nutzungsrecht ist um diese Bestimmungen
> erg�nzt.
Das gilt alles nur dann, wenn Urheberrecht �berhaupt zutrifft. Das ist
nur dann der Fall, wenn ein ausreichend komplexer Programmteil direkt
kopiert wird. Ist er nur sehr simpel, fehlt die n�tige Sch�pfungsh�he.
Und die Reimplementation eines bekannten Algorithmus f�llt nicht unter
Urheberrecht, da Algorithmen nicht gesch�tzt sind.
Und gerade bei Numerical Recipes ist das Reimplementieren sehr sinnvoll.
Man kann die Algorithmen allesamst deutlich sch�ner und in deutlich
besserem C aufschreiben. Den albernen "wir machen die Arrays eins gr��er
bzw. dekrementieren unter Nutzung undefinierten Verhaltens Array-
Adressen, um wie in FORTRAN 1-basiert indizieren zu k�nnen"-Trick muss
man ja nun wirklich nicht nachbauen.
Stefan
Rainer Weikusat
Das 'die Amerikaner' auf den ausgesprochen gloriosen Gedanken
gekommen sind, ihr Aequivalent eines 'Patentamtes' in eine
Organisation umzuwandeln, die ihr Budget durch Verkauf bestimmter
Rechtstitel zu erwirtschaften hat, ist jedenfalls ihr Problem :->.
Es gibt auch einen 'einfachen' Workaround: Patentklagen werden nicht
angestrengt, weil jemandem langweilig ist, oder er boesartige
Absichten hegt, sondern weil dieser jemand sich dadurch einen
finanziellen Zugewinn erhofft. Sonst waere das ja ein
Zuschussgeschaeft. Dh jeder, bei dem diesbezueglich nichts zu holen
ist, kann sich einigermassen sicher fuehlen und im anderen Fall waere
eben ein bestimmter Teil der Gewinne als Notopfer Microsoft nach
Redmond abzufuehren. Laestig, aber machbar.
Markus Wichmann
in die Numerik, C, Compiler und Juristerei. Ich hab einfach mal kurz die
Version, die mein Compiler erzeugt hat, mit meinem Sinus verglichen und
festgestellt, dass trotz des ABI-Wrappers der x87-Code immer noch kürzer
ist als die Apporximierung mit SSE alleine. Und wie ich dank einer
funktionierenden Anwendung herausfand: SSE brächte mir nur dann
Vorteile, wenn ich für den Sinus eine Massenproduktion eröffnen würde.
Aber das tue ich ja nicht: Ein Funktionsaufruf zu sin() berechnet genau
einen Sinus. Und Ende.
Danke nochmals und tschö,
Georg Bauhaus
> Georg Bauhaus wrote:
>> Thomas Richter schrieb:
>>> http://www.nr.com/aboutNR3license.html
>>
>> Na, das kl�rt es doch. "All rights reserved", Urheberrecht
>> und allf�lliges Nutzungsrecht ist um diese Bestimmungen
>> erg�nzt.
>
> Das gilt alles nur dann, wenn Urheberrecht �berhaupt zutrifft.
Auch deswegen oben: "erg�nzt".
Ob die Nutzung der Ideen in X in den USA zu Klagen f�hren
kann oder nicht ist eine umfassendere Frage, als dass eine
Auskunft zum deutschen Urheberrecht zu X sie ersch�pfend
beantworten k�nnte. MP3.
BTW: Wer sich der Auffassung anschlie�t, dass Arbeit
an Algorithmen auf dem Weltmarkt belohnt, vulgo: verkauft
werden k�nnen sollte, wird der Frage nicht ausweichen
k�nnen, ob ihre Ver�ffentlichung dazu beitragen kann,
oder ob im Gegenteil damit faktisch nur ideeler Gewinn
oder suboptimaler Gewinn erzielt werden kann, denn wo
Reimplementierung kostenfrei m�glich ist, kann man u.U.
auf Anstand oder auf strategische H�flichkeit, ausgedr�ckt
in Bezahlung hoffen... Ist ein Algorithmus wie eine Formel
in der pharmazeutischen oder der Agrar-Industrie
(Pflanzenschutzmittel-Neuentwicklung, zugelassene
neue Kartoffelsorten, etc.)? Steckt in dem ver�ffentlichen
Arbeitsprodukt, den Algorithmen, genug Potential f�r
die Dienstleistung Support?
(In den Lobbies wird von beiden Seiten intensiv gearbeitet.
Ich hatte mal das Vergn�gen, die Sache: Paptentierbarkeit von
Algorithmen zwischen u.a. einem Vertreter des FFII (casual),
einem Sohn eines Br�sseler Anwalts, auch Anwalt (gr�ngrauer
Zweireiher) und einem Patentanwalt von Siemens (blauer Anzug)
verhandelt zu sehen.
Ein einziger "Begriffsstreit", in dem juristische und andere
Worte zu Instrumenten wurden, f�r den Versuch, die ohnehin
klar motivierte Einstellung zur Sache rhetorisch zu st�tzen.
Z.B. "Was ist Technizit�t?"
Kann/Soll sie f�r Algorithment ein anwendbares Kriterium
sein? (Siemens: "Ja!" FFII: "Nein!" Der Sohn: "Mein Vater ...!")
Welche Folgen hat so etwas? Der Patentanwalt empfahl, versucht
auch ein, zwei Patente zu erlangen, dann habt ihr was zu
verkaufen/verhandeln... Letzteres sehen manche Leute gern
als veritable Chance, sie hat so einen sch�nen Schein --
vor�bergehender Anerkennung. Oder erlaubt tats�chlich, auf Dauer
bei den Gro�en gern als Partnerunternehmen gesehen zu sein,
gekauft zu werden, etc.)
Google hat gelegentlich �ffentlich Ver�rgerung �ber den
rechtlichen Schutz der Privatsp�hre in Europa ge�u�ert.
Microsoft hat D�nemark mit Abzug (ergo: Ausfall von
Steuern und Fototerminen mit Politikern) gedroht, ebenfalls
wegen der MS nicht gef�lligen rechtlichen Situation.
Auch das sind Hinweise auf juristische Aktivit�t, die
das Urheberrecht nicht beantwortet, die es aber vielleicht
irgendwann doch ber�hren wird. Mal sehen, was das EU-
Parlament jetzt macht.
Rainer Weikusat
> On 12/8/09 6:32 PM, Stefan Reuther wrote:
[...]
> BTW: Wer sich der Auffassung anschlie锟絫, dass Arbeit
> an Algorithmen auf dem Weltmarkt belohnt, vulgo: verkauft
> werden k锟絥nen sollte, wird der Frage nicht ausweichen
> k锟絥nen, ob ihre Ver锟絝fentlichung dazu beitragen kann,
> oder ob im Gegenteil damit faktisch nur ideeler Gewinn
> oder suboptimaler Gewinn erzielt werden kann, denn wo
> Reimplementierung kostenfrei m锟絞lich ist, kann man u.U.
> auf Anstand oder auf strategische H锟絝lichkeit, ausgedr锟絚kt
> in Bezahlung hoffen... Ist ein Algorithmus wie eine Formel
> in der pharmazeutischen oder der Agrar-Industrie
> (Pflanzenschutzmittel-Neuentwicklung, zugelassene
> neue Kartoffelsorten, etc.)?
Wuerdest Du Dich der Ansicht anschliessen, dass jemand, der eine
praktische Abkuerzung eines Fussweges entdeckt hat (eine
Wegbeschreibung ist ein Algorithmus) auf dieser ein Zollhaeuschen
errichten duerfen sollte, damit hinfort niemand mehr diesen Weg ohne
finanzielle Kompensation benutzen kann, oder eher, dass er ja mal
versuchen koennte, sein Wissen um diesen neuen Algorithmus zu Geld zu
machen, dh stehen Leuten, die Stadtplaene veroeffentlichen, Abgaben
von Touristen zu, die diese benutzen, ohne einen Plan gekauft zu
haben? Und wie aendert sich die Ansicht, falls es ploetzlich nicht
mehr um offensichtlich bloedsinnige Beispiele geht, sondern um Richtig
Viel Geld[tm]?
Rainer Weikusat
> On 12/8/09 6:32 PM, Stefan Reuther wrote:
[...]
> BTW: Wer sich der Auffassung anschlie锟絫, dass Arbeit
> an Algorithmen auf dem Weltmarkt belohnt, vulgo: verkauft
> werden k锟絥nen sollte, wird der Frage nicht ausweichen
> k锟絥nen, ob ihre Ver锟絝fentlichung dazu beitragen kann,
> oder ob im Gegenteil damit faktisch nur ideeler Gewinn
> oder suboptimaler Gewinn erzielt werden kann, denn wo
> Reimplementierung kostenfrei m锟絞lich ist, kann man u.U.
> auf Anstand oder auf strategische H锟絝lichkeit, ausgedr锟絚kt
> in Bezahlung hoffen... Ist ein Algorithmus wie eine Formel
> in der pharmazeutischen oder der Agrar-Industrie
> (Pflanzenschutzmittel-Neuentwicklung, zugelassene
> neue Kartoffelsorten, etc.)?
Wuerdest Du Dich der Ansicht anschliessen, dass jemand, der eine
praktische Abkuerzung eines Fussweges entdeckt hat (eine
Wegbeschreibung ist ein Algorithmus) auf dieser ein Zollhaeuschen
errichten duerfen sollte, damit hinfort niemand mehr diesen Weg ohne
finanzielle Kompensation benutzen kann, oder eher, dass er ja mal
versuchen koennte, sein Wissen um diesen neuen Algorithmus zu Geld zu
machen, dh stehen Leuten, die Stadtplaene veroeffentlichen, Abgaben
von Touristen zu, die diese Wege benutzen, ohne einen Plan gekauft zu
Georg Bauhaus
> Wuerdest Du Dich der Ansicht anschliessen, dass jemand, der eine
> praktische Abkuerzung eines Fussweges entdeckt hat (eine
> Wegbeschreibung ist ein Algorithmus) auf dieser ein Zollhaeuschen
> errichten duerfen sollte,
An der Stelle, ob "Algorithmus" mit "Entdeckung eines
Wegs" gleichgestellt werden kann oder soll, und was eine
Entdeckung ist, scheiden sich die Geister. (Auff�llig
erscheint die Trennlinie zwischen stereotypen menschlichen
Z�gen, aber die Hinterb�hne sei nur am Rande erw�hnt.)
Was macht eine technische Erfindung aus und was macht
zwar technisches Ger�t patentierbar, Algorithmen allein in
der EU aber (noch) nicht?
Freilich, bei "abstrakten" Algorithmen ist es folgenreich,
wenn die jeweils ersten Entdecker sich eine Okkupationschance
sichern k�nnen: insofern sie die "�hnlichkeit" von Algorithmen
juristisch ausnutzen k�nn(t)en, blockiert die Ausnutzung
jede Konkurrenz, jede Weiterentwicklung, sofern diese �hnliche
Ideen verfolgt und keinen Okkupationstribut zahlen m�chte.
Im NR-Beispiel w�ren dann (NR-Patente nur mal angenommen) die
Umst�nde zu betrachten, unter denen die NR-Ideen neu und
vielleicht anders programmiert werden k�nnten. Eine juristische
Empfehlung lautet dann, neben den anders programmierten NR
weitere eigene Entwicklungen in die Waagschale werfen zu k�nnen,
wenn mit den "Schutzinhabern" verhandelt werden muss. Etwas zu
haben, ohne das die Entwickler der NR selbst nicht so recht
arbeiten k�nnen, zum Beispiel. Als eigene Bewaffnung. Dann
k�nnen die Juristen wieder Vertr�ge aushandeln. Dazu f�llt
mir "des Esels Schatten" ein.
Als n�chstes gelangt man argumentativ zum Anteil, den eine
okkupierte Idee an den Entwicklungen Dritter haben darf,
um sie blockieren zu k�nnen usw.
Charles Goodyear (der die Gr�ndung der gleichnamigen
Reifen-Firma nicht erlebt hat) wurde von einem etablierten
Industriellen ohne Umschweife seiner Gummi-Idee beraubt
(Stiefel in Goldrausch-Zeiten). Da er aber ein Patent auf
seine jahrelang betriebene Entwicklung durchfechten konnte,
musste der Industrielle Dieb immerhin $5000 zahlen. Nicht, dass
das Gooyears k�rglich lebende Familie endlich reich gemacht h�tte.
Ist ein "frei stehender" Algorithmus ebenso aufw�ndig, ebenso
detailliert, ebenso technisch, ebenso patentierbar wie so
ein Industrieverfahren zur Herstellung von Gummi?
D�rfen neben Bonner Gummib�rchen auch andere Weingummis
hergestellt werden, ohne dass eine Plagiatsklage im
Raum steht? (Nicht dass ich eine Klageabsicht unterstellen
m�chte, nur des bekannten Beispiels wegen.) Sagen wir,
1-Swish statt 1-Click... Wie steht es mit
der Sensortechnik von Apple zur Erkennung von Fingerbewegungen?
Sind da die Algorithmen allein patentiert?
Oder darf ich das Verfahren in Holz mit Linsentechnik
einfach so nachbauen?
Bei Algorithmen w�re es schon im allgemeinen Interesse,
denke ich, wenn "allgemein gangbare Wege" nicht
patentierbar w�ren. - Microsoft hat sich mal wieder irgend
so ein Verfahren f�r algorithmisch parallele
Verarbeitung von queues patentieren lassen, IIRC. Da muss man
also aufpassen und, wie du sagst, sp�testens dann mit der
Rechtsabteilung von Microsoft rechnen, wenn eine hinreichend
�hnliche Eigenentwicklung zuf�llig in deren Patentschutz hinein
geredet werden kann.
�brigens scheint es Algorithmen f�r die moderne Informations-
Verarbeitung zu geben, die so delikat ist, dass ein Produkt
(MS Exchange, die letzte) nur mit einem anderen Produkt (MS
Windows Server, die letzte) richtig genutzt werden kann.
Aus technischen Gr�nden... Aber ja.
Thomas Richter
> Rainer Weikusat schrieb:
>
>> Wuerdest Du Dich der Ansicht anschliessen, dass jemand, der eine
>> praktische Abkuerzung eines Fussweges entdeckt hat (eine
>> Wegbeschreibung ist ein Algorithmus) auf dieser ein Zollhaeuschen
>> errichten duerfen sollte,
>
/* snip */
> Im NR-Beispiel w�ren dann (NR-Patente nur mal angenommen) die
> Umst�nde zu betrachten, unter denen die NR-Ideen neu und
> vielleicht anders programmiert werden k�nnten. Eine juristische
> Empfehlung lautet dann, neben den anders programmierten NR
> weitere eigene Entwicklungen in die Waagschale werfen zu k�nnen,
> wenn mit den "Schutzinhabern" verhandelt werden muss. Etwas zu
> haben, ohne das die Entwickler der NR selbst nicht so recht
> arbeiten k�nnen, zum Beispiel. Als eigene Bewaffnung. Dann
> k�nnen die Juristen wieder Vertr�ge aushandeln.
/* snip */
Vielleicht w�re es durchaus hilfreich, mal in besagtes Buch
hineinzusehen, um sich ein Bild von der Situation zu machen, wie sie
sich wirklich darstellt. Bei Algorithmen verweisen die Autoren -
wissenschaftlich korrekt - auf die jeweiligen Quellen, dementsprechend
kann f�r diese Algorithmen kein Patentschutz von den Autoren
wahrgenommen werden, sondern lediglich ein Urheberrechtsschutz (bzw. im
amerikanischen Recht, welches ein Copyright und keinen Urheberschutz
kennt, ein Copyright) auf die jeweilige Implementierung. Man sollte eine
Implementierung also nicht zu dicht an der Quelle f�hren, sondern eben
darauf achten, nur das Funktionsprinzip und nicht den Code zu verwenden.
Dies ist nat�rlich keine rechtlich wasserdichte Strategie, aber heute
kann schwerlich eine Zeile Code schreiben ohne ein Trivialpatent zu
verletzen...
Gr��e,
Thomas
Georg Bauhaus
> Georg Bauhaus schrieb:
> /* snip */
>
>> Im NR-Beispiel w�ren dann (NR-Patente nur mal angenommen) die
>> Umst�nde zu betrachten, unter denen die NR-Ideen neu und
>> vielleicht anders programmiert werden k�nnten. Eine juristische
>> Empfehlung lautet dann, neben den anders programmierten NR
>> weitere eigene Entwicklungen in die Waagschale werfen zu k�nnen,
>> wenn mit den "Schutzinhabern" verhandelt werden muss. Etwas zu
>> haben, ohne das die Entwickler der NR selbst nicht so recht
>> arbeiten k�nnen, zum Beispiel. Als eigene Bewaffnung. Dann
>> k�nnen die Juristen wieder Vertr�ge aushandeln.
>
> /* snip */
>
> Vielleicht w�re es durchaus hilfreich, mal in besagtes Buch
> hineinzusehen, um sich ein Bild von der Situation zu machen, wie sie
> sich wirklich darstellt.
Die Patentannahme f�r NR habe ich nur hinein fingiert,
um eine sonst durchaus �bliche Strategie am so fingierten
Beispiel sichtbar zu machen. Tats�chlich finde ich es
w�nschenswert, wenn man vor der Benutzung von Algorithmen
klare Informationen �ber Bedingungen und Abh�ngigkeiten
bekommen kann. Und ich w�rde durchaus f�r nicht-triviale
Vorarbeit bezahlen, wo m�glich. (Wenn ich in einem anderen
Leben denn mal Gelegenheit h�tte, US-relevante Software zu
produzieren, die gegen Trivialpatente dann abgesichert
werden m�sste.)
Stefan Reuther
> Das 'die Amerikaner' auf den ausgesprochen gloriosen Gedanken
> gekommen sind, ihr Aequivalent eines 'Patentamtes' in eine
> Organisation umzuwandeln, die ihr Budget durch Verkauf bestimmter
> Rechtstitel zu erwirtschaften hat, ist jedenfalls ihr Problem :->.
> Es gibt auch einen 'einfachen' Workaround: Patentklagen werden nicht
> angestrengt, weil jemandem langweilig ist, oder er boesartige
> Absichten hegt, sondern weil dieser jemand sich dadurch einen
> finanziellen Zugewinn erhofft. Sonst waere das ja ein
> Zuschussgeschaeft. Dh jeder, bei dem diesbezueglich nichts zu holen
> ist, kann sich einigermassen sicher fuehlen
Ganz so einfach ist das nicht. Dem "Rechteinhaber" geht es gerne auch
darum, unliebsame Mitbewerber loszuwerden, egal, ob dort was zu holen
ist oder nicht und ob das ein Zuschussgesch�ft ist oder nicht.
Stefan
Rainer Weikusat
> Thomas Richter schrieb:
>> Georg Bauhaus schrieb:
>
>> /* snip */
>>
>>> Im NR-Beispiel w�ren dann (NR-Patente nur mal angenommen) die
>>> Umst�nde zu betrachten, unter denen die NR-Ideen neu und
>>> vielleicht anders programmiert werden k�nnten. Eine juristische
>>> Empfehlung lautet dann, neben den anders programmierten NR
>>> weitere eigene Entwicklungen in die Waagschale werfen zu k�nnen,
>>> wenn mit den "Schutzinhabern" verhandelt werden muss. Etwas zu
>>> haben, ohne das die Entwickler der NR selbst nicht so recht
>>> arbeiten k�nnen, zum Beispiel. Als eigene Bewaffnung. Dann
>>> k�nnen die Juristen wieder Vertr�ge aushandeln.
>>
>> /* snip */
>>
>> Vielleicht w�re es durchaus hilfreich, mal in besagtes Buch
>> hineinzusehen, um sich ein Bild von der Situation zu machen, wie sie
>> sich wirklich darstellt.
>
> Die Patentannahme f�r NR habe ich nur hinein fingiert,
> um eine sonst durchaus �bliche Strategie am so fingierten
> Beispiel sichtbar zu machen.
Es waere sehr freundlich, wenn Du Strategien, die Du, aus welchen
Gruenden auch immer, fuer 'ueblich' haelst, nicht auf Kosten von
Leuten demonstrieren wuerdest, die vollkommen andere Standpunkte
vertreten. Aber das war ja gerade diese uebliche Strategie in leichter
Abwandlung, nicht wahr?
> Tats�chlich finde ich es w�nschenswert, wenn man vor der Benutzung
> von Algorithmen klare Informationen �ber Bedingungen und
> Abh�ngigkeiten bekommen kann.
Tatsaechlich finde ich es wuenschenswert, wenn es keine
fussballgrossen Eiskoerner hagelt, auch dem wird wohl jeder
zustimmen, und es hat auch keinen Bezug zum Thema.
> Und ich w�rde durchaus f�r nicht-triviale Vorarbeit bezahlen, wo
> m�glich.
Das ist Dir ja auch vollkommen freigestellt. Genauso wie es jedem
freigestellt ist, nicht an einer Universitaet zu lehren, sondern als
Consultant zu arbeiten. Der wird nun allerdings, falls er die
Weitergabe von Wissen verhindern moechte, Patentschutz brauchen, denn
sonst wird man ihm ins Gesicht lachen und einen Konkurrenten
beauftragen. Ausgesprochen laestige Sache das. Man muesste mal was
dagegen tun, dass jeder alles selbst erlernen und nach Massgabe seiner
Kenntnisse und Faehigkeiten kommerziell oder nichtkommerziell zu
verwerten ...
Rainer Weikusat
> Thomas Richter schrieb:
>> Georg Bauhaus schrieb:
>
>> /* snip */
>>
>>> Im NR-Beispiel w�ren dann (NR-Patente nur mal angenommen) die
>>> Umst�nde zu betrachten, unter denen die NR-Ideen neu und
>>> vielleicht anders programmiert werden k�nnten. Eine juristische
>>> Empfehlung lautet dann, neben den anders programmierten NR
>>> weitere eigene Entwicklungen in die Waagschale werfen zu k�nnen,
>>> wenn mit den "Schutzinhabern" verhandelt werden muss. Etwas zu
>>> haben, ohne das die Entwickler der NR selbst nicht so recht
>>> arbeiten k�nnen, zum Beispiel. Als eigene Bewaffnung. Dann
>>> k�nnen die Juristen wieder Vertr�ge aushandeln.
>>
>> /* snip */
>>
>> Vielleicht w�re es durchaus hilfreich, mal in besagtes Buch
>> hineinzusehen, um sich ein Bild von der Situation zu machen, wie sie
>> sich wirklich darstellt.
>
> Die Patentannahme f�r NR habe ich nur hinein fingiert,
> um eine sonst durchaus �bliche Strategie am so fingierten
> Beispiel sichtbar zu machen.
Es waere sehr freundlich, wenn Du Strategien, die Du, aus welchen
Gruenden auch immer, fuer 'ueblich' haelst, nicht auf Kosten von
Leuten demonstrieren wuerdest, die vollkommen andere Standpunkte
vertreten. Aber das war ja gerade diese uebliche Strategie in leichter
Abwandlung, nicht wahr?
> Tats�chlich finde ich es w�nschenswert, wenn man vor der Benutzung
> von Algorithmen klare Informationen �ber Bedingungen und
> Abh�ngigkeiten bekommen kann.
Tatsaechlich finde ich es wuenschenswert, wenn es keine
fussballgrossen Eiskoerner hagelt, auch dem wird wohl jeder
zustimmen, und es hat auch keinen Bezug zum Thema.
> Und ich w�rde durchaus f�r nicht-triviale Vorarbeit bezahlen, wo
> m�glich.
Das ist Dir ja auch vollkommen freigestellt. Genauso wie es jedem
freigestellt ist, nicht an einer Universitaet zu lehren, sondern als
Consultant zu arbeiten. Der wird nun allerdings, falls er die
Weitergabe von Wissen verhindern moechte, Patentschutz brauchen, denn
sonst wird man ihm ins Gesicht lachen und einen Konkurrenten
beauftragen. Ausgesprochen laestige Sache das. Man muesste mal was
dagegen tun, dass jeder alles selbst erlernen und nach Massgabe seiner
Kenntnisse und Faehigkeiten kommerziell oder nichtkommerziell verwerten
darf ...
Georg Bauhaus
> Georg Bauhaus <rm.dash...@futureapps.de> writes:
>> Thomas Richter schrieb:
>>> Georg Bauhaus schrieb:
>>> /* snip */
>>>
>>>> Im NR-Beispiel w�ren dann (NR-Patente nur mal angenommen) die
>>>> Umst�nde zu betrachten, unter denen die NR-Ideen neu und
>>> /* snip */
>>>
>>> Vielleicht w�re es durchaus hilfreich, mal in besagtes Buch
>>> hineinzusehen, um sich ein Bild von der Situation zu machen, wie sie
>>> sich wirklich darstellt.
>> Die Patentannahme f�r NR habe ich nur hinein fingiert,
>> um eine sonst durchaus �bliche Strategie am so fingierten
>> Beispiel sichtbar zu machen.
>
> Es waere sehr freundlich, wenn Du Strategien, die Du, aus welchen
> Gruenden auch immer, fuer 'ueblich' haelst, nicht auf Kosten von
> Leuten demonstrieren wuerdest, die vollkommen andere Standpunkte
> vertreten.
Stimmt, obwohl ich nicht glaube, dass den Autoren durch unser
Gerede -- in dem Konjunktive wohl nicht hinreichen, mit
rhetorischer Verrenkung h�tte ein alternatives Beispiel
eingef�hrt werden m�ssen -- Kosten entstehen.
Vielleicht schauen sich ein paar Menschen die
Nutzungsm�glichkeiten genauer an.
> Tatsaechlich finde ich es wuenschenswert, wenn es keine
> fussballgrossen Eiskoerner hagelt, auch dem wird wohl jeder
> zustimmen, und es hat auch keinen Bezug zum Thema.
In der Tat, ja. Das hat keinen Bezug zum Thema.
> Man muesste mal was
> dagegen tun, dass jeder alles selbst erlernen und nach Massgabe seiner
> Kenntnisse und Faehigkeiten kommerziell oder nichtkommerziell verwerten
> darf ...
Naja, erstens darf ein Mensch Gelerntes nicht �berall kommerziell oder
nichtkommerziell verwerten, das finde ich erw�hnenswert.
Oder jedes Wissen, zum Beispiel algorithmisches, f�r jeden Zweck.
Dir wird nicht entgangen sein, dass f�r die erlaubte Nutzung von
Gelerntem manchmal Vertr�ge erforderlich sein k�nnen. (Oder das
Recht zur Nutzung. Nicht jeder Rechtsstaat hat zu jeder Zeit erlaubt,
Kommunikation nach Belieben zu verschl�sseln, sofern man gelernt
hat, mit welchem Algorithmus das geht. Sinn hin, politischer
Sachverstand her...)
Zweitens sind es gerade die, die du leider nicht nennst, die
"Man m�sste mal verhindern" tats�chlich gesagt haben und sagen, die
das Umfeld mitschaffen, in dem man sich um die Nutzungsbedingungen von
ver�ffentlichtem und nicht ver�ffentlichtem Wissen k�mmern wird.
"The licenses for most software and other practical works are designed
to take away your freedom to share and change the works." (GPL)
Drittens tat man was f�r und gegen die Nutzung von Gelerntem und
wird es immer tun:
f�r Kunst, Kunstf�lschung, Innovatives, Plagiate, u.v.m. muss
man wissen und k�nnen. Fragen der Bewertung und der Zul�ssigkeit
der einzelnen Nutzung werden mal mit daf�r, mal mit dagegen
beantwortet.
Bei den NR hat man Antworten und Anschrift.