obliczenie błędu ułamków n/m w double

0

Obliczam dowolny ułamek, albo i nawet wpisuję stałe typu: 7/15, lub 1/360, itp.
i potem używam to w obliczeniach i wychodzą błędy.

Przykładowo obliczamy wyrażenie typu:

s = -1/3 a + 4/3 b;

i tu wiadomo że dla a = b powinno być: s = a dokładnie,
no ale z obliczeń tak nie będzie;

double c1 = -1/3, c2 = 4/3;

c = c1 + c2,
i wcale nie będzie całe 1, lecz coś chyba mniej nieco.

No i chodzi mi właśnie o wyliczenie tych błędów reprezentacji,
tj. ile jest: e1 = c1 + 1/3 i e2 = c2 - 4/3 ?

Tego nie można obliczyć wprost, czyli tak:

double e1 = c1 - 1/3;

ponieważ to będzie równoważne z tym: c1 - c1 = 0.

Ale przecież błąd jest tu około taki: 0.0000000000000000333333333333333333,
czyli spokojnie wejdzie do double w postaci: 3.333333333333333e-17
bo zasięg jest tu aż do 1e-300, chyba.

No ale jak to obliczyć?

1/360 - fl(1/360) = ?
i inne takie - w zasadzie dowolne liczby.
fl(1e44/5465333333333333333) - 1e44/5465333333333333333 = ?

0

z C w nagłówku float.h zdefiniowane jest makro (stała) DBL_EPSILON reprezentująca różnice między jedynką a najmniejszą liczbą większą od jedynki, której reprezentacje można zapisać w double.
liczba*DBL_EPSILON daje ci dokładność reprezentacji danej liczby.
W momencie wykonywania jakiś operacji możesz w odpowiedni sposób sumować te błędy.

Pamiętaj jednak, że kompilator same obliczenia może przeprowadzać na typie innym niż double, ale na typie zmiennoprzecinkowym, który działa najszybciej w danej architekturze i dokonać konwersji, gdy wynik zapisywany jest w jakiejś zmiennej (co prowadzi do poprawy dokładności). Patrz link powyżej opis stałej FLT_EVAL_METHOD.

0

Jak chcesz to sumować?

A to DBL_EPSILON to po prostu 2^-53, czyli trochę powyżej 1e-16.

To chyba trzeba obliczać metodami numerycznymi,
czyli zrobić w zasadzie dzielenie o precyzji 2 x double.

x2 = (x_hi, x_lo); i to lo trzyma kolejne 53 bity,
po dodaniu wprost: x = x_hi + x_lo, otrzymamy tu x = x_lo, dokładnie - tracimy te 53 bity z lo.

Algorytm dzielenia Newtona podwaja liczbę bitów w każdym kroku, czyli wystarczy jeden krok.

0

albo podwajasz precyzję i tak estymujesz błędy (lepsza ale bardzie kosztowna metoda), albo robisz rachunek błędu. Z rachunek błędów wychodzą takie wzory (zakładając, że obliczenia robione są na double FLT_EVAL_METHO!=2):

s = -1.0/3.0 * a + 4.0/3.0 * b;
sErrMax = (fabs(1.0/3.0 * a)+fabs(4.0/3.0 * b))*2*DBL_EPSILON; // błąd maksymalny
sErrAvrg = sqrt(2.0/9.0*a*a + 32.0/9.0*b*b)*DBL_EPSILON; // błąd średni
0
MarekR22 napisał(a):

Z rachunek błędów wychodzą takie wzory (zakładając, że obliczenia robione są na double FLT_EVAL_METHO!=2):

s = -1.0/3.0 * a + 4.0/3.0 * b;
sErrMax = (fabs(1.0/3.0 * a)+fabs(4.0/3.0 * b))*2*DBL_EPSILON; // błąd maksymalny
sErrAvrg = sqrt(2.0/9.0*a*a + 32.0/9.0*b*b)*DBL_EPSILON; // błąd średni

To są nieprawidłowe wzory na błąd.

Wystarczy wstawić dokładnie reprezentowane liczby i wtedy tu otrzymasz err > 0, a faktycznie będzie zero.

Ty obliczasz błąd operacji zmiennoprzecinkowych, czyli: er = |a op b - fl(a op b)| <= |a op b|*eps;

Np. błąd obliczania wyrażenia: y = ab + cd

fl(y) = (ab(1+eps) + cd(1+eps)((1+eps) <= ab(1+2eps) + cd(1+2eps) = ab + cd + (ab+cd)2eps;
każda operacja osobno: dwa mnożenia i jedno odejmowanie.

Jak widać średni błąd jest tu równy około: er = (ab-cd)2eps, czyli to co tam napisałeś.

A ten drugi wzór jest chyba zupełnie nieprawidłowy -
błąd obliczania wyrażenia: sqrt( (ab)^2^ + (cd)^2^ ) byłby znacznie większy.

dodanie znaczników <code> oraz `` - fp

0

Ogólny wzór na błąd średni:
\Delta f(x) = \sqrt{\sum_i \left(\frac{\partial}{\partial x_i}f(x)\cdot\Delta x_i\right)^2}
Wzór na błąd maksymalny:
\Delta f(x) = \sum_i \left|\frac{\partial}{\partial x_i}f(x)\cdot\Delta x_i\right|
Nie ma siły średni wychodzi mniejszy niż maksymalny.

1 użytkowników online, w tym zalogowanych: 0, gości: 1