$SPAD/src/input contfrac.input

Timothy Daly December 3, 2016

Abstract

Contents

2

1 License

— license —

--Copyright The Numerical Algorithms Group Limited 1991.

———-

— * —

)set break resume
)spool contfrac.output
)set message test on
)set message auto off
)clear all
--S 1 of 40
r1 := 3/4
--R
--R
--R        3
--R   (1)  -
--R        4
--R                                                      Type: Fraction(Integer)
--E 1
--S 2 of 40
r2 := 314159/100000
--R
--R
--R        314159
--R   (2)  ------
--R        100000
--R                                                      Type: Fraction(Integer)
--E 2
--S 3 of 40
c1 := r1 :: ContinuedFraction Integer
--R
--R
--R          1 |     1 |
--R   (3)  +---+ + +---+
--R        | 1     | 3
--R                                             Type: ContinuedFraction(Integer)
--E 3
--S 4 of 40
c2 := r2 :: ContinuedFraction Integer
--R
--R

--R              1 |     1  |     1 |     1  |     1 |     1 |     1 |
--R   (4)  3 + +---+ + +----+ + +---+ + +----+ + +---+ + +---+ + +---+
--R            | 7     | 15     | 1     | 25     | 1     | 7     | 4
--R                                             Type: ContinuedFraction(Integer)
--E 4
-- We can view these in the list notation
--S 5 of 40
partialQuotients c1
--R
--R
--R   (5)  [0,1,3]
--R                                                        Type: Stream(Integer)
--E 5
--S 6 of 40
partialQuotients c2
--R
--R
--R   (6)  [3,7,15,1,25,1,7,4]
--R                                                        Type: Stream(Integer)
--E 6
-- These are algebraic objects, so we can manipulate them accordingly
--S 7 of 40
c1 + c2
--R
--R
--R   (7)
--R         1 |     1 |     1 |     1 |     1 |     1 |     1 |     1  |     1 |
--R   3 + +---+ + +---+ + +---+ + +---+ + +---+ + +---+ + +---+ + +----+ + +---+
--R       | 1     | 8     | 4     | 2     | 5     | 1     | 2     | 32     | 2
--R                                             Type: ContinuedFraction(Integer)
--E 7
--S 8 of 40
c1 * c2
--R
--R
--R   (8)
--R           1 |     1 |     1 |     1 |     1 |     1 |     1  |     1 |     1 |
--R     2 + +---+ + +---+ + +---+ + +---+ + +---+ + +---+ + +----+ + +---+ + +---+
--R         | 2     | 1     | 4     | 5     | 6     | 2     | 13     | 1     | 1
--R   +
--R       1 |
--R     +---+ + ...
--R     | 1
--R                                             Type: ContinuedFraction(Integer)
--E 8

--S 9 of 40
1 / c2
--R
--R
--R          1 |     1 |     1  |     1 |     1  |     1 |     1 |     1 |
--R   (9)  +---+ + +---+ + +----+ + +---+ + +----+ + +---+ + +---+ + +---+
--R        | 3     | 7     | 15     | 1     | 25     | 1     | 7     | 4
--R                                             Type: ContinuedFraction(Integer)
--E 9
--S 10 of 40
c1 - c2
--R
--R
--R   (10)
--R           1 |     1 |     1 |     1 |     1 |     1 |     1 |     1 |      1  |
--R   - 3 + +---+ + +---+ + +---+ + +---+ + +---+ + +---+ + +---+ + +---+ + +-----+
--R         | 1     | 1     | 1     | 1     | 4     | 6     | 2     | 2     | 131
--R                                             Type: ContinuedFraction(Integer)
--E 10
--S 11 of 40
c2 - c1
--R
--R
--R               1 |     1 |     1 |     1 |     1 |     1 |     1 |      1  |
--R   (11)  2 + +---+ + +---+ + +---+ + +---+ + +---+ + +---+ + +---+ + +-----+
--R             | 2     | 1     | 1     | 4     | 6     | 2     | 2     | 131
--R                                             Type: ContinuedFraction(Integer)
--E 11
-- and can convert them back to rational numbers.
--S 12 of 40
convergents %
--R
--R
--R            5 7 12 55 342 739 1820 239159
--R   (12)  [2,-,-,--,--,---,---,----,------]
--R            2 3  5 23 143 309  761 100000
--R                                              Type: Stream(Fraction(Integer))
--E 12
)clear all
-- Continued fractions over other Euclidean domains
--S 13 of 40
a0 := ((-122 + 597* %i)/(4 - 4*%i))
--R
--R

--R          719   475
--R   (1)  - --- + --- %i
--R           8     8
--R                                             Type: Complex(Fraction(Integer))
--E 13
--S 14 of 40
b0 := ((-595 - %i)/(3 - 4*%i))
--R
--R
--R          1781   2383
--R   (2)  - ---- - ---- %i
--R           25     25
--R                                             Type: Complex(Fraction(Integer))
--E 14
--S 15 of 40
a  := continuedFraction(a0)
--R
--R
--R                           1    |         1     |
--R   (3)  - 90 + 59%i + +---------+ + +-----------+
--R                      | 1 - 2%i     | - 1 + 2%i
--R                                    Type: ContinuedFraction(Complex(Integer))
--E 15
--S 16 of 40
b  := continuedFraction(b0)
--R
--R
--R                            1     |      1  |
--R   (4)  - 71 - 95%i + +-----------+ + +-----+
--R                      | - 1 + 2%i     | - 2
--R                                    Type: ContinuedFraction(Complex(Integer))
--E 16
--S 17 of 40
a + b
--R
--R
--R                             1     |      1  |        1     |
--R   (5)  - 161 - 36%i + +-----------+ + +-----+ + +----------+
--R                       | - 7 - 3%i     | 2%i     | - 1 - %i
--R                                    Type: ContinuedFraction(Complex(Integer))
--E 17
--S 18 of 40
convergents %
--R
--R

--R                      - 1020 - 735%i - 2004 - 1631%i 362 - 4655%i
--R   (6)  [- 161 - 36%i,--------------,---------------,------------]
--R                          7 + 3%i        14 + 7%i      4 + 28%i
--R                                     Type: Stream(Fraction(Complex(Integer)))
--E 18
--S 19 of 40
last % - (a0 + b0)
--R
--R
--R   (7)  0
--R                                             Type: Complex(Fraction(Integer))
--E 19
--S 20 of 40
a / b
--R
--R
--R   (8)
--R                 1    |     1 |        1     |        1     |      1  |
--R     - %i + +---------+ + +---+ + +----------+ + +----------+ + +-----+
--R            | 4 - 8%i     | 3     | - 4 - %i     | - 2 - %i     | 3%i
--R   +
--R          1    |
--R     +---------+
--R     | 2 + 4%i
--R                                    Type: ContinuedFraction(Complex(Integer))
--E 20
--S 21 of 40
convergents %
--R
--R
--R   (9)
--R         4 - 7%i 13 - 21%i 69 - 64%i  215 - 80%i 709 - 171%i 2279 - 2022%i
--R   [- %i,-------,---------,---------,-----------,-----------,-------------]
--R         8 + 4%i 24 + 13%i 75 + 72%i 102 + 232%i 234 + 771%i 2376 + 2384%i
--R                                     Type: Stream(Fraction(Complex(Integer)))
--E 21
--S 22 of 40
last % - (a0/b0)
--R
--R
--R   (10)  0
--R                                             Type: Complex(Fraction(Integer))
--E 22
--S 23 of 40
(a = b)::Boolean

--R
--R
--R   (11)  false
--R                                                                Type: Boolean
--E 32
--S 24 of 40
c := continuedFraction(3 + 4*%i, repeating [1 + %i], repeating [5 - %i])
--R
--R
--R   (12)
--R                 1 + %i |     1 + %i |     1 + %i |     1 + %i |     1 + %i |
--R     3 + 4%i + +--------+ + +--------+ + +--------+ + +--------+ + +--------+
--R               | 5 - %i     | 5 - %i     | 5 - %i     | 5 - %i     | 5 - %i
--R   +
--R       1 + %i |     1 + %i |     1 + %i |     1 + %i |     1 + %i |
--R     +--------+ + +--------+ + +--------+ + +--------+ + +--------+ + ...
--R     | 5 - %i     | 5 - %i     | 5 - %i     | 5 - %i     | 5 - %i
--R                                    Type: ContinuedFraction(Complex(Integer))
--E 24
--S 25 of 40
a/c
--R
--R
--R   (13)
--R                        1     |          1      |        1     |         1     |
--R     - 1 + 20%i + +-----------+ + +-------------+ + +----------+ + +-----------+
--R                  | - 1 - 2%i     | - 11 - 16%i     | - 1 + %i     | - 9 - 2%i
--R   +
--R          1    |      1  |        1    |      1  |         1     |        1    |
--R     +---------+ + +-----+ + +---------+ + +-----+ + +-----------+ + +---------+
--R     | 1 + 2%i     | 2%i     | 8 - 2%i     | - 2     | - 1 + 8%i     | 3 - 3%i
--R   +
--R     ...
--R                                    Type: ContinuedFraction(Complex(Integer))
--E 25
-- (a = c)::Boolean -- should give error
--S 26 of 40
d := complete continuedFraction(3+4*%i, repeating [1+%i],[i-%i for i in 1..5])
--R
--R
--R   (14)
--R               1 + %i |     1 + %i |     1 + %i |     1 + %i |     1 + %i |
--R   3 + 4%i + +--------+ + +--------+ + +--------+ + +--------+ + +--------+
--R             | 1 - %i     | 2 - %i     | 3 - %i     | 4 - %i     | 5 - %i
--R                                    Type: ContinuedFraction(Complex(Integer))
--E 26

--S 27 of 40
(a = d)::Boolean
--R
--R
--R   (15)  false
--R                                                                Type: Boolean
--E 27
--S 28 of 40
q : Fraction UnivariatePolynomial('x, Fraction Integer)
--R
--R                                                                   Type: Void
--E 28
--S 29 of 40
q := (2*x**2 - x + 1) / (3*x**3 - x + 8)
--R
--R
--R         2  2   1     1
--R         - x  - - x + -
--R         3      3     3
--R   (17)  --------------
--R           3   1     8
--R          x  - - x + -
--R               3     3
--R                    Type: Fraction(UnivariatePolynomial(x,Fraction(Integer)))
--E 29
--S 30  of 40
c := continuedFraction q
--R
--R
--R              1    |          1      |            1        |
--R   (18)  +---------+ + +-------------+ + +-----------------+
--R         | 3     3     |   8     204     |    343     1421
--R         | - x + -     | - - x - ---     | - ---- x + ----
--R         | 2     4     |   7      49     |   6112     6112
--R           Type: ContinuedFraction(UnivariatePolynomial(x,Fraction(Integer)))
--E 30
--S 31 of 40
d := continuedFraction differentiate q
--R
--R
--R   (19)
--R              1         |          1       |               1            |
--R     +------------------+ + +--------------+ + +------------------------+
--R     |   3  2   3     9     |    1      47     |   69696     6381055032
--R     | - - x  - - x + -     | - -- x + ---     | - ----- x - ----------

--R     |   2      2     4     |   11     264     |   32963     1086559369
--R   +
--R                            1                      |
--R     +---------------------------------------------+
--R     |    35816256480347      79911907817759707445
--R     | - --------------- x + ---------------------
--R     |   218199728406528     307361810626373910528
--R   +
--R                                   1                              |
--R     +------------------------------------------------------------+
--R     |   39359803441398779644674048     9259889268740766802477056
--R     | - -------------------------- x + -------------------------
--R     |    1979914602262093317951025     1979914602262093317951025
--R           Type: ContinuedFraction(UnivariatePolynomial(x,Fraction(Integer)))
--E 31
--S 32 of 40
c/d
--R
--R
--R   (20)
--R           1         1     |              1          |
--R     - x - - + +-----------+ + +---------------------+
--R           2   | 6     219     |  686       83926465
--R               | - x + ---     | ----- x - ---------
--R               | 7      49     | 20131     405257161
--R   +
--R                        1                  |
--R     +-------------------------------------+
--R     |   8158231908091     299222030081511
--R     | - ------------- x - ---------------
--R     |    268162004432     350018456284868
--R   +
--R                              1                        |
--R     +-------------------------------------------------+
--R     |    7309785441053183312      1206863637270224864
--R     | - -------------------- x + --------------------
--R     |   30383172810229285385     30383172810229285385
--R           Type: ContinuedFraction(UnivariatePolynomial(x,Fraction(Integer)))
--E 32
--S 33 of 40
convergents %
--R
--R
--R   (21)
--R                2   40     121     3   14615  2   114691     168376
--R             - x  - -- x - ---  - x  + ----- x  - ------ x - ------
--R          1          7      84         40262      120786      20131
--R   [- x - -, -----------------, -----------------------------------,

--R          2            73               2   17373     307727
--R                   x + --              x  - ----- x + ------
--R                       14                   20131     120786
--R       4    3497  3    791   2   15030     18082
--R    - x  + ----- x  + ----- x  - ----- x + -----
--R           10442      31326       5221     15663
--R    --------------------------------------------,
--R            3   4359  2   48833     25886
--R           x  - ---- x  + ----- x - -----
--R                5221      31326      5221
--R       5   1  4   1  3   17  2   3     4
--R    - x  + - x  - - x  - -- x  + - x - -
--R           2      6       6      2     3
--R    ------------------------------------]
--R          4    3   11  2   16     7
--R         x  - x  + -- x  - -- x + -
--R                    6       3     6
--R            Type: Stream(Fraction(UnivariatePolynomial(x,Fraction(Integer))))
--E 33
--S 34 of 40
q/differentiate q
--R
--R
--R            5   1  4   1  3   17  2   3     4
--R         - x  + - x  - - x  - -- x  + - x - -
--R                2      6       6      2     3
--R   (22)  ------------------------------------
--R               4    3   11  2   16     7
--R              x  - x  + -- x  - -- x + -
--R                         6       3     6
--R                    Type: Fraction(UnivariatePolynomial(x,Fraction(Integer)))
--E 34
)clear all
)set streams calculate 7

———-

This file illustrates continued fractions.

Use the notation

ϕ(ai/bi, i = 1..n)

for continued fractions

a1/(b1 + (a2/b2 + ...(an/bn)...))

1/(e − 1) may be written ϕ(i/i, i = 1..) — * —

--S 35 of 40

s := continuedFraction(0, expand [1..], expand [1..])
--R
--R
--R          1 |     2 |     3 |     4 |     5 |     6 |     7 |
--R   (1)  +---+ + +---+ + +---+ + +---+ + +---+ + +---+ + +---+ + ...
--R        | 1     | 2     | 3     | 4     | 5     | 6     | 7
--R                                             Type: ContinuedFraction(Integer)
--E 35

———-

Euler discovered the relation

(e − 1)/(e + 1) = ϕ(1/(4i − 2), i = 1..)

— * —

--S 36 of 40
t := reducedContinuedFraction(0, [4*i-2 for i in 1..])
--R
--R
--R          1 |     1 |     1  |     1  |     1  |     1  |     1  |
--R   (2)  +---+ + +---+ + +----+ + +----+ + +----+ + +----+ + +----+ + ...
--R        | 2     | 6     | 10     | 14     | 18     | 22     | 26
--R                                             Type: ContinuedFraction(Integer)
--E 36

———-

Arithmetic on infinite continued fractions is supported. The results are given in reduced form. We illustrate by using the values s = 1/(e − 1) and t = (e − 1)/(e + 1) to recover the expansion for e.

— * —

--S 37 of 40
e := 1/(s*t) - 1
--R
--R
--R              1 |     1 |     1 |     1 |     1 |     1 |     1 |
--R   (3)  2 + +---+ + +---+ + +---+ + +---+ + +---+ + +---+ + +---+ + ...
--R            | 1     | 2     | 1     | 1     | 4     | 1     | 1
--R                                             Type: ContinuedFraction(Integer)
--E 37
--S 38 of 40
c := convergents e
--R
--R
--R             8 11 19 87 106
--R   (4)  [2,3,-,--,--,--,---,...]

--R             3  4  7 32  39
--R                                              Type: Stream(Fraction(Integer))
--E 38
--S 39 of 40
for i in 1..15 repeat
  output numeric c.i
--R
--R   2.0
--R   3.0
--R   2.6666666666 666666667
--R   2.75
--R   2.7142857142 857142857
--R   2.71875
--R   2.7179487179 487179487
--R   2.7183098591 549295775
--R   2.7182795698 924731183
--R   2.7182835820 895522388
--R   2.7182817182 817182817
--R   2.7182818352 059925094
--R   2.7182818229 439497119
--R   2.7182818287 356957267
--R   2.7182818284 45401318
--R                                                                   Type: Void
--E 39
--S 40 of 40
(s = t)::Boolean
--R
--R
--R   (6)  false
--R                                                                Type: Boolean
--E 40
)spool
)lisp (bye)

———-

References^{[1] nothing}

14