%I #27 Sep 08 2022 08:45:38
%S 8,24,48,80,168,224,360,440,728,840,1088,1224,1368,1848,2208,2400,
%T 3024,3720,3968,4760,5040,5624,5928,7920,8648,10608,11448,13688,14160,
%U 14640,16128,17160,18224,19320,21024,24024,25920,28560,29928,31328,33488
%N Nonsquarefree numbers such that n-1 is prime and n+1 is square.
%C Also numbers n > 3 such that n-1 is prime and n+1 is square.
%C Sequence gives values x of fundamental solution (x,y) to Pellian x^2 - D*y^2 = 1, with D = n-1 = A049002, corresponding values y being sqrt(n+1) = A028870. (Substituting back into the Pellian we indeed have n^2 - (n-1)(n+1) = 1.) - _Lekraj Beedassy_, Feb 23 2019
%H Charles R Greathouse IV, <a href="/A146980/b146980.txt">Table of n, a(n) for n = 1..10000</a>
%t Select[Range[35000], !SquareFreeQ[#] && PrimeQ[#-1] && IntegerQ[Sqrt[#+1] ] &] (* _G. C. Greubel_, Feb 22 2019 *)
%t Mean/@SequencePosition[Table[Which[PrimeQ[n],1,IntegerQ[Sqrt[ n]],3,!SquareFreeQ[ n],2,True,0],{n,33500}],{1,2,3}] (* Requires Mathematica version 10 or later *) (* _Harvey P. Dale_, Sep 03 2020 *)
%o (Magma) [ n: n in [1..35000] | not IsSquarefree(n) and IsPrime(n-1) and IsSquare(n+1) ]; // _Klaus Brockhaus_, Nov 05 2008
%o (PARI) list(lim)=my(v=List()); forstep(k=3,sqrtint(lim\1+1),2, if(isprime(k^2-2), listput(v,k^2-1))); Vec(v) \\ _Charles R Greathouse IV_, Jun 13 2017
%o (Sage) [n for n in (1..35000) if not is_squarefree(n) and is_prime(n-1) and is_square(n+1)] # _G. C. Greubel_, Feb 22 2019
%Y Cf. A013929, A049002.
%K nonn
%O 1,1
%A _Giovanni Teofilatto_, Nov 04 2008
%E Extended beyond a(6) by _Klaus Brockhaus_, Nov 05 2008