Witam,


Napisałem drobny program w c z użyciem gsl, oblicz silnię(100).
Wynik wydaje się dobry ale ten błąd obliczeń ? Jak to interpretować ?
#include <stdio.h>
#include <gsl/gsl_errno.h>
#include <gsl/gsl_sf.h>


// gcc -I/usr/include/gsl/ -L/usr/local/lib/ -lgsl -lgslcblas -lm f.c
// gcc -I/usr/include/gsl/ -L/usr/lib/ -lgsl -lgslcblas -lm f.c
// factorial
//
http://www.gnu.org/software/gsl/manual/html_node/Spe
cial-Function-Usage.html#Special-Function-Usage
// The error-handling function
//
http://www.gnu.org/software/gsl/manual/html_node/Spe
cial-Functions-Examples.html#Special-Functions-Examp
les
/*




*/

int main(void)
{
unsigned int n = 100;
int status;
gsl_sf_result result;

status = gsl_sf_fact_e (n, &result);


printf("status = %s\n", gsl_strerror(status));
printf("factorial( %2d ) = %f\n", n,result.val);
printf("+/- % .18f\n",result.err);

return 0;
}


./a.out
status = success
factorial( 100 ) =
9332621544394415096564670479595388257840097037318409
8831012889540582227238570431295066113089288327277825
8496640065242705545359762897193828521818658959597240
32.000000
+/-
4144516527479786869998377376409676501474209436767641
6602380879768122591161803228174716423865792481641279
576393802186138583881092629521428905984.000000000000
000000



TIA

Adam

do góry

On 24/07/2013 21:04, Adam Majewski wrote:
> Witam,
>
>
> Napisałem drobny program w c z użyciem gsl, oblicz silnię(100).
> Wynik wydaje się dobry ale ten błąd obliczeń ? Jak to interpretować ?

Precyzja wyniku z tym błedem obliczeń wynosi dokładnie 51 bitów. Double
ma precyzję 52 bitów. Nie wiem dokładnie dlaczego akurat jeden bit ci
przycięło, ale sądzę, że to nie przypadek.

do góry

W dniu 2013-07-24 22:04, Adam Majewski pisze:
> Witam,
>
>
> Napisałem drobny program w c z użyciem gsl, oblicz silnię(100).
> Wynik wydaje się dobry ale ten błąd obliczeń ? Jak to interpretować ?
> #include <stdio.h> #include <gsl/gsl_errno.h> #include
> <gsl/gsl_sf.h>
>
>
> // gcc -I/usr/include/gsl/ -L/usr/local/lib/ -lgsl -lgslcblas -lm f.c
> // gcc -I/usr/include/gsl/ -L/usr/lib/ -lgsl -lgslcblas -lm f.c //
> factorial //
> http://www.gnu.org/software/gsl/manual/html_node/Spe
cial-Function-Usage.html#Special-Function-Usage
>
>
>
> // The error-handling function //
> http://www.gnu.org/software/gsl/manual/html_node/Spe
cial-Functions-Examples.html#Special-Functions-Examp
les
>
>
>
> /*
>
>
> */
>
> int main(void) { unsigned int n = 100; int status; gsl_sf_result
> result;
>
> status = gsl_sf_fact_e (n, &result);
>
>
> printf("status = %s\n", gsl_strerror(status)); printf("factorial(
> %2d ) = %f\n", n,result.val); printf("+/- % .18f\n",result.err);
>
> return 0; }
>
>
> ./a.out status = success factorial( 100 ) =
> 9332621544394415096564670479595388257840097037318409
8831012889540582227238570431295066113089288327277825
8496640065242705545359762897193828521818658959597240
32.000000
>
>
>
> +/-
> 4144516527479786869998377376409676501474209436767641
6602380879768122591161803228174716423865792481641279
576393802186138583881092629521428905984.000000000000
000000
>
>
>


gsl_sf_result to para _double_'i. Siłą rzeczy dokładność jest
na poziomie 10^-15..10^-16. Nie bardzo wiem, jak udało Ci się
wypisać te wszystkie znaki, ale większość z nich nie ma znaczenia.

Wynik powinien raczej być zapisany jako

9.332621544394415*10^(157)+-4.14452*10^(142)

Wartość / błąd to 4.44089*10^(-16)
Dość rozsądnie:)


Jeśli chcesz liczby naturalnej, musisz użyć zmienych naturalnych
dowolnej precyzji. Nie widzę, by gsl coś takiego miał.

Jest oczywiście superkombajn gmp (pod windowsem podobno lepiej
_podobno_ użyć forka MPIR).

W c++ masz bardzo wygodną lekką ttmath i nawet w boscie coś
jest (i możę działać jako frontend dla gmp:))

pzdr
bartekltg

do góry