%I #29 Feb 09 2025 03:40:25
%S 2,2,5,5,13,13,5,13,89,89,233,233,89,233,1597,1597,5,1597,1597,13,
%T 28657,28657,13,28657,28657,5,514229,514229,2,514229,514229,89,13,89,
%U 89,13,233,233,0,233,433494437,433494437,5,433494437,2971215073,2971215073,13,2971215073
%N Greatest prime Fibonacci divisor of F(n)^2 + 1 where F(n) is the n-th Fibonacci number, or 0 if no such prime factor exists.
%C a(n) = 0 for n = 39, 60, 69, 72, ... .
%C a(5385) has 1126 decimal digits. - _Chai Wah Wu_, Nov 19 2020
%H Chai Wah Wu, <a href="/A338762/b338762.txt">Table of n, a(n) for n = 1..5384</a> (n = 1..1000 from Alois P. Heinz)
%e a(6) = 13 because F(6)^2 + 1 = 8^2 + 1 = 65 = 5*13 and 13 is the greatest prime Fibonacci divisor.
%p a:= proc(n) local F, m, t; F, m, t:=
%p [1, 2], 0, (<<0|1>, <1|1>>^n)[2, 1]^2+1;
%p while F[2]<=t do if isprime(F[2]) and irem(t, F[2])=0
%p then m:=F[2] fi; F:= [F[2], F[1]+F[2]]
%p od; m
%p end:
%p seq(a(n), n=1..50); # _Alois P. Heinz_, Nov 07 2020
%t a[n_] := Module[{F, m, t}, F = {1, 2}; m = 0; t = MatrixPower[{{0, 1}, {1, 1}}, n][[2, 1]]^2 + 1; While[F[[2]] <= t, If[PrimeQ[F[[2]]] && Mod[t, F[[2]]] == 0, m = F[[2]]]; F = {F[[2]], F[[1]] + F[[2]]}]; m];
%t Table[a[n], {n, 1, 50}] (* _Jean-François Alcover_, Feb 09 2025, after _Alois P. Heinz_ *)
%o (PARI) isfib(n) = my(k=n^2); k+=(k+1)<<2; issquare(k) || (n>0 && issquare(k-8));
%o a(n) = my(f=factor(fibonacci(n)^2+1)[,1]~, v=select(x->isfib(x), f)); if (#v, vecmax(v), 0); \\ _Michel Marcus_, Nov 07 2020
%Y Cf. A000045, A005478, A245306.
%K nonn
%O 1,1
%A _Michel Lagneau_, Nov 07 2020