Integers and rationals

In an earlier post, I asked readers to implement the bijection rational->integer and integer->rational. John Cowan suggested the Calkin-Wilf tree as a starting point. The Calkin-Wilf tree is a rooted binary tree where the nodes (or vertices, if you like) are labeled with positive rational numbers. It is infinite and complete: every node has two children. The Calkin-Wilf tree is constructed so that every rational number is assigned a unique node. Every positive rational number appears once and exactly once in the tree. The path from the root node to any selected rational is unique and can be encoded (in binary) as an integer.

1 ]=> (rational->integer 355/113)

;Value: 67107847

1 ]=> (integer->rational 67107847)

;Value: 355/113

1 ]=> (cwt/value *the-calkin-wilf-tree*)

;Value: 1

1 ]=> (cwt/value (cwt/left *the-calkin-wilf-tree*))

;Value: 1/2

1 ]=> (cwt/value (cwt/right *the-calkin-wilf-tree*))

;Value: 2

1 ]=> (cwt/value (cwt/left (cwt/right (cwt/left (cwt/left *the-calkin-wilf-tree*)))))

;Value: 4/7

Ho hum. We've all seen this sort of thing before.

Here's the unusual part:

1 ]=> cwt/left

;Value 1236: #[linear-fractional-transform 1236 x/(x + 1)]

1 ]=> cwt/right

;Value 1237: #[linear-fractional-transform 1237 (x + 1)]

So I can write

1 ]=> (cwt/value ((compose cwt/left cwt/right cwt/left cwt/left) *the-calkin-wilf-tree*))

;Value: 4/7

1 ]=> (lft/compose cwt/left cwt/right cwt/left cwt/left)

;Value 1260: #[linear-fractional-transform 1260 (3x + 1)/(5x + 2)]

(See "Playing with Linear Fractional Transforms")

---

The dyadic fractions are those rational numbers whose denominator is a power of 2. Numbers like 1/4, 3/8, or 11/32. These are the divisions you'd find on a ruler (in the US). Floating point numbers are usually implemented as dyadic fractions.
You can put the dyadic fractions into a binary tree as follows:

(define *the-dyadic-fraction-tree* 1)

(define (dft/left node)
  (/ (- (* (numerator node) 2) 1)
     (* (denominator node) 2)))

(define (dft/right node)
  (/ (+ (* (numerator node) 2) 1)
     (* (denominator node) 2)))

1 ]=> *the-dyadic-fraction-tree*

;Value: 1

1 ]=> (dft/left *the-dyadic-fraction-tree*)

;Value: 1/2

1 ]=> (dft/right *the-dyadic-fraction-tree*)

;Value: 3/2

1 ]=> (dft/left (dft/left (dft/right (dft/left *the-dyadic-fraction-tree*))))

;Value: 9/16

The next question is, what happens if I use a path derived from the Calkin-Wilf tree and use it on the dyadic fraction tree? Yes, this is a fairly random thing to try, but the trees are the same (that is, have the same structure) even if the values at the nodes are not. Either set of fractions is in a one-to-one mapping with the tree, so there is a one-to-one mapping between rational numbers and dyadic fractions.

This is Minkowski's ? (question mark) function. It maps the rational numbers on the X axis to the dyadic fraction on the Y axis. It has a number of weird properties. For example, it is strictly increasing and continuous, but it is not absolutely continuous. The function does not have a derivative in the usual sense.

A textbook case

This is priceless. The url is http://imgur.com/tDSX24E. It says

This is in a high school math textbook in Texas.

Example 2 : Is there a one-to-one and onto correspondence between integers and rational numbers? ... -4 -3 -2 -1 0 1 2 3 4 ... ... -1/4 -1/3 -1/2 -1/1 -2/4 -2/3 -2/2 -2/1 -3/4 -3/3 -3/2 -3/1 No matter how you arrange the sets, there is never one unique integer for each rational number. There is not a one-to-one and onto correspondence.

---

The challenge: implement rational->integer that returns "one unique integer for each rational number", and its inverse, integer->rational.

Currying and confusion

A friend of mine recently said to me "Don't know anything about currying except for food". I'm sure that nearly everyone who reads this blog is familiar with currying functions (and has probably curried a function within the last few hours), but it makes a good blog topic anyway.

"Currying" a function (as opposed to an entrée) is named after Haskell Curry, but he was inspired by a paper by Moses Schönfinkel, and it appears that Gottlob Frege came up with idea. So much for the name.

Pick your favorite binary function. I like "multiply", but any binary function will do. It need not be associative or commutative. As an example, imagine a print-to function that takes a document and a device.

Now consider these unary functions:

(define (multiply-by-five n) (* 5 n))
(define (multiply-by-negative-one n) (* -1 n))
(define (multiply-by-thirty-seven n) (* 37 n))
etc.

There is an obvious pattern here that we can abstract out:

(define (make-multiply-by x) (lambda (n) (* x n))

(define multiply-by-five (make-multiply-by 5))
(define multiply-by-negative-one (make-multiply-by -1))
(define multiply-by-thirty-seven (make-multiply-by 37))

Or, if we use print-to

(define (make-print-to-device device) (lambda (doc) (print-to doc device)))

(define print-to-inkjet (make-print-to-device the-inkjet-printer))
(define print-to-daisywheel (make-print-to-device the-daisy-wheel-printer))

Note the similarity between make-multiply-by and make-print-to-device.

(define (make-multiply-by x) (lambda (n) (* x n))
(define (make-print-to-device device) (lambda (doc) (print-to doc device)))

(define (make-<unary> an-argument)
  (lambda (another-argument) (<binary> an-argument another-argument)))

We can abstract this operation:

(define (make-unary-maker binary-operation)
  (define (make-unary-operation an-argument)
    (lambda (other-argument)
      (binary-operation an-argument other-argument)))
  make-unary-operation)

We have a pile of functions, all similar because they were created with the make- function. And the make- functions are all similar because they were created with make-unary-maker.

This is the meta-operation we call "currying". We take a function of several arguments, decide on some fixed values for some subset of arguments, and return a new function of the remaining, unfixed arguments.

---

Captain Obvious has a few things to say about functions.

A function won't return a value unless you call it.

The cardinality of the set of return values cannot be larger than the cardinality of the set of valid arguments. Functions don't make up values. There can be fewer possible return values than possible valid arguments, but never more.

If two sets are different sizes, then if you try to pair up elements from each set, you'll have some left over.

If we compose two functions and the set of possible output from the first is different in size from the set of possible valid input to the second, then there will either be output from the first that the second cannot accept, or possible inputs to the second that the first cannot produce. The latter case isn't a problem, but in the former case, you'll get an error.

If you "invert" a function, its output becomes input and vice versa.

Since inverse functions (when they exist, that is) are functions, The same constraints on the cardinality of input and output apply, but in the other direction. If you invert a composition, then the valid arguments and results swap roles. If the cardinalities match, we're cool, but if there is a mismatch, we're in the opposite situation of the non-inverse, and what was not an error (i.e. accepting an input that can never be produced) becomes one (i.e. producing a value that cannot be accepted).

If the appropriate cardinalities never increase, you can compose two functions. If the appropriate cardinalities match, you can invert the composition.

---

What does this have to do with currying?

I'm reading a bit about group theory and I chanced upon some exposition that was thick with comments and arguments surrounding the cardinality of various sets and subsets within the group. Once the exposition introduced a number of sets with equal cardinalities, the author pointed out that you can define an invertible function between any two sets with the same cardinality. Since the set of integers and the set of multiply-by- functions have the same cardinality, you could define a function that maps integers to multiply-by- functions. This isn't much of a revelation because that's what we did to get those functions in the first place. After slogging through this, I finally worked out what they were trying to say: "You can curry the binary operation." But they completely avoided the word "curry", and made arguments about set cardinality.

The next paragraph in the exposition was worse. You can curry the other argument to the binary operation. This leads to a completely different set of curried functions over a different set of arguments. Now it ultimately doesn't matter which argument you curry, when all is said and done, both arguments show up and are passed to the binary function. But all the intermediate values are drawn from different sets than if you curried the other way. It takes a very complicated argument about cardinalities to derive that either way of currying will work.

Something less abstract

To recap, a linear fractional transform (lft) is a function of this form:

(lambda (x)
  (/ (+ (* A x) B)
     (+ (* C x) D)))

Using an object-oriented extension to Scheme, I've given these functions a special print method so they look pretty. (I'm using MIT/Gnu Scheme with SOS, a CLOS-like object system. I haven't included all the code in my posts, but it should be easy for an intermediate programmer to replicate this.)

1 ]=> (make-linear-fractional-transform 1 2 3 4)

;Value 13: #[linear-fractional-transform 13 (x + 2)/(3x + 4)]

But otherwise they are normal one argument functions. Since they take and return objects of the same type, you can chain them together with function composition.

(define foo (make-linear-fractional-transform 1 2 3 4))
(define bar (make-linear-fractional-transform 1 4 0 2))

1 ]=> (compose foo bar)

;Value 14: #[compiled-closure 14 (lambda "global" #x35) #x84 #x25de16c #x3a5ea78]

1 ]=> ((compose foo bar) 1/10)

;Value: 81/203

lft/compose knows to explicitly construct an lft rather than a generic closure,

1 ]=> (lft/compose foo bar)

;Value 15: #[linear-fractional-transform 15 (x + 8)/(3x + 20)]

1 ]=> ((lft/compose foo bar) 1/10)

;Value: 81/203

but there is no difference when you apply the composed function.

I've been distracted by group theory, so let's do something less abstract. Here is an infinite stream of linear fractional transforms.

(define lft-stream
  (cons-stream
   (make-linear-fractional-transform 0 4 1 0)
   (stream-map
    (lambda (n d)
      (make-linear-fractional-transform n 1 d 0))
    (odds)
    (cons-stream 1 (stream-map square (naturals))))))

1 ]=> (pp (stream-head lft-stream 10))
(#[linear-fractional-transform 17 4/x]
 #[linear-fractional-transform 26 (x + 1)/x]
 #[linear-fractional-transform 25 (3x + 1)/x]
 #[linear-fractional-transform 34 (5x + 1)/4x]
 #[linear-fractional-transform 33 (7x + 1)/9x]
 #[linear-fractional-transform 32 (9x + 1)/16x]
 #[linear-fractional-transform 31 (11x + 1)/25x]
 #[linear-fractional-transform 30 (13x + 1)/36x]
 #[linear-fractional-transform 29 (15x + 1)/49x]
 #[linear-fractional-transform 28 (17x + 1)/64x])

This is a stream of functions, so let's compose them. Imagine our function stream is (f g h ...). We want to generate the stream (f (compose f g) (compose f g h) ...)

1 ]=> (fold-stream lft/compose lft:identity lft-stream)

;Value 36: {#[linear-fractional-transform 37 x] ...}

1 ]=> (pp (stream-head (fold-stream lft/compose lft:identity lft-stream) 10))
(#[linear-fractional-transform 37 x]
 #[linear-fractional-transform 17 4/x]
 #[linear-fractional-transform 45 4x/(x + 1)]
 #[linear-fractional-transform 44 (12x + 4)/(4x + 1)]
 #[linear-fractional-transform 43 (19x + 3)/(6x + 1)]
 #[linear-fractional-transform 42 (160x + 19)/(51x + 6)]
 #[linear-fractional-transform 41 (1744x + 160)/(555x + 51)]
 #[linear-fractional-transform 40 (23184x + 1744)/(7380x + 555)]
 #[linear-fractional-transform 39 (10116x + 644)/(3220x + 205)]
 #[linear-fractional-transform 38 (183296x + 10116)/(58345x + 3220)])

Let's apply these transforms to a number. 42 is random enough.

1 ]=> (pp (stream-head 
     (stream-map 
      (lambda (transform) (transform 42.0)) 
      (fold-stream lft/compose lft:identity lft-stream))
     10))
  (42. 
   .09523809523809523
   3.9069767441860463
   3.0059171597633134
   3.16600790513834
   3.137337057728119
   3.1423312358203845
   3.141464985588458
   3.1416146775443905
   3.1415888593191537)

That looks familiar. We'll pick a point further on in the stream.

1 ]=> (stream-ref (fold-stream lft/compose lft:identity lft-stream) 20)

;Value 48: #[linear-fractional-transform 48 (2166457145737216x + 48501417558016)/(689604727481670x + 15438480702645)]

1 ]=> ((stream-ref (fold-stream lft/compose lft:identity lft-stream) 20) 1.0)

;Value: 3.1415926535898056

This stream of functions are functions that approximate π. Each element in the stream is a function that is a better approximation of π than the last.

1 ]=> (stream-ref (fold-stream lft/compose lft:identity lft-stream) 10)

;Value 51: #[linear-fractional-transform 51 (3763456x + 183296)/(1197945x + 58345)]

1 ]=> (define foo (stream-ref (fold-stream lft/compose lft:identity lft-stream) 10))

1 ]=> (foo 0.0)

;Value: 3.1415888250921244

1 ]=> (foo 1.0)

;Value: 3.141593103503172

1 ]=> (exact->inexact (foo 'infinity))

;Value: 3.1415933118799275

Ok, let's subtract off that three.

(define (make-lft-subtract n)
  (make-linear-fractional-transform 1 (- n) 0 1))

1 ]=> (lft/compose (make-lft-subtract 3) foo)

;Value 54: #[linear-fractional-transform 54 (169621x + 8261)/(1197945x + 58345)]

1 ]=> ((lft/compose (make-lft-subtract 3) foo) 42.)

;Value: .14159330668296408

And multiply by ten

(define (make-lft-multiply n)
  (make-linear-fractional-transform n 0 0 1))

1 ]=> (lft/compose (make-lft-multiply 10) (make-lft-subtract 3) foo)

;Value 55: #[linear-fractional-transform 55 (339242x + 16522)/(239589x + 11669)]

1 ]=> ((lft/compose (make-lft-multiply 10) (make-lft-subtract 3) foo) 42.)

;Value: 1.4159330668296406

and now the next decimal digit is in the ones place. If we repeat this process, we'll generate the decimal digits of our approximation to π one at a time.

1 ]=> ((lft/compose (make-lft-multiply 10) (make-lft-subtract 3) foo) 42.)

;Value: 1.4159330668296406

1 ]=> ((lft/compose (make-lft-multiply 10)
                    (make-lft-subtract 1)
                    (make-lft-multiply 10)
                    (make-lft-subtract 3) foo) 42.)

;Value: 4.159330668296407

1 ]=> ((lft/compose (make-lft-multiply 10)
                    (make-lft-subtract 4)
                    (make-lft-multiply 10)
                    (make-lft-subtract 1)
                    (make-lft-multiply 10)
                    (make-lft-subtract 3) foo) 42.)

;Value: 1.5933066829640692

We'll define a composition as a data structure with an lft and a promise (delayed evaluation) to provide more lfts to compose. Operations on a composition will work on the lft and try to avoid forcing the promise if possible. If the lft isn't a good enough approximation, we force the promise, get a new lft, and create a new composition.

(define-class (<lft-composition>
        (constructor make-lft-composition 
                            (linear-fractional-transform promise)))
  ()
  (linear-fractional-transform define accessor accessor composition-transform)
  (promise                     define accessor accessor composition-promise))

(define (composition/refine lft-composition)
  (let ((lft (composition-transform lft-composition))
        (lft-stream (force (composition-promise lft-composition))))
    (make-lft-composition (lft/compose lft (head lft-stream)) (cdr lft-stream))))

Getting the integer part of a composition involves getting the integer part of the transform, if possible, and if not, refining the composition until it is possible. composition/integer-and-fraction returns two values: the integer extracted, and the final composition remaining after any refinement steps and after subtracting off the integer.

(define (composition/%integer-part lft-composition)
  (lft/integer-part (composition-transform lft-composition)))

(define (composition/integer-and-fraction lft-composition)
  (let ((integer (composition/%integer-part lft-composition)))
    (if (not integer)
        (composition/integer-and-fraction (composition/refine lft-composition))
        (values integer (make-lft-composition
                         (lft/compose 
                            (make-linear-fractional-transform 1 (- integer) 0 1)
                            (composition-transform lft-composition))
                         (composition-promise lft-composition))))))

We'll create a stream of digits by peeling off the integer parts and multiplying the fraction by 10.

(define (composition->digit-stream c radix)
  (receive (ipart fpart) (composition/integer-and-fraction c)
    (cons-stream 
     ipart 
     (composition->digit-stream
      (make-lft-composition
       (lft/compose (make-linear-fractional-transform radix 0 0 1)
      (composition-transform (force fpart)))
       (composition-promise (force fpart)))
      radix))))

(define foocomp (make-lft-composition (head lft-stream) (cdr lft-stream)))

1 ]=> (stream-head (composition->digit-stream foocomp 10) 50)

;Value 204: (3 1 4 1 5 9 2 6 5 3 5 8 9 7 9 3 2 3 8 4 6 2 6 4 3 3 8 3 2 7 9 5 0 2 8 8 4 1 9 7 1 6 9 3 9 9 3 7 5 1)

Flipping the number line

Here is a number line. It looks funny because I did two things: first, I stretched and compressed different parts so I could fit all the numbers in. Second, I rolled it into a circle.

All the numbers are in there. Most of them are really crammed near the top.

What happens to the number line if you apply a linear fractional transform? Here I apply 1/x:

1/x flips the entire diagram vertically. (What happens if you apply -x?)

In this example, I apply x + 1:

The numbers in the circle rotate clockwise, but not exactly linearly. The negative numbers are compressed and the positive numbers expanded.

There are, of course, an unlimited number of linear fractional transformations I can apply to my number line, but they all act by flipping or rotating the number line in some manner. We can derive every transform by a sequence of flips and rotates. Can we find a basic set of operations from which we can construct all the rest?

We can enumerate all the linear fractional transforms the same way we enumerate all rational numbers. We conceptually make a table of all of them, but we walk the table diagonally. There are four integers per linear fractional transform, rather than the two in a rational number, but that makes only a minor difference. It seems likely that our basic set of operations will be found in the set of linear fractional transforms with small integers for coefficients. If we use the positive and negative integers with absolute value of four or less, we'll have a good sized set to start with. Just from the combinatorics, we'd expect (expt 9 4) possible ways to list them, but many of these are duplicates. Nonetheless, there are still 2736 unique linear fractional transforms with coefficients that have an absolute value of four or less.

(I'm tired of typing "linear fractional transform". I'm going to abbreviate it lft.)

What happens if you compose an lft with itself? Obviously constants and identity remain unchanged, but the other ones are more interesting. These lfts are self inverses. If you compose them with themselves, you get the identity:

#[linear-fractional-transform 20 -x]
#[linear-fractional-transform 14 1/x]
#[linear-fractional-transform 13 -1/x]
#[linear-fractional-transform 227 x/(x - 1)]
#[linear-fractional-transform 226 -x/(x + 1)]
#[linear-fractional-transform 225 -(x + 1)]
#[linear-fractional-transform 224 (1 - x)]
#[linear-fractional-transform 223 (x + 1)/(x - 1)]
#[linear-fractional-transform 222 (1 - x)/(x + 1)]

These ones are not identities when you compose them twice (lft/compose x x), but they are identities if you compose them three times (lft/compose x x x)

#[linear-fractional-transform 522 1/(1 - x)]
#[linear-fractional-transform 521 -1/(x + 1)]
#[linear-fractional-transform 397 -(x + 1)/x]
#[linear-fractional-transform 396 (x - 1)/x]

These two are examples that require four composes

#[linear-fractional-transform 528 (x + 1)/(1 - x)]
#[linear-fractional-transform 527 (x - 1)/(x + 1)]

And these examples require six:

#[linear-fractional-transform 59 (x + 1)/(2 - x)]
#[linear-fractional-transform 58 (x - 1)/(x + 2)]
#[linear-fractional-transform 57 (2x + 1)/(1 - x)]
#[linear-fractional-transform 56 (2x - 1)/(x + 1)]
#[linear-fractional-transform 55 1/(3 - 3x)]

I have not seen an lft that gives the identity when composed with itself five times, nor seven, or eight.

Some lfts never self compose to identities. Obviously repeated self composition of x+1 will never equal the identity.

Linear fractional transforms and permutations

So what does group theory have to do with linear fractional transforms? I'm glad you asked. The answer is pretty complicated, though.

There's no good place to start here, so let's just dive in. This function computes the cross-ratio of four numbers:

(define (cross-ratio A B C D)
  (/ (* (- a b) (- c d))
     (* (- a d) (- c b))))

1 ]=> (cross-ratio 1 9 5 13)

;Value: 4/3

There's a geometric interpretation of the cross ratio, but for the moment, just think of it as a function that takes any four numbers and produces a value.

1 ]=> (cross-ratio 3 2 1 9)

;Value: -4/3

1 ]=> (cross-ratio 5 1 8 4)

;Value: 16/7

The cross ratio is preserved under linear fractional transforms:

1 ]=> (define lft2 (make-linear-fractional-transform 3 -1 2 7))

;Value: lft2

1 ]=> (cross-ratio (lft2 3) (lft2 2) (lft2 1) (lft2 9))

;Value: -4/3

1 ]=> (cross-ratio (lft 3) (lft 2) (lft 1) (lft 9))

;Value: -4/3

1 ]=> (cross-ratio 3 2 1 9)

;Value: -4/3

The cross ratio is also preserved under some permutations of its arguments.

1 ]=> (cross-ratio 3 2 1 9)

;Value: -4/3

1 ]=> (cross-ratio 1 9 3 2)

;Value: -4/3

1 ]=> (cross-ratio 2 3 9 1)

;Value: -4/3

but not all

1 ]=> (cross-ratio 3 2 9 1)

;Value: 4/7

1 ]=> (cross-ratio 2 9 3 1)

;Value: 7/3

Now here's the interesting part. There's a linear fractional transform that will "undo" the permutation of the arguments to cross-ratio.

1 ]=> ((make-linear-fractional-transform -1 1 0 1) 7/3)

;Value: -4/3

1 ]=> ((make-linear-fractional-transform 1 0 1 -1) 4/7)

;Value: -4/3

There are 24 permutations of the argument list to cross-ratio, and here are the equivalence classes:

((-4/3 (1 9 3 2) (2 3 9 1) (3 2 1 9) (9 1 2 3))
 (-3/4 (1 2 3 9) (2 1 9 3) (3 9 1 2) (9 3 2 1))
  (3/7 (1 2 9 3) (2 1 3 9) (3 9 2 1) (9 3 1 2))
  (4/7 (1 9 2 3) (2 3 1 9) (3 2 9 1) (9 1 3 2))
  (7/4 (1 3 2 9) (2 9 1 3) (3 1 9 2) (9 2 3 1))
  (7/3 (1 3 9 2) (2 9 3 1) (3 1 2 9) (9 2 1 3)))

And these are the linear fractional transforms that permute among these:

(#[linear-fractional-transform 89333 x]
#[linear-fractional-transform 89334 1/x]
#[linear-fractional-transform 89335 (1 - x)]
#[linear-fractional-transform 89403 1/(1 - x)]
#[linear-fractional-transform 89370 x/(x - 1)]
#[linear-fractional-transform 89393 (x - 1)/x])

This is cool. Instead of thinking of a linear fractional transform as a continuous function that traces out a hyperbola, or as an interval on the real number line, we can think of a linear fractional transform as a function that computes a permutation.

---

Here is a function called d3 that is defined on the symbols a, b, c, d, e, and f.

1 ]=> (d3 'a 'b)

;Value: d

With six symbols, there's only 36 possible outcomes, so we can enumerate them:

    a b c d e f
a ((e d f b a c) 
b  (f e d c b a) 
c  (d f e a c b) 
d  (c a b f d e) 
e  (a b c d e f) 
f  (b c a e f d))

Here's a hack. The identity element can be seen to be 'e, so list that row and column first:

    e a b c d f
e ((e a b c d f) 
a  (a e d f b c) 
b  (b f e d c a) 
c  (c d f e a b) 
d  (d c a b f e) 
f  (f b c a e d))

Now the table headers are exactly the same as the first entries in the respective rows or columns, so you can do without them.

((e a b c d f) 
 (a e d f b c) 
 (b f e d c a) 
 (c d f e a b) 
 (d c a b f e) 
 (f b c a e d))

Another weird thing is that you can permute any rows but the first, or permute any columns but the first, and still have an equivalent table.

Anyway, d3 is the "symmetry group" of a triangle. Here the operations are

e = leave it alone
f = rotate clockwise 120 degrees
d = rotate counter-clockwise 120 degrees
a = hold vertex A in place, and flip the triangle over swapping vertex B and C
b = hold vertex B in place, and flip the triangle over swapping vertex A and C
c = hold vertex C in place, and flip the triangle over swapping vertex A and B

Now consider this,

e = #[linear-fractional-transform 89333 x]
f = #[linear-fractional-transform 89393 (x - 1)/x]
d = #[linear-fractional-transform 89403 1/(1 - x)]
a = #[linear-fractional-transform 89370 x/(x - 1)]
b = #[linear-fractional-transform 89335 (1 - x)]
c = #[linear-fractional-transform 89334 1/x]

Just as d3 permutes a triangle, these linear fractional transforms permute among the equivalence classes of cross ratios.

n-ary functions and argument lists

Suppose we have a binary operation Op2 = (lambda (left right) ...). If it is closed over some type, you can make expression trees.

(Op2 (Op2 <arg1> <arg2>)
     (Op2
          (Op2 <arg3> <arg4>)
          <arg5>))

If Op2 is associative as well, these are equal:

(Op2 (Op2 <arg1> <arg2>)
     (Op2 <arg3>
         (Op2 <arg4> <arg5>)))

(Op2 <arg1> 
     (Op2 (Op2 <arg2> <arg3>)
          (Op2 <arg4> <arg5>)))

This makes the tree structure irrelevant so long as the fringe of the tree stays in order, so we can flatten the tree by making an N-ary version of the binary operation:

(define (binary->nary Op2)
  (lambda (a1 a2 . an)
    (fold-left Op2 (Op2 a1 a2) an)))

((binary->n-ary Op2) <arg1> <arg2> <arg3> <arg4> <arg5> <arg6> etc.)

The value of this expression naturally depends upon the values of the arguments. Changing the argument list is highly likely to change the value of the entire expression. However, we can make certain changes to the argument list without changing the value of the entire expression. If we know what those changes are, we can manipulate the argument list before invoking the operator, and still get the same answer. Naturally, most of the changes we can make to the argument list depend on the specifics of Op2, but it turns out that some interesting changes are possible without knowing any specifics, only knowing a few high-level properties of Op2.

Obviously, the first thing we can do is reduce the argument list through evaluation. Simply replace the first two arguments with the value of (Op2 <arg1> <arg2>)

((binary->n-ary Op2) <arg1> <arg2> <arg3> <arg4> <arg5> <arg6> etc.) =

((binary->n-ary Op2) (Op2 <arg1> <arg2>) <arg3> <arg4> <arg5> <arg6> etc.) =

((binary->n-ary Op2) <result> <arg3> <arg4> <arg5> <arg6> etc.)

Since Op2 is associative, we can replace any 2 adjacent arguments with their combination.

Now suppose there is an identity element among the arguments we can give to Op2.

(Op2 <arg> id) = <arg>  and
(Op2 id <arg>) = <arg>

We can do this:

(define (binary->nary Op2)
  (lambda an
    (fold-left Op2 id an)))

(define Op (binary->nary Op2))

Which is cleaner than the original. We also get a new way to manipulate the argument list to Op. We can add the identity element anywhere we wish, or we can delete the identity element wherever we find one.

(Op <arg1> <arg2> <arg3> Id <arg4> <arg5> <arg6> etc.) =

(Op <arg1> <arg2> <arg3> <arg4> <arg5> <arg6> etc.) =

(Op <arg1> Id <arg2> <arg3> <arg4> <arg5> Id <arg6> etc.)

One more restriction. We want Op2 to be invertible. Suppose (Op2 <arg1> <arg2>) = <result>. Op2 is invertible if, given any two of <arg1>, <arg2>, and <result>, the third can be uniquely determined. If you have one <arg> and a <result>, you can run things backwards and get the other <arg>.

Requiring Op2 to be invertible has many consequences, some of them quite non-obvious. An obvious consequence, though, is that we can define inverse elements. If (Op2 <argA> <argB>) = Id, then we say that <argB> is the inverse of <argA> (and vice versa). We find the inverse of an argument by fixing the output as the identity element and running Op2 backwards to find the other argument.

This gives us the final way to manipulate the argument list. If an element appears next to its inverse, both can be removed:

(Op <arg1> <arg2> <arg3> (inverse <arg3>) <arg5> <arg6> etc.) =
(Op <arg1> <arg2> (Op2 <arg3> (inverse <arg3>)) <arg5> <arg6> etc.) =
(Op <arg1> <arg2> Id <arg5> <arg6> etc.) =
(Op <arg1> <arg2> <arg5> <arg6> etc.)

So here are all the restrictions on Op2:

• Closed over an set of arguments
• Associative
• Has an identity argument
• Invertible

If Op2 has these properties (and a lot of binary operations do), then we can define an n-ary Op and play with its argument list. If you do this, you might notice that it looks kind of familiar:

(op f g (inverse g) j id h) = (op f j id h) = (op f j h)

The argument list sort of looks like a function pipeline. The allowed manipulations of the argument list are compatible with a function pipeline, too. In fact, it could be a function pipeline if Op is the compose operator, and f, g, and h are appropriate invertible unary functions. But whatever it is, the point is that it looks enough like a function pipeline that we can pretend that it is.

Playing with linear fractional transforms

I wanted to play with continued fractions and linear fractional transforms, so I wrote some code to make it easier. A linear fractional transform (also called a homographic function or Mobius transform) is a function like this:

In MIT/GNU Scheme:

;; x => (3x + 1)/(x + 4)

1 ]=> (define foo (make-linear-fractional-transform 3 1 1 4))

;Value: foo

1 ]=> (foo 2)

;Value: 7/6

I used an entity object so, in addition to invoking it on a number, there are two more ways to manipulate a linear fractional transform:

;; A predicate
1 ]=> (linear-fractional-transform? foo)

;Value: #t

;; And a CPS accessor
1 ]=> (lft/spread-coefficients foo (lambda (A B C D) (list A B C D)))

;Value 307: (3 1 1 4)

I also added a print method:

1 ]=> foo

;Value 308: #[linear-fractional-transform 308 (3x + 1)/(x + 4)]

As I mentioned in a prior post, you can partly apply a linear fractional transform:

1 ]=> foo

;Value 308: #[linear-fractional-transform 308 (3x + 1)/(x + 4)]

1 ]=> (lft/partly-apply foo 2)

;Value 315: #[linear-fractional-transform 315 (7x + 3)/(6x + 1)]

Since I want to reason about applying a linear fractional transform to an argument, I wrote an abstraction for that:

;; Apply LFT foo to continued fraction phi.
1 ]=> (make-lft-application foo phi)

;Value 311: #[lft-application 311 (3x + 1)/(x + 4) {1 ...}]

So now we can write a procedure that takes an application, peels off the first term in the continued fraction, feeds it to the linear fractional transform, and creates a new application:

(define (application/step lft-application)
  (let ((lft (application-function lft-application))
 (cf  (application-continued-fraction lft-application)))
    (make-lft-application
     (lft/partly-apply lft (head cf))
     (tail cf))))

1 ]=> (define appl (make-lft-application lft/identity sqrt-two))

;Value: appl

1 ]=> appl

;Value 317: #[lft-application 317 x {1 2 2 2 ...}]

1 ]=> (application/step appl)

;Value 318: #[lft-application 318 (x + 1)/x {2 2 2 ...}]

1 ]=> (application/step (application/step appl))

;Value 319: #[lft-application 319 (3x + 1)/(2x + 1) {2 2 ...}]

1 ]=> (application/step (application/step (application/step appl)))

;Value 320: #[lft-application 320 (7x + 3)/(5x + 2) {2 ...}]

All these lft-application objects should be (numerically) equal.

In an earlier post I showed how a linear fractional transform can be partly evaluated by determining the integer-part of the transform. The integer-part of an application is the integer-part of the application-function. You can get the fractional part by subtracting the integer-part.

---

A digression

If you apply a linear fractional transform to zero, it's obvious the answer is B/D. On the other hand, if you apply a transform to a sufficiently large x, you'll get as close as you want to A/C.

If the denominator of a linear fractional transform is zero for some value of x, there should be a vertical asymptote at that point. That's the pole of the transform. The pole is at (- D)/C. The pole will be at zero if D is zero. It will be at a negative number if D and C are the same sign and at a positive number if D and C differ in sign.

If you take a linear fractional transform with a pole at a negative number, and you sweep the input from 0 up on to infinity, the output will vary smoothly and monotonically from B/D approaching A/C and staying between both values at all times.

1 ]=> lft1

;Value 675: #[linear-fractional-transform 675 (3x + 1)/(4x + 2)]

1 ]=> (lft1 0)

;Value: 1/2

1 ]=> (lft1 1000000)

;Value: 3000001/4000002

1 ]=> (exact->inexact 3000001/4000002)

;Value: .7499998750000625

(On the other hand, if the pole is at a positive number, as you sweep the input from 0 up to infinity, the output starts at B/D, but flees away from A/C until the input gets to the pole. Then the output approaches A/C, but from the opposite direction. In any case, if the pole is positive, then the output will vary from B/D and eventually approach A/C, but never being between them.)

We can represent intervals as linear fractional transforms. The endpoints of the interval are A/C and B/D.

To get the width of the interval, just subtract the endpoints: A/C - B/D = (A*D - B*C)/(C * D)

---

Imagine you are performing some calculation with continued fractions. Since there may be an infinite number of terms, the calculation will proceed incrementally, using up terms as needed and generating other terms. So you can think of a more complex calculation as a tree, where a node in the tree is a linear fractional transform and the continued fraction terms flow between the nodes.

When we do an application/step, we move a term from the continued fraction into the linear fractional transform. Now consider a term as an element of information. We've moved this information out of the continued fraction and into the linear fractional transform. The information is apparently "stored" in the linear fractional transform until it is extracted as an output term for the next stage in the computation. But if you think about it, the "format" of the information is different depending upon whether it is flowing between nodes, where it is a series of continued fraction terms, or if it is stored in a linear fractional transform, where it is encoded somehow in the values of the coefficients. The act of partly-evaluating a linear fractional transform is in effect "encoding" some information as a continued fraction term. Partly applying a linear fractional transform is in effect "decoding" the continued fraction term (presumably generated by an earlier computation). Why not change to a more efficient communication channel?

When a node sends information to another node, instead of converting the information to several continued fraction terms to be assembled at the other end, we'll send the information as a single linear fractional transform. It contains the desired information already in the right "format". (See Peter Potts's work.)

A digression

What happens if we compose two linear fractional transforms?

(compose (lambda (x)
           (/ (+ (* A x) B)
              (+ (* C x) D)))
         (lambda (y)
           (/ (+ (* p y) q)
              (+ (* r y) s))))

We get

(lambda (x)
   (/ (+ (* A (/ (+ (* p x) q)
                 (+ (* r x) s))) B)
      (+ (* C (/ (+ (* p x) q)
                 (+ (* r x) s))) D)))

Which, after some algebra, turns into this:

(lambda (x)
   (/ (+ (* (+ (* A p) (* B r)) x) (+ (* A q) (* B s)))
      (+ (* (+ (* C p) (* D r)) x) (+ (* C q) (* D s)))))

Which is equivalent to this:

(lambda (x)
  (let ((E (+ (* A p) (* B r)))
        (F (+ (* A q) (* B s)))
 (G (+ (* C p) (* D r)))
 (H (+ (* C q) (* D s))))

   (/ (+ (* E x) F)
      (+ (* G x) H))))

Which you can see is another linear fractional transform.

If we have a linear fractional transform

(lambda (x)
  (/ (+ (* A x) B)
     (+ (* C x) D)))

It's inverse (if it has one) is:

(lambda (x)
  (/ (+ (* D x) (- B))
     (+ (* (- C) x) A))))

Which is yet another linear fractional transform. These things are everywhere.

Let's see, if we have a binary operation binop that is

1. Closed over some set, i.e. given any two elements of the set, the operation applied to the elements produces another element of the set. In other words, binop takes two arguments, returns one value, and the type of both arguments and return value are the same.
2. Associative, i.e. (binop a (binop b c)) = (binop (binop a b) c)
3. Has an identity argument. A "left identity" is an argument such that (binop left-identity x) = x. A "right identity" is an argument such that (binop x right-identity) = x. An "identity" argument works as a left or right identity.
4. Is invertible, i.e. for any objects a and b, there is a unique object x such that (binop a x) = b and a unique object y such that (binop y b) = a
   
   

then we have a group.

The compose function is a binary operation. When you compose a linear fractional transform with another, you get a third linear fractional transform.

1 ]=> (define lft1 (make-linear-fractional-transform 3 1 4 2))

;Value: lft1

1 ]=> (define lft2 (make-linear-fractional-transform 5 1 1 0))

;Value: lft2

1 ]=> (lft/compose lft1 lft2)

;Value 662: #[linear-fractional-transform 662 (16x + 3)/(22x + 4)]

Linear fractional transforms are associative.

1 ]=> (define lft3 (make-linear-fractional-transform 7 2 1 3))

;Value: lft3

1 ]=> (lft/compose lft1 (lft/compose lft2 lft3))

;Value 663: #[linear-fractional-transform 663 (115x + 41)/(158x + 56)]

1 ]=> (lft/compose (lft/compose lft1 lft2) lft3)

;Value 664: #[linear-fractional-transform 664 (115x + 41)/(158x + 56)]

The linear fractional transform (make-linear-fractional-transform 1 0 0 1) is both a left and right identity when composed with another linear fractional transform.

1 ]=> (define lft/identity (make-linear-fractional-transform 1 0 0 1))

;Value: lft/identity

1 ]=> (lft/compose lft/identity lft1)

;Value 665: #[linear-fractional-transform 665 (3x + 1)/(4x + 2)]

1 ]=> (lft/compose lft1 lft/identity)

;Value 666: #[linear-fractional-transform 666 (3x + 1)/(4x + 2)]

Given lft1 and lft2, there is a unique linear fractional transform x such that (compose lft1 x) = lft2, and a unique linear fractional transform y such that (compose y lft1) = lft2. x should be (compose (inverse lft1) lft2), and y should be (compose lft2 (inverse lft1))

1 ]=> lft1

;Value 675: #[linear-fractional-transform 675 (3x + 1)/(4x + 2)]

1 ]=> lft2

;Value 687: #[linear-fractional-transform 687 (5x + 1)/x]

1 ]=> (define x (lft/compose (lft/inverse lft1) lft2)))

;Value: x

1 ]=> (lft/compose lft1 x)

;Value 690: #[linear-fractional-transform 690 (5x + 1)/x]

1 ]=> (define y (lft/compose lft2 (lft/inverse lft1)))

;Value: y

1 ]=> (lft/compose y lft1)

;Value 691: #[linear-fractional-transform 691 (5x + 1)/x]

It sure looks like linear fractional transforms form a group under function composition.
I guess it's time to learn a little group theory.