%I #37 Mar 07 2020 11:42:55
%S 1,1,11,173,22931,1319183,233526463,2673857519,39959591850371,
%T 8797116290975003,4872532317019728133,1657631603843299234219,
%U 247098748783812523360613,77729277912104164732573547,1503342018433974345747514544039
%N A sequence related to the period T of a simple gravity pendulum for arbitrary amplitudes.
%C For small angles the period T of a simple gravity pendulum obeys Christiaan Huygens’s law, i.e. T = 2*Pi*sqrt(L/g) with L the length of the pendulum and g the acceleration due to gravity. For arbitrary amplitudes the period T is given below, see Wikipedia. The Taylor series expansion of T as a function of the angular displacement phi leads for the numerators of the even powers of phi to the sequence given above and for the denominators to A223068.
%D C. D. Andriesse and Sally Miedema, Huygens: The Man Behind the Principle, Ch. 8, 2005.
%H Wikipedia, <a href="http://en.wikipedia.org/wiki/Pendulum">Pendulum</a> and <a href="http://en.wikipedia.org/wiki/Pendulum_(mathematics)">Pendulum mathematics</a>.
%F T = 2*Pi*sqrt(L/g)*(2/Pi)*K(sin(phi/2)) with K(k) the complete elliptic integral of the first kind.
%F T = 2*Pi*sqrt(L/g)/M(1,cos(phi/2)) where M(x,y) = (Pi/4)*((x+y)/(K((x-y)/(x+y)) is the arithmetic-geometric mean of x and y. - _Johannes W. Meijer_, Dec 28 2016
%F Let S = Sum_{n>=0} (-1)^n*euler(2*n)*x^n/(2*n) then a(n) = numerator(1/(2*n)! * [x^n] exp(S)). - _Peter Luschny_, Jan 05 2017
%e T = 2*Pi*sqrt(L/g) * (1 + (1/16)*phi^2 + (11/3072)*phi^4 + (173/737280)*phi^6 + … ).
%p nmax:=14: f := series((2/Pi)*EllipticK(sin(phi/2)), phi, 2*nmax+1): for n from 0 to nmax do a(n):= numer(coeff(f, phi, 2*n)) od: seq(a(n), n=0..nmax); # End first program.
%p nmax:=14: f := series(1/((Pi/4)*(1+cos(phi/2))/EllipticK((1-cos(phi/2))/(1+cos(phi/2)))), phi, 2*nmax+1): for n from 0 to nmax do a(n):= numer(coeff(f, phi, 2*n)) od: seq(a(n), n=0..nmax); # End second program. - _Johannes W. Meijer_, Dec 28 2016
%t s = Series[EllipticK[Sin[t/2]^2 ], {t, 0, 60}]; CoefficientList[s/Pi, t^2] // Numerator (* _Jean-François Alcover_, Oct 07 2014 *)
%o (Sage)
%o def A223067_list(prec):
%o P.<x> = PowerSeriesRing(QQ, default_prec=2*prec)
%o g = lambda x: exp(sum((-1)^k*euler_number(2*k)*x^k/(2*k) for k in (1..prec+1)))
%o R = P(g(x)).coefficients()
%o return [numerator(R[n]/factorial(2*n)) for n in (0..prec)]
%o print(A223067_list(14)) # _Peter Luschny_, Jan 05 2017
%Y Cf. A223068 (denominators), A019692 (2*Pi).
%Y Cf. A280442, A280443
%K nonn,easy,frac
%O 0,3
%A _Johannes W. Meijer_, Mar 14 2013