Discussion:
die Rechengenauigkeit bei floating-point kann zum Problem werden ...
(zu alt für eine Antwort)
Walter H.
2018-03-14 21:22:53 UTC
Permalink
Raw Message
warum man manchmal grübelt ...
und sich wundert, daß der Schluß, weil es bei kleinen
Werten passt, auch bei großen Werten passe, trügerisch sein kann ...

void printfact( unsigned long num )
{
unsigned long i;
long double factlog = 0.0;

for ( i = 2; i <= num; i++ )
factlog += log10l( (long double) i );

printf( "fact(%lu) \367 %8.6lfE+%010lu.\n\n", num,
powl( 1.0, factlog - floorl( factlog ) ),
(unsigned long) floorl( factlog ) );
}


so klappt es dann ...

void printfact( unsigned long num )
{
unsigned long i;
unsigned long factloghlp = 0;
long double factlog = 0.0;

for ( i = 2; i <= num; i++ )
{
while ( factlog > 10.0 ) // kann auch 'if' sein ...
{
factlog -= 10.0;
factloghlp += 10;
}
factlog += log10l( (long double) i );
}

printf( "fact(%lu) \367 %8.6lfE+%010lu.\n\n", num,
powl( 1.0, factlog - floorl( factlog ) ),
(unsigned long) floorl( factlog ) + factloghlp );
}
Thomas Jahns
2018-03-20 12:24:10 UTC
Permalink
Raw Message
Post by Walter H.
warum man manchmal grübelt ...
und sich wundert, daß der Schluß, weil es bei kleinen
Werten passt, auch bei großen Werten passe, trügerisch sein kann ...
void printfact( unsigned long num )
{
    unsigned long i;
    long double factlog = 0.0;
    for ( i = 2; i <= num; i++ )
        factlog += log10l( (long double) i );
    printf( "fact(%lu) \367 %8.6lfE+%010lu.\n\n", num,
        powl( 1.0, factlog - floorl( factlog ) ),
        (unsigned long) floorl( factlog ) );
}
ich glaube, Du willst Dich in Kahan-Summation[1] einlesen. Ansonsten gilt der
Hinweis, dass es am einfachsten ist, numerische Verfahren von Leuten
implementieren zu lassen, die sich damit auskennen.

Thomas

[1] <https://en.wikipedia.org/wiki/Kahan_summation_algorithm>
Thomas Koenig
2018-03-25 12:01:25 UTC
Permalink
Raw Message
Post by Thomas Jahns
Post by Walter H.
warum man manchmal grübelt ...
und sich wundert, daß der Schluß, weil es bei kleinen
Werten passt, auch bei großen Werten passe, trügerisch sein kann ...
void printfact( unsigned long num )
{
    unsigned long i;
    long double factlog = 0.0;
    for ( i = 2; i <= num; i++ )
        factlog += log10l( (long double) i );
    printf( "fact(%lu) \367 %8.6lfE+%010lu.\n\n", num,
        powl( 1.0, factlog - floorl( factlog ) ),
        (unsigned long) floorl( factlog ) );
}
ich glaube, Du willst Dich in Kahan-Summation[1] einlesen. Ansonsten gilt der
Hinweis, dass es am einfachsten ist, numerische Verfahren von Leuten
implementieren zu lassen, die sich damit auskennen.
Man kann aber auch dafür sorgen, dass man sich selber auskennt :-)
Goldberg ist da Pflichtlektüre: "What Every Computer Scientist
Should Know Abut Floating-Point Arithmetic", findet man problemlos
im Netz.

Peter J. Holzer
2018-03-24 10:51:58 UTC
Permalink
Raw Message
Post by Walter H.
warum man manchmal grübelt ...
und sich wundert, daß der Schluß, weil es bei kleinen
Werten passt, auch bei großen Werten passe, trügerisch sein kann ...
void printfact( unsigned long num )
{
unsigned long i;
long double factlog = 0.0;
for ( i = 2; i <= num; i++ )
factlog += log10l( (long double) i );
printf( "fact(%lu) \367 %8.6lfE+%010lu.\n\n", num,
powl( 1.0, factlog - floorl( factlog ) ),
(unsigned long) floorl( factlog ) );
}
so klappt es dann ...
void printfact( unsigned long num )
{
unsigned long i;
unsigned long factloghlp = 0;
long double factlog = 0.0;
for ( i = 2; i <= num; i++ )
{
while ( factlog > 10.0 ) // kann auch 'if' sein ...
{
factlog -= 10.0;
factloghlp += 10;
}
factlog += log10l( (long double) i );
}
printf( "fact(%lu) \367 %8.6lfE+%010lu.\n\n", num,
powl( 1.0, factlog - floorl( factlog ) ),
(unsigned long) floorl( factlog ) + factloghlp );
}
Aber auch nur zufällig. Du implementierst hier im Wesentlichen
Fixed-Point-Arithmetik (mit ca. 61 Bit hinter dem Komma (bei x87
Arithmetik).

Das Grundproblem, dass der Logarithmus einer natürlichen Zahl eine
irrationale Zahl ist, die in endlich viel Speicher nicht exakt
dargestellt werden kann und daher gerundet werden muss, bleibt bestehen.
Du addierst dann gerundete Zahlen auf, und der Rundungsfehler
akkumuliert sich. Und ganz zum Schluss verwendest Du eine Funktion, die
Sprungstellen einfügt. Dass da manchmal der wahre Wert knapp an einer
Seite der Sprungstelle liegt, der berechnete aber auf der anderen Seite,
ist unvermeidlich.

hp

PS: Bei welchen Werten liefern die zwei Funktionen eigentlich wirklich
unterschiedliche Resultate? Ich bin jetzt bei 100000 angelangt (das Ding
ist langsam) und bis noch sind beide gleich.
--
_ | Peter J. Holzer | Fluch der elektronischen Textverarbeitung:
|_|_) | | Man feilt solange an seinen Text um, bis
| | | ***@hjp.at | die Satzbestandteile des Satzes nicht mehr
__/ | http://www.hjp.at/ | zusammenpaßt. -- Ralph Babel
Loading...