|
%I #57 Nov 29 2023 19:00:08
%S 1,5,14,26,47,67,100,136,177,221,286,338,415,491,568,652,761,861,990,
%T 1098,1219,1351,1520,1652,1813,1985,2162,2334,2555,2727,2960,3172,
%U 3397,3641,3878,4098,4383,4659,4944,5204,5537,5805,6158,6466,6775,7123,7524,7848
%N Total number of steps of Euclid's GCD algorithm to calculate gcd(x,y) for all pairs x,y in the range 1 <= x,y <= n.
%C The number of steps to calculate gcd(x,y) is A107435(x,y) for x<=y, and is the length of the continued fraction expansion of x/y (including 0 for the integer part).
%C A018806(n) <= a(n) <= A064951(n).
%H Project Euler, <a href="https://projecteuler.net/problem=433">Problem 433: Steps in Euclid's Algorithm</a>
%F a(n) = a(n-1) + n + A049826(n) * 2 and a(1) = 1.
%F a(n) = a(n-1) + n + (Sum_{j=1..n-1} A107435(n,j)) * 2 and a(1) = 1.
%F a(n) = Sum_{x=1..n} Sum_{y=1..n} A107435(x,y). - _Michel Marcus_, Nov 28 2023
%p g:= proc(x, y) option remember;
%p `if`(y=0, 0, 1+g(y, irem(x, y)))
%p end:
%p a:= proc(n) option remember; `if`(n=0, 0,
%p a(n-1)+n+2*add(g(n, j), j=1..n-1))
%p end:
%p seq(a(n), n=1..100); # _Alois P. Heinz_, Nov 27 2023
%o (Python)
%o from functools import cache
%o A107435 = lambda x,y: 0 if y == 0 else 1 + A107435(y, x % y)
%o A049826 = lambda n: sum(A107435(n,j) for j in range(1, n))
%o @cache
%o def a(n):
%o # Code after _Alois P. Heinz_
%o if n == 0: return 0
%o return a(n-1) + n + A049826(n) * 2
%o print([a(n) for n in range(1,49)])
%o (PARI)
%o A107435(n, k) = if (k == 0, 0, 1 + A107435(k, n % k));
%o a(n) = sum(x=1, n, sum(y=1, n, A107435(x,y)));
%o print(vector(49,n,a(n)));
%Y Cf. A018806, A049826, A064951, A107435.
%K nonn
%O 1,2
%A _Darío Clavijo_, Nov 26 2023
|
