(* From Primoz Pirnat,Jun 15 2024 *) (* the program is written in Mathematica 11.3. *) Fn[n_] := Fibonacci[n]; A005653[num_] := Module[{n = 2 num, st, stm, cif, fl = 1, tip, odm1, odm2, d = {0, 0}, rd = 0, j, str = ""}, stm = Quotient[Floor[Log[GoldenRatio, n + 3/(10 n)] + 5.675], 6]*6; odm1 = {{-1, 1}, {0, 0}, {-1, 1}, {0, 0}, {-1, 1}}; odm2 = {{1, -1}, {-1, 1}, {1, -1}}; For[st = stm, st >= 6, st -= 6, cif = {Fn[st]/2, Fn[st], 3 Fn[st]/2, 2 Fn[st], 5 Fn[st]/2}; If[cif[[1]] > n, Goto[forward]]; tip = 5; While[tip > 0, If[n >= cif[[tip]], Break[]]; tip-- ]; d += odm1[[tip]]*fl; n -= cif[[tip]]; If[MemberQ[{1, 3, 5}, tip], fl *= -1]; Label[forward]; cif = {Fn[st - 4], Fn[st - 4] + Fn[st - 3]/2, Fn[st - 2]}; If[st == 6 && n == 0, d -= {1, -1}*fl]; If[cif[[1]] > n, Continue[]]; tip = 3; While[tip > 0, If[n >= cif[[tip]], Break[]]; tip-- ]; If[st == 6, odm2 = {{0, 0}, {0, 0}, {0, 0}}]; d += odm2[[tip]]*fl; n -= cif[[tip]]; If[tip == 2, fl *= -1]; ]; Block[{$MaxExtraPrecision = Max[50, Round[N[Log[10, 2 num]]]]}, For[j = 0, j <= Abs[d[[1]]], j++, If[FractionalPart[ N[(2 num + j*Sign[d[[1]]]) GoldenRatio, $MaxExtraPrecision + 10]] > 1/2, str = str <> "\[FilledCircle]", str = str <> "\[EmptyCircle]"] ] ]; If[d[[1]] != 0, If[d[[1]] > 0, rd = StringCount[str, "\[FilledCircle]\[FilledCircle]"] - StringCount[str, "\[EmptyCircle]\[EmptyCircle]"], rd = StringCount[str, "\[EmptyCircle]\[EmptyCircle]"] - StringCount[str, "\[FilledCircle]\[FilledCircle]"] ] ]; Return[2 num + d[[1]] + rd] ] SetAttributes[A005653, Listable];