[Axiom-math] continued fractions

[Martin Rubey]

Does anybody here know how the arithmetic for continued fractions in CONTFRAC
is supposed to work?

More precisely, it does work for CONTFRAC INT, but not for

CONTFRAC UnivariatePolynomial(x, Fraction Integer)

(in the following: CONTFRAC UP(x, FRAC INT)

but I'm not even sure, whether it *should* work.

Stephen: since you are mentioned as author, is there anything where I can read
about it? Maybe the following hint helps to refresh memory (after all, it's
from 1991):

the main routine is (I believe)

genFromSequence: Stream Q -> %

which takes a stream of approximants and returns (I believe: ) a reduced
continued fraction:

eucWhole(a: Q): R == numer a quo denom a

eucWhole0(a: Q): R ==
isOrdered =>
n := numer a
d := denom a
q := n quo d
r := n - q*d
if r < 0 then q := q - 1
q
eucWhole a

genFromSequence apps ==
lo := first apps; apps := rst apps
hi := first apps; apps := rst apps
while eucWhole0 lo ^= eucWhole0 hi repeat
lo := first apps; apps := rst apps
hi := first apps; apps := rst apps
wh := eucWhole0 lo
[[wh, genReducedForm(wh::Q, apps, moebius(1,0,0,1))], canReduce?]

genReducedForm(wh0, apps, mt) ==
lo: Q := first apps - wh0; apps := rst apps
hi: Q := first apps - wh0; apps := rst apps
lo = hi and zero? eval(mt, lo) => empty()
mt := recip mt
wlo := eucWhole eval(mt, lo)
whi := eucWhole eval(mt, hi)
while wlo ^= whi repeat
wlo := eucWhole eval(mt, first apps - wh0); apps := rst apps
whi := eucWhole eval(mt, first apps - wh0); apps := rst apps
concat([1,wlo], delay genReducedForm(wh0, apps, shift(mt, -wlo::Q)))

In the introduction it says

++ Description: \spadtype{ContinuedFraction} implements general
++ continued fractions. This version is not restricted to simple,
++ finite fractions and uses the \spadtype{Stream} as a
++ representation. The arithmetic functions assume that the
++ approximants alternate below/above the convergence point.
++ This is enforced by ensuring the partial numerators and partial
++ denominators are greater than 0 in the Euclidean domain view of \spad{R}
++ (i.e. \spad{sizeLess?(0, x)}).

but that would make any reduced continued fraction in UP(x, FRAC INT) illegal:
the Euclidean size of any rational number, in particular of 1, is zero in that
domain.



Help would be greatly appreciated!

Martin

[Martin Rubey]

The continued fraction domain is, however, described (but not in depth) in
the paper
Infinite Structures in Scratchpad II, W.H. Burge and S.M. Watt,
pp. 138-148, Proc. 1987 European Conference on Computer Algebra, (EUROCAL
87), June 2-5 1987, Leipzig, German Democratic Republic, LNCS 378, Springer
1989.
I've attached the preprint version from 1987.

Thank you for the paper. Unfortunately, it only (roughly) describes what
CONTFRAC can do, but not how it does it or where the algorithm would be
described.

I guess I'll have to do some reverse-engineering. If you have time though, I'd
appreciate some help.

Have a good stay,

Martin

[Martin Rubey]

Dear Martin,
The quick answer is that I don't remember off the top of my head what
restrictions exist on the domain parameter.
I'm on the road now (just 10 minutes ago checked into my hotel) and so
don't have access to a couple of things that would refresh my memory.

Stephen, do you happen to remember whether you have been aware of Gosper's
algorithms for continued fractions? I cannot be sure, but it seems that
CONTFRAC is not using them, is it?

Martin

[Stephen Watt]

The continued fraction domain is using Moebius transforms internally
(z -> (az+b)/(cz+d)) like Gosper. I wasn't specifically aware of his
work at the time, though.

The main issue is not in converting a particular value to a continued
fraction, which is easy, but rather to construct a c.f. out of a stream
of convergents, e.g. obtained as the result of arithmetic on two other c.f.s.

-- Stephen

Post by Martin Rubey

Dear Martin,
The quick answer is that I don't remember off the top of my head what
restrictions exist on the domain parameter.
I'm on the road now (just 10 minutes ago checked into my hotel) and so
don't have access to a couple of things that would refresh my memory.

Stephen, do you happen to remember whether you have been aware of Gosper's
algorithms for continued fractions? I cannot be sure, but it seems that
CONTFRAC is not using them, is it?
Martin

[Stephen Watt]

P.S. I agree it would be interesting to work through the details of using
Gosper's z(x,y) = (a x y + b x + c y + d) / (e x y + f x + g y + h)
compositions in a categorical framework.

-- Stephen

Post by Martin Rubey

Dear Martin,
The quick answer is that I don't remember off the top of my head what
restrictions exist on the domain parameter.
I'm on the road now (just 10 minutes ago checked into my hotel) and so
don't have access to a couple of things that would refresh my memory.

Stephen, do you happen to remember whether you have been aware of Gosper's
algorithms for continued fractions? I cannot be sure, but it seems that
CONTFRAC is not using them, is it?
Martin

[Stephen Watt]

Dear Martin,

The quick answer is that I don't remember off the top of my head what
restrictions exist on the domain parameter.

I'm on the road now (just 10 minutes ago checked into my hotel) and so
don't have access to a couple of things that would refresh my memory.

The continued fraction domain is, however, described (but not in depth) in
the paper

Infinite Structures in Scratchpad II, W.H. Burge and S.M. Watt, pp. 138-148,
Proc. 1987 European Conference on Computer Algebra, (EUROCAL 87),
June 2-5 1987, Leipzig, German Democratic Republic, LNCS 378, Springer 1989.

I've attached the preprint version from 1987.

In any case, the construction (if not this incarnation) should work for
continued fractions with Q[x] entries, as shown on page 8 of the preprint.

-- Stephen (in Hong Kong)

Post by Martin Rubey
Does anybody here know how the arithmetic for continued fractions in CONTFRAC
is supposed to work?
More precisely, it does work for CONTFRAC INT, but not for
CONTFRAC UnivariatePolynomial(x, Fraction Integer)
(in the following: CONTFRAC UP(x, FRAC INT)
but I'm not even sure, whether it *should* work.
Stephen: since you are mentioned as author, is there anything where I can read
about it? Maybe the following hint helps to refresh memory (after all, it's
the main routine is (I believe)
genFromSequence: Stream Q -> %
which takes a stream of approximants and returns (I believe: ) a reduced
eucWhole(a: Q): R == numer a quo denom a
eucWhole0(a: Q): R ==
isOrdered =>
n := numer a
d := denom a
q := n quo d
r := n - q*d
if r < 0 then q := q - 1
q
eucWhole a
genFromSequence apps ==
lo := first apps; apps := rst apps
hi := first apps; apps := rst apps
while eucWhole0 lo ^= eucWhole0 hi repeat
lo := first apps; apps := rst apps
hi := first apps; apps := rst apps
wh := eucWhole0 lo
[[wh, genReducedForm(wh::Q, apps, moebius(1,0,0,1))], canReduce?]
genReducedForm(wh0, apps, mt) ==
lo: Q := first apps - wh0; apps := rst apps
hi: Q := first apps - wh0; apps := rst apps
lo = hi and zero? eval(mt, lo) => empty()
mt := recip mt
wlo := eucWhole eval(mt, lo)
whi := eucWhole eval(mt, hi)
while wlo ^= whi repeat
wlo := eucWhole eval(mt, first apps - wh0); apps := rst apps
whi := eucWhole eval(mt, first apps - wh0); apps := rst apps
concat([1,wlo], delay genReducedForm(wh0, apps, shift(mt, -wlo::Q)))
In the introduction it says
++ Description: \spadtype{ContinuedFraction} implements general
++ continued fractions. This version is not restricted to simple,
++ finite fractions and uses the \spadtype{Stream} as a
++ representation. The arithmetic functions assume that the
++ approximants alternate below/above the convergence point.
++ This is enforced by ensuring the partial numerators and partial
++ denominators are greater than 0 in the Euclidean domain view of \spad{R}
++ (i.e. \spad{sizeLess?(0, x)}).
the Euclidean size of any rational number, in particular of 1, is zero in that
domain.
Help would be greatly appreciated!
Martin