%I #18 May 13 2022 05:49:51
%S 1,1,1,1,5,21,109,2621,46409,1290665,143900249,10827502009,
%T 1650671059149,826986635160509,267834927512726725,
%U 226491424023176449909,497347467521206399078801,706331193237960728379843409
%N a(n) = (4*a(n-1)*a(n-3) + a(n-2)^2)/a(n-4), with a(1)=...=a(4)=1.
%C This is a (4,1) generalized Somos-4 sequence. - _Michael Somos_, May 12 2022
%H G. C. Greubel, <a href="/A121881/b121881.txt">Table of n, a(n) for n = 1..96</a>
%H A. N. W. Hone, <a href="http://arXiv.org/abs/0807.2538">Algebraic curves, integer sequences and a discrete Painlevé transcendent</a>, arXiv:0807.2538 [nlin.SI], 2008; Proceedings of SIDE 6, Helsinki, Finland, 2004.
%F a(n) = a(5-n) = (-a(n-1)*a(n-4) + 26*a(n-2)*a(n-3))/a(n-5) for all n in Z. - _Michael Somos_, May 12 2022
%p a:= proc(n) option remember;
%p if n < 5 then 1
%p else (4*a(n-1)*a(n-3) + a(n-2)^2)/a(n-4)
%p fi;
%p end proc;
%p seq(a(n), n = 1..30); # _G. C. Greubel_, Oct 08 2019
%t a[n_]:= a[n]= If[n<5, 1, (4*a[n-1]a[n-3] + a[n-2]^2)/a[n-4]]; Table[a[n], {n, 30}]
%o (PARI) my(m=30, v=concat([1,1,1,1], vector(m-4))); for(n=5, m, v[n] = (4*v[n-1]*v[n-3] + v[n-2]^2)/v[n-4]); v \\ _G. C. Greubel_, Oct 08 2019
%o (Magma) [n lt 5 select 1 else (4*Self(n-1)*Self(n-3) + Self(n-2)^2 )/Self(n-4): n in [1..30]]; // _G. C. Greubel_, Oct 08 2019
%o (Sage)
%o @CachedFunction
%o def a(n):
%o if (n<5): return 1
%o else: return (4*a(n-1)*a(n-3) + a(n-2)^2)/a(n-4)
%o [a(n) for n in (1..30)] # _G. C. Greubel_, Oct 08 2019
%o (GAP) a:=[1,1,1,1];; for n in [5..30] do a[n]:=(4*a[n-1]*a[n-3] + a[n-2]^2)/a[n-4]; od; a; # _G. C. Greubel_, Oct 08 2019
%Y Cf. Somos-4 sequence A006720, A014125, A095708.
%Y Cf. A121883.
%K nonn
%O 1,5
%A _Roger L. Bagula_, Sep 09 2006
%E Edited by _N. J. A. Sloane_, Sep 15 2006