|
%I #27 Sep 08 2022 08:45:31
%S 1,1926,850137,2564464982,1132138928657,3415133918621062,
%T 1507685261236261801,4547981651299964079126,2007805569980855008712097,
%U 6056610836775865229750164742,2673822786819976661810784866297,8065673443881586606920210924732502
%N a(n) is the n-th pentagonal number which is the sum of two consecutive pentagonal numbers.
%C We solve the equation P(p) = P(r) + P(r+1) with unknowns p and r, equivalent to (6*p-1)^2 = 2*(6*r+2)^2 + 17. The Diophantine equation X^2 = 2*Y^2 + 17 whose solutions give p and r are obtained by (x(n), y(n)) such that:
%C x(1)=5, x(2)=215, x(3)=4517, x(4)=248087 and the same recurrence relation on the odd and even indices x(n+2) = 1154*x(n+1) - x(n);
%C y(1)=2, y(2)=152, y(3)=3194, y(4)=175424 and the same recurrence relation on the odd and even indices y(n+2) = 1154*y(n+1) - y(n).
%C The solutions (p,r) are given by the (u(n),v(n)) such that
%C u(1)=1, u(2)=36, u(3)=753, u(4)=41348 and the same recurrence relation on the odd and even indices u(n+2) = 1154*u(n+1) - u(n) - 192 or u(n+1) = 577*u(n) - 96 + 68*sqrt(72*u(n)^2 - 24*u(n) - 32);
%C v(1)=0, v(2)=25, v(3)=532, v(4)=29237 and the same recurrence relation on the odd and even indices v(n+2) = 1154*v(n+1) - v(n) + 384 or v(n+1) = 577*v(n) + 192 + 68*sqrt(72*u(n)^2 + 48*u(n) + 15).
%H Colin Barker, <a href="/A133301/b133301.txt">Table of n, a(n) for n = 1..300</a>
%H <a href="/index/Rec#order_05">Index entries for linear recurrences with constant coefficients</a>, signature (1,1331714,-1331714,-1,1).
%F For odd and even indices respectively: a(n+2) = 1331714*a(n+1) - a(n) - 416160; on the odd and the even indices respectively we also have a(n+1) = 665857*a(n) - 208080 + 19618*sqrt(1152*a(n)^2 - 720*a(n) - 32).
%F The g.f., h(z), is such that h(z) = a(1)*z + a(2)*z^2 + ... and is given by h(z) = z*(1 + 1925*z - 483503*z^2 + 65395*z^3 + 22*z^4)/((1-z)*(1 - 1331714*z^2 + z^4)).
%e With P(m) = m*(3*m-1)/2,
%e a(1)=1 because a(1) = P(1) = P(0) + P(1);
%e a(2)=1926 because P(36) = 1926 = P(25) + P(26) = 925 + 1001;
%e a(3)=850137 because P(753) = 850137 = P(532) + P(533) = 424270 + 425867 ...
%p a:=proc(m) if type (sqrt(18*m^2-6*m-8)/6-1/3), integer=true then m*(3*m-1)/2 else fi end : seq(a(m),m=1..100000)od; # _Emeric Deutsch_
%t # (3 # - 1)/2 &@ Select[Range[10^6], IntegerQ[Sqrt[18 #^2 - 6 # - 8]/6 - 1/3] &] (* or *)
%t Rest@ CoefficientList[Series[x*(1+1925*x-483503*x^2+65395*x^3+22*x^4)/( (1-x)*(1 -1154*x + x^2)*(1 +1154*x + x^2)), {x, 0, 12}], x] (* _Michael De Vlieger_, Jul 14 2016 *)
%o (PARI) Vec(x*(1+1925*x-483503*x^2+65395*x^3+22*x^4)/((1-x)*(1 -1154*x + x^2)*(1 +1154*x + x^2)) + O(x^20)) \\ _Colin Barker_, Oct 20 2014
%o (Magma) R<x>:=PowerSeriesRing(Integers(), 20); Coefficients(R!( x*(1+1925*x-483503*x^2+65395*x^3+22*x^4)/((1-x)*(1 -1154*x + x^2)*(1 +1154*x + x^2)) )); // _G. C. Greubel_, Mar 16 2019
%o (Sage) a=(x*(1+1925*x-483503*x^2+65395*x^3+22*x^4)/((1-x)*(1 -1154*x + x^2)*(1 +1154*x + x^2))).series(x, 20).coefficients(x, sparse=False); a[1:] # _G. C. Greubel_, Mar 16 2019
%Y Cf. A046173, A144796 (x(n)).
%K nonn,easy
%O 1,2
%A _Richard Choulet_, Dec 20 2007
%E Fixed typo in g.f. in formula, and more terms from _Colin Barker_, Oct 20 2014
| 