No.343 Fractional Sequences
問題．
任意の正の整数 \(k\) に対して，分数 \(y_i/x_i\) で表される数 \(a_i\) の有限列を次のように定義する：
\(a_1 = 1/k\)
\(a_i = (x_{i-1}+1)/(y_i-1)\quad i>1\)
(但し最も簡単な形に約分する．)
\(a_i\) がなんらかの整数になったとき，数列の終わりとする．(つまり，\(y_i=1\) のとき．)
\(f(k)=n\) と定義する．
例えば，\(k=20\) のとき
\(1/20\to 2/19\to 3/18=1/6\to2/5\to3/4\to4/3\to5/2\to6/1=6\)
なので，\(f(20)=6\).
また，\(f(1)=1, f(2)=2, f(3)=1, \sum_{1\leq k\leq 100}f(k^3)=118937\).
では，\(\sum_{1\leq k\leq 2\times 10^6}f(k^3)\) を求めよ．
考え方とプログラム例．(Python3)
まず，実際にすべての数について規則通りに数列を構成していくのは大変なので，この関数 \(f\) について少し考えてみます．
分数の分母を 1 減らして分子を 1 増やすという操作において，分母と分子の和は変わりません．
約分されるのはどのようなときかを考えると，もちろん，分母と分子が 1 でない共通の約数をもつときですが，その約数は分母と分子の和の約数でもあります．
つまり，最小の素因数に分子がなったときに約分されることが繰り返されていきます．\(k+1\) の最大の素因数が \(p\) だとすると，\(1/(p-1)\) になったあとは \((p-1)/1=p-1\) まで約分がおこることはなく，終了します．よって，\(f(k)\) は \(k+1\) の最大の素因数 -1 を表すことが分かります．
(実際，\(k=20\) のとき \(20+1=3\cdot 7\) の最大素因数は 7 なので \(f(20)=7-1=6\). )
そこで，\(1\leq k\leq 2\times 10^6\) について，(\(k^3+1\) の最大の素因数-1) の和を取ればよいと分かります．
\(k^3+1\) は因数分解できて，
\begin{align*}
k^3+1 = (k+1)(k^2-k+1)
\end{align*}
となることを利用します．
\(k^2-k+1<k^2\) から，\(k^2-k+1\) は \(k\) 以上の素因数をたかだか 1 個しか持たない
ことがわかるので，はじめに \(k+1\) 以下の素数のリストを作っておいて小さい素数から順に割っていきます．
今回はPython3 で351秒もかかってしまいました．1 分以内で解ける方法があるはずですが．
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 |
import time, math start = time.clock() N = 2 * 10**6 num = [0] * (N+1) num[0] = 1 prime = [] for i in range(N+1): if num[i] == 0: prime.append(i+1) num[i] = i+1 n = 2*(i+1) while n <= N+1: num[n-1] = i+1 n += (i+1) if i+1 > math.sqrt(N+1): for (j, k) in enumerate(num[i:]): if k == 0: prime.append(j+i+1) num[j+i] = j+i+1 n = 2*(j+i+1) while n <= N+1: num[n-1] = j+i+1 n += (j+i+1) break Sum = 0 for i in range(1, N+1): K = i * i - i + 1 if K <= N: Sum = Sum + max(num[i], num[K-1]) - 1 else: for j in prime: if j * j > K: break while K % j == 0: K = K / j if K == 1: Sum = Sum + max(num[i], j) - 1 break if K == 1: break if K != 1: Sum = Sum + max(num[i], K) - 1 print(Sum) print(str(time.clock()-start)+"sec") |
答えはこの下にあります．反転してください．
269533451410884183
