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.