Числовая ошибка Matlab и как получить правильный ответ

Я каким-то образом получаю следующее выражение в Matlab (R2014a на W7, 64b)

1/1034591578977116160000*prod(1:19)*(29576428208904825-17729494921579950*k - 20479697577410832*k^2 + 13867226524449248*k^3 - 836937224095392*k^4 - 869194297188672*k^5 + 163710902234880*k^6 + 2589894827520*k^7 - 2476912838400*k^8 + 144848704000*k^9)

где k изначально является символьной переменной. Затем я устанавливаю k=10 и получаю результат 370371188237528, используя выходной формат LONGG. Но если я применю это же выражение к Mathematica (заменив prod(1:19) на 19!), я получу 370371188237525, что, как мне кажется, является правильным ответом. Кажется, это ошибка округления, описанная много раз на этом сайте (правильно ли это?). Как избежать этого с символическим набором инструментов Matlab или без него?


matlab rounding-error symbolic-computation
person user155002    schedule 28.11.2014    source источник

Ответы (1)


arrow_upward
1
arrow_downward

Существует класс High Precision Float. После добавления его к пути вы можете просто написать k = hpf(10); перед оценкой вашего выражения. Результат будет

370371188237525.0106290979251578332118698122510380699168308638036

С Symbolic Math Toolbox я бы написал

syms k
expr = 1/1034591578977116160000*prod(1:19)*(29576428208904825-17729494921579950*k - 20479697577410832*k^2 + 13867226524449248*k^3 - 836937224095392*k^4 - 869194297188672*k^5 + 163710902234880*k^6 + 2589894827520*k^7 - 2476912838400*k^8 + 144848704000*k^9);
subs(expr, k, 10);

который оценивается как 3150006955960150124/8505 = 370371188237525.

person Jommy    schedule 28.11.2014
comment
О боже! Я не знал, что это существует! Я бы хотел +10. - person Ander Biguri; 28.11.2014
comment
Пока я не перезапустил Matlab, я не мог повторить результат (у меня было 2962969505900221/8 от subs(expr, k, 10)). Интересно, почему.. Но теперь все в порядке. Спасибо. - person user155002; 28.11.2014