computing the infinite sum
15 次查看(过去 30 天)
显示 更早的评论

I am trying to solve this problem but I get an error when I use the code below. Anyone have any ideas?
s=(1.1:.2:5.9);
syms k;
symsum(1/k^s,k,1,Inf);
0 个评论
采纳的回答
Star Strider
2016-4-24
The bsxfun function makes this straightforward:
k = 1:2E+6;
s = 1.1:0.2:5.9;
Zser = bsxfun(@power, k, -s'); % Generate Matrix
Precision = diff(Zser(:,end-1:end),[],2); % Check Precision
Z = sum(Zser, 2); % Calculate Sum
figure(1)
plot(s, Z)
grid
title('\zeta(s) = \Sigma k^{-s}')
xlabel('s')
ylabel('\zeta(s)')
This meets the criteria for accuracy, as the ‘Precision’ vector demonstrates.
2 个评论
Matt J
2016-4-24
It looks to me that "Precision" is merely calculating the difference in the final 2 terms of the series. I think they are expected to show that whole truncated tail of the series falls below the 14-digit precision threshold.
更多回答(1 个)
Roger Stafford
2016-4-24
It is my opinion that this problem would require something other than adding the given terms out sufficiently far for 14-digit accuracy. The following inequality is obviously true:
int(1/x^s,n,n+1) > 1/(n+1)^s > int(1/x^s,n+1,n+2)
which means that the infinite sum
1/(n+1)^s + 1/(n+2)^s + 1/(n+3)^s +...
lies between int(1/x^s,n,inf) and int(1/x^s,n+1,inf). These two quantities can be solved as (s-1)/n^(s-1) and (s-1)/(n+1)^(s-1). With the s = 1.1 value, to ensure that the sum which is left out,
1/(n+1)^s + 1/(n+2)^s + 1/(n+3)^s +...
is not more than 10^(-14) it we would have to have
(s-1)/(n+1)^(s-1) < 10^(-14),
or substituting s = 1.1
.1/(n+1)^.1 < 10^(-14)
10^13 < (n+1)^.1
10^130 < n+1
The number 10^130 is horribly large and the consequent cumulative round off error would be enormous, though the universe would suffer a heat death before the computation was completed. However, if we get smart, we can use the above integral to complete the summation. Choose n sufficiently large that the difference between the two integrals is less than 10^(-14):
(s-1)*(1/n^(s-1)-1/(n+1)^(s-1)) < 10^(-14)
1/n^.1-1/(n+1)^.1 < 10^(-13)
The answer to that is n = 81114515936, which above 81 billion. The plan would be then, for s = 1.1, to add the first n = 81114515936 terms, but in reverse order to minimize accumulated rounding error, and then add on the value .1/81114515936^.1 for the integral.
This is a far cry from stopping after a couple million terms. However, in any case, the cumulative rounding error from 81 billion additions is likely to greatly exceed 10^(-14) which indicates in my mind, a hastily conceived homework problem.
2 个评论
Roger Stafford
2016-4-25
If I had been a little smarter in my answer, the infinite trailing series
1/(n+1)^s + 1/(n+2)^s + 1/(n+3)^s + ...
can better be approximated (in true 'trapz' spirit) by
int(1/x^s,n+1/2,n+3/2)
so that it would not be necessary to have 81 billion odd additions but some much smaller number (which I am too lazy at the moment to determine.)
Roger Stafford
2016-4-25
Just for the record, I erred in the integration above. The integral, int(1/x^s,n,inf), is 1/(s-1)*1/n^(s-1), and not (s-1)/n^(s-1). However, this just makes matters worse. One would have to do 10^150 terms to achieve 14-digit accuracy.
另请参阅
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!