Pathological floating point problems: Difference between revisions
Content added Content deleted
Thundergnat (talk | contribs) m (→{{header|C}}: Remove vanity tags) |
|||
Line 1,876: | Line 1,876: | ||
3rd: Rump's example: f(77617.0, 33096.0) = -0.827396059946821 |
3rd: Rump's example: f(77617.0, 33096.0) = -0.827396059946821 |
||
</pre> |
|||
=={{header|Phix}}== |
|||
Standard maths with the IEEE-754 hardware floats fails predictably (that code left in as comments), so roll out the bigatoms.<br> |
|||
Task1: needs at least 120 decimal places to avoid serious divergance, and 196 to be accurate to 78 digits<br> |
|||
Task2: needs at least 41 decimal places (and 42 passed as the first argument to ba_euler)<br> |
|||
Task3: apparently only needs just 15 decimal places, then again ba_scale() defines the minimum, and it may (as in is permitted to) use more. |
|||
<lang Phix>include builtins\bigatom.e |
|||
puts(1,"Task 1\n") |
|||
constant {fns,fmts} = columnize({{3,1},{4,6},{5,6},{6,6},{7,6},{8,7},{20,16},{30,24},{50,40},{100,78}}) |
|||
{} = ba_scale(196) |
|||
sequence v = {2,-4} |
|||
for n=3 to 100 do |
|||
-- v = append(v,111 - 1130/v[n-1] + 3000/(v[n-1]*v[n-2])) |
|||
v = append(v,ba_add(ba_sub(111,ba_divide(1130,v[n-1])),ba_divide(3000,ba_multiply(v[n-1],v[n-2])))) |
|||
if n<9 or find(n,{20,30,50,100}) then |
|||
-- printf(1,"n = %-3d %20.16f\n", {n, v[n]}) |
|||
string fmt = sprintf("%%.%dB",fmts[find(n,fns)]) |
|||
printf(1,"n = %-3d %s\n", {n, ba_sprintf(fmt,v[n])}) |
|||
end if |
|||
end for |
|||
puts(1,"\nTask 2\n") |
|||
--atom balance = exp(1)-1 |
|||
--for i=1 to 25 do balance = balance*i-1 end for |
|||
--printf(1,"\nTask 2\nBalance after 25 years: $%12.10f", balance) |
|||
{} = ba_scale(41) |
|||
bigatom balance = ba_sub(ba_euler(42,true),1) |
|||
for i=1 to 25 do balance = ba_sub(ba_multiply(balance,i),1) end for |
|||
ba_printf(1,"Balance after 25 years: $%.16B\n\n", balance) |
|||
puts(1,"Task 3\n") |
|||
{} = ba_scale(15) -- fine! |
|||
integer a = 77617, |
|||
b = 33096 |
|||
--atom pa2 = power(a,2), |
|||
-- pb2a211 = 11*pa2*power(b,2), |
|||
-- pb4121 = 121*power(b,4), |
|||
-- pb6 = power(b,6), |
|||
-- pb855 = 5.5*power(b,8), |
|||
-- f_ab = 333.75 * pb6 + pa2 * (pb2a211 - pb6 - pb4121 - 2) + pb855 + a/(2*b) |
|||
--printf(1,"f(%d, %d) = %.15f\n\n", {a, b, f_ab}) |
|||
bigatom pa2 = ba_power(a,2), |
|||
pb2a211 = ba_multiply(11,ba_multiply(pa2,ba_power(b,2))), |
|||
pb4121 = ba_multiply(121,ba_power(b,4)), |
|||
pb6 = ba_power(b,6), |
|||
pa2mid = ba_multiply(pa2,ba_sub(ba_sub(ba_sub(pb2a211,pb6),pb4121),2)), |
|||
pb633375 = ba_multiply(333.75,pb6), |
|||
pb855 = ba_multiply(5.5,ba_power(b,8)), |
|||
f_ab = ba_add(ba_add(ba_add(pb633375,pa2mid),pb855),ba_divide(a,ba_multiply(2,b))) |
|||
printf(1,"f(%d, %d) = %s\n", {a, b, ba_sprintf("%.15B",f_ab)})</lang> |
|||
Aside: obviously you don't ''have'' to use lots of temps like that, but I find that style often makes things much easier to debug. |
|||
{{out}} |
|||
<pre> |
|||
Task 1 |
|||
n = 3 18.5 |
|||
n = 4 9.378378 |
|||
n = 5 7.801153 |
|||
n = 6 7.154414 |
|||
n = 7 6.806785 |
|||
n = 8 6.5926328 |
|||
n = 20 6.0435521101892689 |
|||
n = 30 6.006786093031205758530554 |
|||
n = 50 6.0001758466271871889456140207471954695237 |
|||
n = 100 6.000000019319477929104086803403585715024350675436952458072592750856521767230266 |
|||
Task 2 |
|||
Balance after 25 years: $0.0399387296732302 |
|||
Task 3 |
|||
f(77617, 33096) = -0.827396059946821 |
|||
</pre> |
</pre> |
||