A320773 - OEIS

%I #38 Jan 25 2020 02:32:20

%S 2,3,5,6,8,10,11,12,15,17,18,20,24,26,27,30,35,37,38,39,40,42,48,50,

%T 51,56,63,65,66,68,72,80,82,83,84,87,90,99,101,102,104,105,110,120,

%U 122,123,132,143,145,146,147,148,150,152,156,168,170,171,182,195,197,198,200

%N Numbers (excluding squares) whose square root has a continued fraction with a period < 3.

%C The Heron sequence of every number a(n) has the following relationship: numerator(h(k))^2 - a(n)*denominator(h(k))^2 = 1 for k > 1.

%C The Heron sequence of every number a(n) has the following relationship with the continued fraction f(s) convergent to sqrt(a(n)): h(k) = f(2^k-1).

%C From _Gerhard Kirchner_, Jan 17 2020: (Start)

%C Numbers k = m^2 + r with m > 0 and 0 < r <= 2m such that r is a divisor of 2m.

%C Continued fraction: k = [m; 2m/r, 2m, 2m/r, 2m, ...].

%C The number of terms that are between m^2 and (m+1)^2  is equal to the number of divisors of 2m, which is A099777(m).

%C Proof see link. The Maxima code below demonstrates the divisor property. Note that there is no divisor of 2m between m and 2m.

%C (End)

%H Gerhard Kirchner, <a href="/A320773/a320773_1.pdf">Divisor property of a(n)</a>

%e The continued fraction of sqrt(6) = 2 + 1/(2 + 1/(4 + 1/(2 + 1/(4 + 1/(2 + 1/(4 + ...)))))) = [2; 2, 4, 2, 4, 2, 4, ...] has repeating portion (2, 4) with period 2, so 6 is a term.

%p Digits:=40: nr:=0:

%p for z from 2 to 200 do

%p   test:=true: c:=sqrt(z):

%p   if (c=floor(c)) the test:=false: end if:

%p   while (test=true) do

%p     b[0]:=floor(c):

%p     r[0]:=c:

%p     for k from 1 to 2 do

%p       r[k]:=evalf(1/(r[k-1]-b[k-1])):

%p       b[k]:=floor(r[k]):

%p     end do:

%p     if (b[1]=2*b[0]) or (b[2]=2*b[0]) then nr:=nr+1: a[nr]:=z: printf("%4d",z): end if:

%p     test:=false:

%p   end do:

%p end do:

%t Select[Range[200], !IntegerQ[Sqrt[#]] && Length@ContinuedFraction[Sqrt[#]][[-1]]<3 &] (* _Amiram Eldar_, Nov 01 2018 *)

%o (Maxima)

%o block([n: 2, m: 0, r: 0, k: 0, kmax: 10,v: ""],

%o   while k<kmax do

%o     (m: floor(sqrt(n)), r: n-m^2,

%o     if mod(2*m,r)=0 then (k: k+1, print(n)),

%o     if r= m then n: n+m else (if r= 2*m then n: n+2 else n: n+1)));

%o /* _Gerhard Kirchner_, Jan 17 2020 */

%K nonn

%O 1,1

%A _Paul Weisenhorn_, Oct 21 2018

%E Edited by _Jon E. Schoenfield_, Oct 19 2019