OFFSET
1,2
COMMENTS
Note that a(n)=n if n is a square. The square root of a squarefree integer n has a continued fraction of the form [e(0);[e(1),...,e(p)]] with e(p)=2e(0) and e(i)=e(p-i) for 0 < i < p, see reference. The symmetric part [e(1),...,e(p-1)] of the continued fraction [m;[e(1),...,e(p-1), 2m]] will be called the pattern of n. 2 has the empty pattern (sqrt(2)=[1,[2]]), 3 has the pattern [1] (sqrt(3)=[1,[1,2]]) and so on. In this sense, the description of the sequence can be simplified as "Least number greater than n with the same pattern".
It can be can proved (see link) that integers with the same pattern are terms of a quadratic sequence.
An ambiguity has to be fixed: sqrt(2)=[1,[2]] = [1,[2,2]] = [1,[2,2,2]] and so on. We define that the shortest pattern is correct, here it is empty. Comment on the third subsequence (2),6,12,... below: The second term 6 has the pattern [2], but the first term 2 in brackets has the "wrong" pattern, after fixing the ambiguity.
REFERENCES
Kenneth H. Rosen, Elementary number theory and its applications, Addison-Wesley, 3rd ed. 1993, page 428.
LINKS
Chai Wah Wu, Table of n, a(n) for n = 1..138
Gerhard Kirchner, Continued fractions with the same pattern
EXAMPLE
1) p=1: f(1)=2, f(2)=a(2)=5, f(3)=a(5)=10, f(4)=a(10)=17,..
sqrt(2)=[1,[2]], sqrt(5)=[2,[4]], sqrt(10)=[3,[6]], sqrt(17)=[4,[8]],..
2) p=2: f(1)=3, f(2)=a(3)=8, f(3)=a(8)=15, f(4)=a(15)=24,..
sqrt(3)=[1,[1,2]], sqrt(8)=[2,[1,4]], sqrt(15)=[3,[1,6]], sqrt(24)=[4,[1,8]],..
3) p=3: f(1)=41, f(2)=a(41)=130, f(3)=a(130)=269,..
sqrt(41)=[6,[2,2,12]], sqrt(130)=[11,[2,2,121]], sqrt(269)=[16,[2,2,256]],..
4) p=4: f(1)=33, f(2)=a(33)=60, f(3)=a(60)=95,..
sqrt(33)=[5,[1,2,1,10]], sqrt(60)=[7,[1,2,1,49]], sqrt(95)=[9,[1,2,1,81]],..
Several subsequences f(k) with f(k+1)=a(f(k)).
k>1 if first term in brackets, k>0 otherwise.
First terms Period Formula Example
1) 2,5,10,17 1 A002522(k)=k^2+1 1
2) 3,8,15,24 2 A005563(k)=(k+1)^2-1 2
3)(2),6,12 2 A002378(k)=k*(k+1)
4) 7,32,75 4 A013656(k)=k*(9*k-2)
5) 11,40,87 2 A147296(k)=k*(9*k+2)
6) 13,74,185 5 A154357(k)=25*k^2-14*k+2
7) (3),14,33 4 A033991(k)=k*(4*k-1) 4
8) (5),18,39 2 A007742(k)=k*(4*k+1)
9) 21,112,275 6 A157265(k)=36*k^2-17*k+2
10)23,96,219 4 A154376(k)=25*k^2-2*k
11)27,104,231 2 A154377(k)=25*k^2+2*k
12)28,299,858 4 A156711(k)=144*k^2-161*k+45
13)29,338,985 5 A156640(k)=169*k^2+140*k+29
14)(8),34,78 4 A154516(k)=9*k^2-k
15)(10),38,84 2 A154517(k)=9*k^2+k
16)(2),41,130 3 A154355(k)=25*k^2-36*k+13 3
17)47,192,435 4 A157362(k)=49*k^2-2*k
PROG
(Maxima)
block([nmax: 100],
/*saves the first nmax terms in the current directory*/
algebraic: true, local(coeff), showtime: true,
fl: openw(sconcat("terms", nmax, ".txt")),
coeff(w, m):=
block(a: m, p: 0, s: w, vv:[],
while a<2*m do
(p: p+1, s: ratsimp(1/(s-floor(s))), a: floor(s),
if a<2*m then vv: append(vv, [a])),
j: floor((p-1)/2),
if mod(p, 2)=0 then v: [1, 0, vv[j+1]] else v: [0, 1, 1],
for i from j thru 1 step(-1) do
(h: vv[i], u: [v[1]+h*v[3], v[3], 2*h*v[1]+v[2]+h^2*v[3]], v: u),
return(v)),
for n from 1 thru nmax do
(w: sqrt(n), m: floor(w),
if w=m then b: n else
(v: coeff(w, m), x: v[1], y: v[2], z: v[3], q: mod(z, 2),
if q=0 then (z: z/2, y: y/2) else x: 2*x,
fr: (x*m+y)/z, m: m+z, fr: fr+x, b: m^2+fr),
printf( fl, "~d, ", b)),
close(fl));
(Python)
from sympy import floor, S, sqrt
def coeff(w, m):
a, p, s, vv = m, 0, w, []
while a < 2*m:
p += 1
s = S.One/(s-floor(s))
a = floor(s)
if a < 2*m:
vv.append(a)
j = (p-1)//2
v = [0, 1, 1] if p % 2 else [1, 0, vv[j]]
for i in range(j-1, -1, -1):
h = vv[i]
v = [v[0]+h*v[2], v[2], 2*h*v[0]+v[1]+h**2*v[2]]
return v
def A334116(n):
w = sqrt(n)
m = floor(w)
if w == m:
return n
else:
x, y, z = coeff(w, m)
if z % 2:
x *= 2
else:
z //= 2
y //= 2
return (m+z)**2+x+(x*m+y)//z # Chai Wah Wu, Sep 30 2021, after Maxima code
CROSSREFS
KEYWORD
nonn
AUTHOR
Gerhard Kirchner, Apr 14 2020
STATUS
approved