Project Euler Problem 368 modified harmonic series
Problem 36822 January 2012
The harmonic series 1 + 1 /2
+ 1 /3
+ 1 /4
+ ... is well known to be divergent.
If we however omit from this series every term where the denominator has a 9 in it, the series remarkably enough converges to approximately 22.9206766193.
This modified harmonic series is called the Kempner series.
Let us now consider another modified harmonic series by omitting from the harmonic series every term where the denominator has 3 or more equal consecutive digits. One can verify that out of the first 1200 terms of the harmonic series, only 20 terms will be omitted.
These 20 omitted terms are:
1/111
, 1/222
, 1 /333
, 1 /444
, 1 /555
, 1 /666
, 1 /777
, 1 /888
, 1 /999
, 1 /1000
, 1 /1110
, 1/1111
, 1 /1112
, 1 /1113
, 1 /1114
, 1 /1115
, 1 /1116
, 1 /1117
, 1 /1118
and1/1119
.
This series converges as well.
Find the value the series converges to.
Give your answer rounded to 10 digits behind the decimal point. 这个应该可以用PolyGamma函数的展开式帮助计算:
http://bbs.emath.ac.cn/viewthread.php?tid=2390&page=2&fromuid=20#pid28839 网上找到一个代码，谁能解释下
function [ s,error ] = Problem42
% Thomas Schmelzer, November 2008
% computes the limit of the Kempner series \sum_{n=1}^\infty 1/n
% where n does not contain the substring "42" --- a problem that has been
% proposed by Folkmar Bornemann. For details see
% Schmelzer, Thomas; Robert Baillie (June朖uly 2008). "Summing a Curious,
% Slowly Convergent Series". American Mathematical Monthly 115 (6): 525?40
% exact solution: 228.44630 41592 30813 25414 80861 26250 58957 81629 ...
%% given the recurrence matrix T, compute f,A and B.
T = ;
= find(T); f = zeros(2,2,10); A = zeros(2,2);
for s = 1:length(I)
f(T(I(s),J(s)),I(s),J(s))=1;
end
for s = 0:9
A = A + f(:,:,s+1)/10;
end
I = eye(size(A)); B = inv(I-A)-I;
%% define the S_1 and S_2 explicitly for integers with up to 5 digits
S = cell(5,2); S{1,1} = ; S{1,2} = 4;
for i = 2:size(S,1)
for m = 0:1:9 % possible digits to attend
= find(f(:,:,m+1));
for s = 1:length(J)
S{i,J(s)} = ;
end
end
end
%% define the numerical parameters
K = 20;% extrapolation
P = 30;% power cutoff
Ignore = -22;
Psi = zeros(K,2,P);
% compute the sums Psi_{i} explicitly for up to 5 digits
for i = 1:size(S,1)
for k = 1:P
Psi(i,1,k) = sum(S{i,1}.^(-k));
Psi(i,2,k) = sum(S{i,2}.^(-k));
end
end
warning('off','MATLAB:log:logOfZero');
for i = size(S,1)+1:1:K
for k = 1:P
if (Psi(i-1,1,k)>0)
for m = 0:9
= find(f(:,:,m+1));
for s = 1:length(J)
for w = 0:P-k
if log10(abs(aCoeff(k,w,m))) + log10(abs(Psi(i-1,L(s),k+w))) > Ignore;
Psi(i,J(s),k) = Psi(i,J(s),k)+aCoeff(k,w,m)*Psi(i-1,L(s),k+w);
end
end
end
end
end
end
end
warning('on','MATLAB:log:logOfZero');
% Extrapolation
s = sum(B*Psi(K,:,1)') + sum(sum(Psi(:,:,1)));
error = 228.44630415923081325414 - s;
function [ a ] = aCoeff( k,w,m )
if (m+w==0)
mw = 1;
else
mw = m^w;
end
a = 10^(-k-w)*mw*(-1)^w*nchoosek(k+w-1,w);
还有这个，但不知怎么改造成368的答案
页:
[1]