Problem mit sqrt

Philipp schrieb:
> Hallo ich führe folgende Berechnung aus:
>

1

> #include <math.h>

2

> int main (void)

3

> {

4

>     signed int x = sqrt((long)((2*3.14159/24)*2*10000000000)/ 375);

5

>     return 0;

6

> }

7

>

> Das Ergebnis von x ist 2393
> sollte laut Taschenrechner aber 3736 sein.
> Wo ist der Fehler???

Da entstehen überall Rundungsfehler. Divisionen solltest du immer 
zuletzt durchführen, genauso wie Umwandlungen in Gleitkommaformate.
Außerdem muss der Wertebereich stimmen.
Probiere mal folgendes:

 sqrt((long long)((2*2*10000000000)*3.14159)/ (24*375))

10000000000, also 10*1e9 ist zu groß für ein 32bit-long. unsigned long 
geht nur bis ca. 4*1e9. Da man leider an deinen Typdefinitionen nicht 
erkennen kann, in welchem Datenformat du rechnest, und ob du 
vorzeichenbehaftet arbeitest, kann man sich nicht sicher sein, ob nicht 
auch dabei schon ein Fehler auftritt.
Verwende am besten Typen aus stdint.h, z.B. uint64_t. Damit siehst du 
auch gleich, welchen Wertebereich du hast. Ein int64_t würde auch nicht 
reichen, da der positive Wertebereich nur 32Bit umfasst. Du brauchst 
hier aber mind. ~37 Bit.

Um die Ungenauigkeit bei der Umrechnung in Gleitkommaformate zu 
vermeiden, wäre aber folgendes höchst sinnvoll:

 sqrt((long long)((2*2*10000000000)*314159)/ (100000*24*375))

Da du mit einem beschränkten Wert für Pi rechnest, ist die 
Genauigkeitsanforderung auf diese 5 für dich signifikanten 
Nachkommastellen beschränkt. Da ist Festkomma besser als Gleitkomma.
Und dann entsprechend die Vereinfachung:

 sqrt((long long)((2*2*100000)*314159)/ (24*375))

Wenn du jetzt noch ein wenig intelligent kürzt (z.b. 2*2/24), dann 
kannst du zusätzliche Bits sparen.

Manuel schrieb:
> Probier mal:
>
>

1

#include <math.h>

2

> int main (void)

3

> {

4

>     signed int x = sqrt((long)((2.0*3.14159/24.0)*2*10000000000.0)/ 

5

> 375.0);

6

>     return 0;

7

> }

8

>

Das ist Unsinn.
Die Ungenauigkeit durch die Verwendung sehr großer (10e9) und 
gleichzeitig sehr kleiner (Pi) Gleitkommazahlen macht die Sache nur noch 
schlimmer. Damit treten massive Rundungsfehler auf.

Übrigens:

ich bekomme schon hiermit das korrekte Ergebnis:

int32_t x = sqrt((uint64_t)((2*3.14159/24)*2*10000000000)/ 375);

Also liegt es vor allem an deinen Typecasts. long (32bit) ist einfach zu 
klein für die Gleichung.

Philipp schrieb:
> imonbln schrieb:
>> Aber da du uns keine Information zu deinen Mirkocontroller, Kompiler
>> oder ähnlichen gibts
>
> Jup das habe ich im Eifer des Gefechts vergessen
> Ich nutze einen STM32F100
> Compiler ist ARM GCC 4.7

Dann binde stdint.h ein und nimm die Zeile:

int32_t x = sqrt((uint64_t)((2*3.14159/24)*2*10000000000)/ 375);
Auf einem F103 läuft die genauso und bringt das richtige Ergebnis.

>
> Manuel schrieb:
>> Probier mal:
> Führt leider zum gleichen Ergebnis

Du solltest dich wirklich intensiv damit beschäftigen, wie Zahlen 
dargestellt werden. Dann erkennst du, dass bei bestimmten 
Rechenoperationen in Verbindung mit bestimmten Datentypen automatisch 
Rundungsfehler entstehen. Ein Beispiel dafür ist die Division mit 
integer-Datentypen. Ein anderes ist die Verwendung von Gleitkommatypen, 
dabei tritt mit jeder Rechenoperation ein Rundungsfehler auf. Die Frage 
ist eben nur, ob man diesen Rundungsfehler akzeptieren kann.

proier mal:

 #include <math.h>
 int main (void)
 {
     signed int x = sqrt((long)((3.14159*40000000)/ 9));
     return 0;
 }

Da ist das gleiche nur mit vorher manuell gekürztem Bruch. Das sieht so 
konstruiert aus. Ist das eine Übungsaufgabe? Oder warum castest du 
((2*3.14159/24)*2*10000000000) nach long und teilst danach durch 375?

Philipp schrieb:
> Hallo Frank Bär
>
> VIELEN DANK FÜR DIE SCHNELLE UND SEHR AUSFÜHRLICHE HILFE!!!
>
> Frank Bär schrieb:
>> ich bekomme schon hiermit das korrekte Ergebnis:
>>
>> int32_t x = sqrt((uint64_t)((2*3.14159/24)*2*10000000000)/ 375);
>
> Ja damit komme ich nun auch zum gewümschten Ergebnis.
>
> Das mit den Datentypen und Darstellungen muss ich mir echt nochmal zu
> Gemüte führen, da hast du recht!!!

auch ein wenig Mathe schadet nicht.

sqrt(a*b) <==> sqrt(a) * sqrt(b)

d.h. du kannst dein Ergebnis schon mal enorm pimpen, indem du den 
Ausdruck

   sqrt( 2 * 3141592654 / 24 * 10000000000 ...

zerlegst in

   sqrt( 2 * 3.141592654 ) * sqrt( 10000000000 )

was du damit erreicht hast: Du hast vermieden, dass du kleine Zahlen (2 
* 3.141592654) mit sehr grossen Zahlen (10000000000) verrechnen musst, 
was dir deine Genauigkeit in den Keller treibt.


Formeln einfach naiv vom Papier in eine Anweisung einer 
Programmiersprache zu übertragen, ohne den Faktor "Was passiert mit 
meiner Genaugikeit bzw. kann ich mit dem Datentyp die notwendige 
Auflösung überhaupt erreichen" zu berücksichtigen, geht oft ins Auge. 
Floating Point ist NICHT die Lösung aller Probleme, die man einfach naiv 
einsetzen kann und alles ist gut.

Karl Heinz Buchegger schrieb:
> was du damit erreicht hast: Du hast vermieden, dass du kleine Zahlen (2
> * 3.141592654) mit sehr grossen Zahlen (10000000000) verrechnen musst,
> was dir deine Genauigkeit in den Keller treibt.

Weil ich das jetzt schon zum zweiten Mal lese: das ist in diesem Kontext 
Unsinn. Das Problem tritt nur bei der Addition derartiger Zahlen auf. 
Hier haben wir aber eine Multiplikation


XL

Irgendwo sollte doch auch mal eine Variable ins Spiel kommen, oder?
Sonst kann man noch weiter kürzen:

1

signed int x = 3736;

:-)

Karl Heinz Buchegger schrieb:
> auch ein wenig Mathe schadet nicht.

Als Freund von ganzzahligen Berechnungen würde ich noch über die 
geforderte Genauigkeit verschärft nachdenken. Manchmal reicht für PI 
auch 22/7 ...

Axel Schwenke schrieb:
> Karl Heinz Buchegger schrieb:
>> was du damit erreicht hast: Du hast vermieden, dass du kleine Zahlen (2
>> * 3.141592654) mit sehr grossen Zahlen (10000000000) verrechnen musst,
>> was dir deine Genauigkeit in den Keller treibt.
>
> Weil ich das jetzt schon zum zweiten Mal lese: das ist in diesem Kontext
> Unsinn. Das Problem tritt nur bei der Addition derartiger Zahlen auf.
> Hier haben wir aber eine *Multiplikation*
>
>
> XL

Willst du jetzt wirklich sagen, dass Gleitkommamultiplikation keine 
Rundungsfehler erzeugt? Klar entstehen dabei Fehler. Vor allem bedingt 
dadurch, dass hohe Genauigkeit und große Zahlen in Kombination schlicht 
nicht funktioniert. Im Beispiel gibt es fünf signifikante 
Nachkommastellen.
Ein Double-Precision-Float stellt für Pi mit 5 Nachkommastellen aber 15 
Nachkommastellen dar.
Jetzt multipliziere ich dieses gerundete Pi mit 10e9 und zaubere nicht 
nur meine 5 signifikanten Nachkommastellen, sondern auch noch 5 weitere 
in das am Ende ganzzahlige Ergebnis.
Gleitkommazahlen und Rechenoperationen vertragen sich einfach nicht, 
wenn man bestimmte Ansprüche hat. Hohe Genauigkeit mag gehen, dann aber 
nur mit kleinen Zahlen. Große Zahlen gehen, dann aber mit geringer 
Genauigkeit. Große Zahlen mit vielen signifikanten Stellen funktioniert 
nicht, egal ob Addition oder Multiplikation.

Frank Bär schrieb:
> Axel Schwenke schrieb:
>> Karl Heinz Buchegger schrieb:
>>> was du damit erreicht hast: Du hast vermieden, dass du kleine Zahlen (2
>>> * 3.141592654) mit sehr grossen Zahlen (10000000000) verrechnen musst,
>>> was dir deine Genauigkeit in den Keller treibt.
>>
>> Weil ich das jetzt schon zum zweiten Mal lese: das ist in diesem Kontext
>> Unsinn. Das Problem tritt nur bei der Addition derartiger Zahlen auf.
>> Hier haben wir aber eine *Multiplikation*
>
> Willst du jetzt wirklich sagen, dass Gleitkommamultiplikation keine
> Rundungsfehler erzeugt?

Das habe ich natürlich nicht gesagt. Und das ist auch gar nicht der 
Punkt.

Der Punkt ist, daß der bei der Gleitkommamultiplikation (dito Division) 
entstehende Fehler ausschließlich von der Mantisse abhängt. Der Exponent 
(vulgo: der Größenunterschied zwischen den Operanden) beeinflußt den 
Fehler überhaupt nicht.
Zumindest so lange nicht, wie kein Exponenten-Überlauf eintritt.

Ganz anders ist die Situation bei der Gleitkommaaddition (dito: 
Subtraktion). Wenn da die Operanden verschiedene Exponenten haben, dann 
müssen die Mantissen vor der Addition gegeneinander verschoben werden. 
Dabei fallen signifikante Bits des kleineren Operanden auf der rechten 
Seite aus dem Mantissen-"Fenster" heraus. Das ist der zusätzliche 
Fehler den man sich bei Gleitkomma-Rechnungen mit stark 
unterschiedlichen Zahlen einhandelt.


XL

Axel Schwenke schrieb:
> Der Punkt ist, daß der bei der Gleitkommamultiplikation (dito Division)
> entstehende Fehler ausschließlich von der Mantisse abhängt. Der Exponent
> (vulgo: der Größenunterschied zwischen den Operanden) beeinflußt den
> Fehler überhaupt nicht.

Hierzu vielleicht mal ein Beispiel:

1

float x = 1.5;          /* 3/2 */

2

float y = 0.1875;       /* 3/16 */

3

float z = 0.0029296875; /* 3/1024 */

Da der Nenner aller drei Brüche jeweils eine Zweierpotenz ist, besteht 
die Mantisse in Gleitkommadarstellung jeweils aus den gleichen Bits wie 
der Zähler. Für x, y und z ist die Mantisse 11000... (2 signifikante 
Bits). [Einschub: das hidden bit in der Mantisse ignoriere ich für diese 
Betrachtung; man stelle sich einfach vor die Mantisse wäre 1 Bit länger 
und enthalte das hidden bit]

Weil die Mantissen rechts so viele Nullen haben, ist hier die 
Multiplikation sogar fehlerfrei. Sowohl x*y als auch x*z haben die 
Mantisse 1001... (4 signifikantes Bits). Der Unterschied in der 
Skalierung wird durch die Addition der Exponenten aufgefangen. Und diese 
Operation ist im Normalfall (ohne Überlauf) exakt.

Die gleiche Betrachtung gilt auch, wenn die Operanden Mantissen mit mehr 
signifikanten Bits haben. Ein Fehler tritt bei der Multiplikation dann 
auf, wenn das Produkt der beiden Mantissen mehr Bits hat als die 
Mantisse lang ist. Der Fehler ist aber unabhängig von den Exponenten.


XL

Wozu ist überhaupt der long-Cast da? Soll da irgendetwas abgerundet
werden? Wenn nicht:

1

    signed int x = sqrt(2 * 3.14159 / 24 * 2 * 10000000000 / 375);

oder lesbarer:

1

    signed int x = sqrt(2 * 3.14159 / 24 * 2 * 1e10 / 375);

-> 3736