cp's OEIS Frontend

This is a front-end for the Online Encyclopedia of Integer Sequences, made by Christian Perfect. The idea is to provide OEIS entries in non-ancient HTML, and then to think about how they're presented visually. The source code is on GitHub.

Showing 1-2 of 2 results.

A228929 Optimal ascending continued fraction expansion of Pi - 3.

Original entry on oeis.org

7, -113, 4739, -46804, 134370, -614063, 1669512, -15474114, -86232481, 1080357006, -8574121305, -24144614592, 133884333083, -2239330253016, -6347915250018, 14541933941298, -42301908155404, -298013673554972, 5177473084279656, -46709468571434452, 1201667304102142095, -68508286025632748778, 850084640720511629243, -2458418086834560217354
Offset: 1

Views

Author

Giovanni Artico, Sep 08 2013

Keywords

Comments

Definition of the expansion: for a positive real number x, there is always a unique sequence of signed integers with increasing absolute value |a(i)|>|a(i-1)| such that x =floor(x)+ 1/a(1) + 1/a(1)/a(2) + 1/a(1)/a(2)/a(3) + 1/a(1)/a(2)/a(3)/a(4) ... or equivalently x=floor(x)+1/a(1)*(1+1/a(2)*(1+1/a(3)*(1+1/a(4)*(1+...)))) giving the fastest converging series with this representation. This formula can be represented as a regular ascending continued fraction. The expansion is similar to Engel and Pierce expansions, but the sign of the terms is not predefined and determined by the algorithm for optimizing the convergence.
For x rational number the expansion has a finite number of terms, for x irrational an infinite number. Empirically the sequence doesn't show any evident regularity except in some interesting cases.

Examples

			Pi = 3 + 1/7*(1 - 1/113*(1 + 1/4739*(1 - 1/46804*(1 + 1/134370*(1 - 1/614063*(1 + 1/1669512*(1 + ...))))))).
		

Crossrefs

Programs

  • Derive
    ArticoExp(x, n) := VECTOR(ROUND(1, ABS(k))*SIGN(k), k, ITERATES(ROUND(1, ABS(u))*ABS(u) - 1, u, MOD(x), n))
    Precision:=Mixed
    PrecisionDigits:=10000
    ArticoExp(PI,20)
  • Maple
    # Slow procedure valid for every number
    ArticoExp := proc (n, q::posint)::list; local L, i, z; Digits := 50000; L := []; z := n-floor(n); for i to q+1 do if z = 0 then break end if; L := [op(L), round(1/abs(z))*sign(evalf(z))]; z := abs(z)*round(1/abs(z))-1 end do; return L end proc
    # Fast procedure, not suited for rational numbers
    ArticoExp := proc (n, q::posint)::list; local L, i, z; Digits := 50000; L := []; z := frac(evalf(n)); for i to q+1 do if z = 0 then break end if; L := [op(L), round(1/abs(z))*sign(z)]; z := abs(z)*round(1/abs(z))-1 end do; return L end proc
    # List the first 20 terms of the expansion of Pi
    ArticoExp(Pi,20)
  • Mathematica
    ArticoExp[x_, n_] :=  Round[1/#] & /@ NestList[Round[1/Abs[#]]*Abs[#] - 1 &, FractionalPart[x], n]; Block[{$MaxExtraPrecision = 50000}, ArticoExp[Pi, 20]]

Formula

Given a positive real number x, let z(0)=x-floor(x) and z(k+1)=abs(z(k))*round(1/abs(z(k)))-1 ; then a(n)=sign(z(n))*round(1/abs(z(n))) for n>0.

A319749 a(n) is the numerator of the Heron sequence with h(0)=3.

Original entry on oeis.org

3, 11, 119, 14159, 200477279, 40191139395243839, 1615327685887921300502934267457919, 2609283532796026943395592527806764363779539144932833602430435810559
Offset: 0

Views

Author

Paul Weisenhorn, Sep 27 2018

Keywords

Comments

The denominator of the Heron sequence is in A319750.
The following relationship holds between the numerator of the Heron sequence and the numerator of the continued fraction A041018(n)/A041019(n) convergent to sqrt(13).
n even: a(n)=A041018((5*2^n-5)/3).
n odd: a(n)=A041018((5*2^n-1)/3).
More generally, all numbers c(n)=A078370(n)=(2n+1)^2+4 have the same relationship between the numerator of the Heron sequence and the numerator of the continued fraction convergent to 2n+1.
sqrt(c(n)) has the continued fraction 2n+1; n,1,1,n,4n+2.
hn(n)^2-c(n)*hd(n)^2=4 for n>1.
From Peter Bala, Mar 29 2022: (Start)
Applying Heron's method (sometimes called the Babylonian method) to approximate the square root of the function x^2 + 4, starting with a guess equal to x, produces the sequence of rational functions [x, 2*T(1,(x^2+2)/2)/x, 2*T(2,(x^2+2)/2)/( 2*x*T(1,(x^2+2)/2) ), 2*T(4,(x^2+2)/2)/( 4*x*T(1,(x^2+2)/2)*T(2,(x^2+2)/2) ), 2*T(8,(x^2+2)/2)/( 8*x*T(1,(x^2+2)/2)*T(2,(x^2+2)/2)*T(4,(x^2+2)/2) ), ...], where T(n,x) denotes the n-th Chebyshev polynomial of the first kind. The present sequence is the case x = 3. Cf. A001566 and A058635 (case x = 1), A081459 and A081460 (essentially the case x = 4). (End)

Examples

			A078370(2)=29.
hn(0)=A041046(0)=5; hn(1)=A041046(3)=27; hn(2)=A041046(5)=727;
hn(3)=A041046(13)=528527.
		

Crossrefs

2*T(2^n,x/2) modulo differences of offset: A001566 (x = 3 and x = 7), A003010 (x = 4), A003487 (x = 5), A003423 (x = 6), A346625 (x = 8), A135927 (x = 10), A228933 (x = 18).

Programs

  • Maple
    hn[0]:=3:  hd[0]:=1:
    for n from 1 to 6 do
    hn[n]:=(hn[n-1]^2+13*hd[n-1]^2)/2:
    hd[n]:=hn[n-1]*hd[n-1]:
       printf("%5d%40d%40d\n", n, hn[n], hd[n]):
    end do:
    #alternative program
    a := n -> if n = 0 then 3 else simplify( 2*ChebyshevT(2^(n-1), 11/2) ) end if:
    seq(a(n), n = 0..7); # Peter Bala, Mar 16 2022
  • Python
    def aupton(nn):
        hn, hd, alst = 3, 1, [3]
        for n in range(nn):
            hn, hd = (hn**2 + 13*hd**2)//2, hn*hd
            alst.append(hn)
        return alst
    print(aupton(7)) # Michael S. Branicky, Mar 16 2022

Formula

h(n) = hn(n)/hd(n); hn(0)=3; hd(0)=1.
hn(n+1) = (hn(n)^2+13*hd(n)^2)/2.
hd(n+1) = hn(n)*hd(n).
A041018(n) = A010122(n)*A041018(n-1) + A041018(n-2).
A041019(n) = A010122(n)*A041019(n-1) + A041019(n-2).
From Peter Bala, Mar 16 2022: (Start)
a(n) = 2*T(2^(n-1),11/2) for n >= 1, where T(n,x) denotes the n-th Chebyshev polynomial of the first kind.
a(n) = 2*T(2^n, 3*sqrt(-1)/2) for n >= 2.
a(n) = ((11 + 3*sqrt(13))/2)^(2^(n-1)) + ((11 - 3*sqrt(13))/2)^(2^(n-1)) for n >= 1.
a(n+1) = a(n)^2 - 2 for n >= 1.
a(n) = A057076(2^(n-1)) for n >= 1.
Engel expansion of (1/6)*(13 - 3*sqrt(13)); that is, (1/6)*(13 - 3*sqrt(13)) = 1/3 + 1/(3*11) + 1/(3*11*119) + .... (Define L(n) = (1/2)*(n - sqrt(n^2 - 4)) for n >= 2 and show L(n) = 1/n + L(n^2-2)/n. Iterate this relation with n = 11. See also Liardet and Stambul, Section 4.)
sqrt(13) = 6*Product_{n >= 0} (1 - 1/a(n)).
sqrt(13) = (9/5)*Product_{n >= 0} (1 + 2/a(n)). See A001566. (End)

Extensions

a(6) and a(7) added by Peter Bala, Mar 16 2022
Showing 1-2 of 2 results.