%I #25 Aug 05 2023 04:12:06
%S 12,20,40,70,88,104,464,650,1504,1888,1952,4030,5830,7192,7912,8925,
%T 9555,10792,13736,17272,30555,30592,32128,32445,78975,130304,442365,
%U 521728,522752,1713592,1848964,4526272,8353792,8378368,8382464,9928792,11547352,17999992
%N Admirable numbers such that the subtracted divisor is a Fibonacci number.
%C Subsequence of A111592.
%C The corresponding Fibonacci numbers are given by the sequence {b(n)} = 2, 1, 5, 2, 2, 1, 1, 1, 8, 2, 1, 2, 2, 8, 8, 3, 21, 8, 34, 8, 21, 8, 2, 3, 13, 1, 3, 2, 1, ....
%H Amiram Eldar, <a href="/A282754/b282754.txt">Table of n, a(n) for n = 1..46</a>
%H Terry Trotter, <a href="https://web.archive.org/web/20180129200429/http://www.trottermath.net/numthry/admirabl.html">Admirable Numbers</a>, 2009. [Wayback Machine copy] [Warning: As of March 2018 this site appears to have been hacked. Proceed with great caution. The original content should be retrieved from the Wayback machine and added here. - _N. J. A. Sloane_, Mar 29 2018]
%e 40 is in the sequence because sigma(40) - 2*5 = 90 - 10 = 80 = 2*40, where 5 is a Fibonacci number, or 1 + 2 + 4 + 8 + 10 + 20 - 5 = 40 where the subtracted divisor is 5.
%p with(numtheory):
%p for n from 1 to 20000 do:
%p x:=divisors(n):n0:=nops(x):
%p for i from 1 to n0 do:
%p u:=sqrt(5*x[i]^2-4):v:=sqrt(5*x[i]^2+4):
%p if (floor(u)=u or floor(v)=v) and sigma(n)-2*x[i]=2*n
%p then
%p printf(`%d %d \n`,n, x[i]):
%p else
%p fi:
%p od:
%p od:
%t With[{nn = 10^6}, Function[s, Flatten@ Position[#, 1] &@ Table[Total@ Boole@ Map[MemberQ[s, #] &, Select[Most@ Divisors@ n, Function[d, DivisorSigma[1, n] - 2 d == 2 n]]], {n, nn}]]@ Fibonacci@ Range[2 + Floor@ Log[GoldenRatio, nn]]] (* _Michael De Vlieger_, Feb 24 2017 *) (* or *)
%t fibQ[n_] := IntegerQ@ Sqrt[5 n^2 + 4] || IntegerQ@ Sqrt[5 n^2 - 4]; ok[n_] := Block[{d = DivisorSigma[1, n] - 2 n}, d>0 && EvenQ@d && Mod[n, d/2] == 0 && fibQ[d/2]]; Select[Range[10^6], ok] (* faster, _Giovanni Resta_, Mar 10 2017 *)
%o (PARI) isadmirable(n)=if(issquare(n)||issquare(n/2), 0, my(d=sigma(n)/2-n); (d>0 && d!=n && n%d==0)*d);
%o isfib(n) = my(k=n^2); k+=(k+1)<<2; issquare(k) || (n>0 && issquare(k-8))
%o isok(n) = (d=isadmirable(n)) && isfib(d); \\ _Michel Marcus_, Mar 10 2017
%Y Cf. A000045, A010056, A111592.
%K nonn
%O 1,1
%A _Michel Lagneau_, Feb 21 2017
%E More terms from _Michel Marcus_, Mar 10 2017