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/203but 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.1415926535898056This 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.1415933118799275Ok, 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: .14159330668296408And 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.4159330668296406and 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)
No comments:
Post a Comment