Lambdas All the Way Down
2024-08-06
This project is based on the work of u/DaVinci103
Since its inception, the bruijn programming language had support for positive and negative integers.
One thing I love about lambda calculus encodings for numbers is that they are usually infinitely precise – only limited by your computer’s memory. Still, lazy (call-by-need) or optimal reducers can sometimes do magical optimizations automatically that would have to be hard-coded in other languages.
It’s natural to want to extend the integers to arbitrary-precision fractions, and later real and complex numbers. Encoding rational numbers, it turns out, is as easy as writing two integers as a (Church) pair. The relevant arithmetic operations are fairly trivial to implement as well.
However, for the real numbers I couldn’t come up with any elegant solution:
This is, until I saw a great post on reddit about this exact topic.
So, in this article I want to show how the new arbitrary-precision arithmetic of bruijn is implemented, including some “more mathematical” real/complex formulas and ways to approximate their result. I also want to use this opportunity to introduce you to coding in bruijn and lambda calculus in general.
Note that all of the definitions (including library imports) have a 1:1 correspondence to pure lambda calculus – bruijn is just syntactic sugar (syntax explanation). We’re literally crafting complex formulas from the bottom up. You can click on any identifier to jump to its definition.
Natural numbers in lambda calculus are most commonly encoded as Church numerals. As in the meme above, they can be extended to integers by writing them as a pair.
While the implementations are quite small and elegant, any number now takes space:
:import std/Pair .
# syntactic sugar with 'u' suffix for unary/Church encoding
:test ((+3u) : (+5u)) ([[1 (1 (1 0))]] : [[1 (1 (1 (1 (1 0))))]])
Instead, I use balanced ternary numbers. Theoretically one could use any base, but an investigation by Torben Mogensen showed that bigger bases result in much more complicated implementations without giving significant performance benefits.
The encoding of a number in arbitrary base :
Balanced ternary numbers provide a very elegant way to denote negative numbers by assigning one of its three digits to a negative multiplier. I use the second de Bruijn index as the negative trit, the zeroth index as the zero trit, and the first index as the positive trit:
# from left to right:
:test ((+2)) ([[[[2 (1 3)]]]])
:test ((-2)) ([[[[1 (2 3)]]]])
Now we only need space!
You can find details about definitions of the arithmetic operations
in std/Number/Ternary
and my Rosetta Code
entry. All of the implementations below import integer maths into
the N
namespace.
The rational numbers represent the ratio between two integers:
The obvious implementation is therefore to use a pair of two balanced ternary numbers! Since we have the restriction of a non-zero , we subtract from the denominator in the encoding:
:import std/Pair .
:import std/Math N
# new syntactic sugar with 'q' suffix!
:test ((+5.0q)) ((+5) : (+0))
:test ((+3.14q)) ((+157) : (+49))
# converts a balanced ternary number to a rational number
number→rational [0 : (+0)] ⧗ Number → Rational
:test (number→rational (+5)) ((+5.0q))
Now, to check equality of two rational numbers
and
,
we only need to compare the product of
’s
denominator and
’s
numerator and vice versa. In the following snippets I uncurry the pairs
using the thrush combinator &‣
.
:import std/List L
:import std/Combinator .
:import std/Logic/Binary .
# returns true if two rational numbers are equal
eq? &[[&[[N.eq? (N.mul 3 N.++0) (N.mul N.++2 1)]]]] ⧗ Rational → Rational → Boolean
…=?… eq?
:test (((+1) : (+3)) =? ((+2) : (+7))) (true)
:test ((+0.25q) =? (+0.25q)) (true)
# returns true if a rational number is greater than another
gt? &[[&[[N.gt? (N.mul 3 N.++0) (N.mul N.++2 1)]]]] ⧗ Rational → Rational → Boolean
Other basic operators can be implemented just as easily:
# negates a rational number
# this could also be done in the denominator
negate &[[N.-1 : 0]] ⧗ Rational → Rational
-‣ negate
# adds two rational numbers
add &[[&[[p : q]]]] ⧗ Rational → Rational → Rational
p N.add (N.mul 3 N.++0) (N.mul N.++2 1)
q N.add (N.mul 2 0) (N.add 2 0)
…+… add
# subtracts two rational numbers
sub &[[&[[p : q]]]] ⧗ Rational → Rational → Rational
p N.sub (N.mul 3 N.++0) (N.mul N.++2 1)
q N.add (N.mul 2 0) (N.add 2 0)
…-… sub
# multiplies two rational numbers
mul &[[&[[p : q]]]] ⧗ Rational → Rational → Rational
p N.mul 3 1
q N.add (N.mul 2 0) (N.add 2 0)
…⋅… mul
# calculates the multiplicative inverse of a rational number
invert &[[N.compare-case eq gt lt 1 (+0)]] ⧗ Rational → Rational
eq Ω
gt N.++0 : N.--1
lt N.-(N.++0) : N.--(N.-1)
~‣ invert
# divides two rational numbers
div [[1 ⋅ ~0]] ⧗ Rational → Rational → Rational
…/… div
# approximates a rational number to n digits
approximate &[[[&[[int : float]] (N.quot-rem 2 N.++1)]]] ⧗ Rational → Number → (List Number)
int 1
float N.number→list (N.div (N.mul 0 (N.pow (+10) 2)) N.++3)
This is nothing special really and appears similar in various
libraries like Haskell’s
Rational
(although somewhat obfuscated).
Such libraries often do a reduction step that shortens fractions similar to how we do intuitively. These steps require a (greatest common divisor) operation, which is currently quite expensive in bruijn1. Still, one could do something like this:
# finds the smallest equivalent representation of a rational number
compress &[[[(N.div 2 0) : N.--(N.div N.++1 0)] (N.gcd 1 N.++0)]] ⧗ Rational → Rational
%‣ compress
# equivalent of Haskell's Fractional reduce function
# similar to compress, but without uncurrying
reduce [[[(N.div 2 0) : N.--(N.div 1 0)] (N.gcd 1 0)]] ⧗ Number → Number → Rational
# shortened (and much slower) add operation
add' &[[&[[reduce p q]]]] ⧗ Rational → Rational → Rational
p N.add (N.mul 3 N.++0) (N.mul 1 N.++2)
q N.mul N.++0 N.++2
All of this could already be used to do some pretty advanced calculations.
λ pow-n = [L.nth-iterate (mul 0) (+1.0q)]
λ n = L.nth-iterate [pow-n (0 + (+0.1q)) (+2)] (+4.2q) (+4)
2179006444280423418821440444946289062500000000000000/152587890625000000000000000000000000000000
(~14280336633.23618300, time: 5.7e-2s, length: 1018)λ approximate n (+5)
{(+14280336633), (+2), (+3), (+6), (+1), (+8)}
Note: All REPL examples ran on my ThinkPad T14G1 AMD using bruijn’s HigherOrder reduction mode (ergo Haskell).
The real magic, of course, happens with real numbers!
The missing insight I got from reddit is the fact that approximations of real numbers can be seen as the infinite limit of a function that maps natural to rational numbers. For example, for specific there exists , such that
This would typically be defined as a sequence that converges to the desired real number2.
Therefore we represent a real number as an abstraction over the rational numbers:
:import std/Pair .
:import std/Math/Rational Q
# if n is not irrational, just abstract the rational directly
# new syntactic sugar with 'r' suffix!
:test ((+5.0r)) ([(+5) : (+0)])
# converts a balanced ternary number to a real number
number→real [[Q.number→rational 1]] ⧗ Number → Real
:test (number→real (+5)) ((+5.0r))
Approximating a real number to an arbitrarily precise rational number can then be done by applying some natural number. In theory we could use any base here, yet I still use balanced ternary. This wastes some space for negative numbers, but is best supported in bruijn’s standard library.
:test ((+4.2r) (+42)) ((+4.2q))
# returns true if two real numbers are approximately equal
# one could further approximate by comparing the difference with some threshold
approx-eq? [[[Q.eq? (1 2) (0 2)]]] ⧗ Number → Real → Real → Boolean
# initialized with some large enough value
…≈?… approx-eq? (+42) ⧗ Real → Real → Boolean
:test ((+4.2r) ≈? (+4.2r)) (true)
This is one of the drawbacks of the implementation: It’s not possible to be 100% certain that two real numbers are equivalent! This is okay though, usecases for this are rare anyway.
For other binary operators I use a small combinator trick that “joins”/“lifts” the argument of both sides such that, inductively, we only ever need to apply a single number to any real number to get its approximation:
:import std/Combinator .
# to remind you:
# phoenix combinator / liftM2
φ [[[[3 (2 0) (1 0)]]]] ⧗ (b → c → d) → (a → b) → (a → c) → (a → d)
# φ (Rational → Rational → Rational) (Number → Rational) (Number → Rational)
# ~> (Number → Rational) ~> Real
# adds two real numbers
add φ Q.add ⧗ Real → Real → Real
…+… add
# subtracts two real numbers
sub φ Q.sub ⧗ Real → Real → Real
…-… sub
# multiplies two real numbers
mul φ Q.mul ⧗ Real → Real → Real
…⋅… mul
# divides two real numbers
div φ Q.div ⧗ Real → Real → Real
…/… div
We can do something similar for the missing unary operators:
# to remind you:
# bluebird combinator / composition
b [[[2 (1 0)]]] ⧗ (b → c) → (a → b) → (a → c)
# b (Rational → Rational) (Number → Rational)
# ~> (Number → Rational) ~> Real
# negates a real number
negate b Q.negate ⧗ Real → Real
-‣ negate
# inverts a real number
invert b Q.invert ⧗ Real → Real
~‣ invert
# finds smallest equivalent representation of a real number
compress b Q.compress ⧗ Real → Real
%‣ compress
Like with the rationals before, we can now go calculate some things!
:import std/List L
# power function: real^integer
pow-n [L.nth-iterate (mul 0) (+1.0r)] ⧗ Real → Number → Real
...
This is boring though. How can we calculate ? The trick is that for positive ,
Therefore we only need to implement exp
and ln
.
# simple exp using infinite limit
exp-lim [[pow-n ++(1 / (number→real 0)) 0]] ⧗ Real → Real
exp-lim
converges very slowly. Since we want to combine exp
with other
(faster converging) formulas, we need to use a different technique using
Taylor series. This is not as trivial:
# exp using Taylor series
# the factorial and power is embedded into iteration
# such that we don't always calculate them redundantly
exp [[L.nth-iterate &[[[[op]]]] start 0 [[[[1]]]] 0]] ⧗ Real → Real
start [0 (+1.0r) (+1.0r) (+1.0r) (+1)]
op [[[2 1 0 (4 + (1 / 0)) N.++3]] pow fac]
pow 6 ⋅ 4
fac 3 ⋅ (number→real 1)
The iteration can be described by a recursive sequence:
The [[[[1]]]]
then extracts the third value,
,
of the 4-tuple start
.
ln
is
similarly complicated:
# natural logarithm using Taylor series
# again, the power is embedded into iteration (see start)
ln [[[L.nth-iterate &[[[op]]] start 1] (--1 / ++1 0) [[[1]]]]] ⧗ Real → Real
start [0 1 (+0.0q) (+0)]
op [0 pow (Q.add 2 go) N.++1]
pow Q.mul 3 (Q.mul 4 4)
go Q.mul ((+2) : (N.mul (+2) 1)) 3
# real power function: real^real
pow [[exp (0 ⋅ (ln 1))]] ⧗ Real → Real → Real
…**… pow
You can already see a common problem: Different approximations converge with different speed (here notated as error in iteration step ). Combining them into different terms can yield unbalanced ratios of convergence, in turn impacting performance of more precise calculations. I could not find a solution yet, so I just try to make most formulas converge roughly with the same speed.
Some example REPL usage, with precision :
λ exp (+1.0r) (+10)
18099969098565397826764800000/6658606584104736522240000000
(~2.7182818011463845, error of ~1.73e-8, time: ~0.03s, 632bit)λ ln (+2.0r) (+10)
233890488187863915541701543725460185763128364510093286440/337432647424663101659020986764900182676573608430256879075
(~0.6931471805498117, error of ~1.01e-11, time: ~0.03s, length: 1051bit)λ ((+2.0r) ** (+5.0r)) (+10)
598079134110276147827564883538659081286604734987944789387834372340947075110943453837234618373552693185584716635684189972636644241652171376916943525486461022728058976842778463558150783565917373041554240375907383228706705686070473570707758765671846250321276278714131843233090967960950100066733563643600273605634465180162835902914836231104763335304562987169953913286974135123478137546638315618992875001966578776051314748446370950239017032303780045954608544272508524030682516609446564296943397047151770455802591948913236771550241140628206362089337306080729456884605657473736006900750062127183174575149591063756525402897186089732707298165462679140041446762061630523298302176019472132682528624760986533490881586211470148319049682387677590759715636301654414721058038524158625347730839117574417904657884025887070370668002870317851865610637261303431591432712751211868675681918441438181009236310163455557218485924474140212254727259954387270719434271208143299910616186908085375881165357827436854844415894159309606791014794031897082269022230707092560145786402326155316644252639477330086223550047352406410530249264625980127611616014612744110563965082672579163660755936886822130116400763646980101508854615001526358232971088757782801868634665475194968733501724384816076327046316070654670769403132424875856535131094386148667653955301369622336790381866119818316339256930323110721023022331415394093223237411495404394416682243161354096950645896692264278922066316269661717458816226041291616605553523205391825489419079585537561636737655730488488440777060126875525392906188334680281882510800224189036509975301509573328453196523002298919394859926114954914787468282500028398858785002799476737665131442186918408368665754947375031323023994205757187950886526437405495638105657127387780793248754925231819696930561066681593613019051902903759963710577600178126551532374230654018917071744634056944518994438021779780944427419530474781506626479456134411453643022372837732547369414691303401174003771686180236856810234420872033503572191499213696774228469644141715532392472118731282252789583127278677420087233498337120738433880137357508449801704906102451094822409824236222734766813101302782137436899378522610042953094252716889434203754695406049119679351289357796953897381003734393866725685601294639966545608927832846146076337329964777111812830183537916681463188543158146621245507172569192113334571731610610293877626003372555010124576440160375294328387964074220385933549089782641666077677368723655537628973398985260911246931635785708900043886771153987974930564246498556798449367380515140417080137186355779417955734278686914202430144848198772372684160394117644408417711866459815541733128505818698901992840753593067270178536796957546332531748299524902599475784658804350813184544953176440746444798346737996269153523375711529045417264644913503743443882289016496822187252992674242720739998853511069087517775444431645507009913242118489351376976411894580661365924273926755289179739713122277279765888522997890201654099882914400318594082414541185666356171548192583342325223014342107882723075838412735160971173131727596228301990777254104614257812500000000000000000000000000000000000000/74760779145070733867908932761413155762971594365211605875643688715795658071296489021978852498729185241962654371126198985856785158579859308675958564026242002487212659492193987786684537196431214280757975925264758942883634825715655210754259881854161850075555156727957799415478911599087787338759644279029501022468741585284712250870563608807772216626319992978677936038223391819445771394867916047740401921773908140195393144179369141320178364174064953698343981949644450414632898184304796423731597697465375566339947566525590816014253540526958504467447037286922629637764004487538446911081259943320301928694733654797538871623891971353728913964342443511026866560197170972549628045705440407603312893994442950853177415705720888410553006515572923881893088120078667252896080200863364262936348551786574482904929604768379106016250412338548222219088828801332785047621040579051084935367004882760987761820574788206148525328752068539097310214079835762720917029086377377340872933055196076185798222540307196361313399629934793229668735088328473360784158080026141115666746992410622458258295577175907820466408926620132964962486225515990672744837581085043326597011264706252949090969847322322495952154924286412491550689193904192548701376743857571103786655274607787837390201233431227051287558269947139643643354748534743176311469903335733573895699231869252769229007514872672526870460166704215242799424899976518534714924877197462712550853495447879709422302180307953400983205499497427437609334221946874236442601011590334150274580045185317079861788907309774171436652200887227409217065183012520489092932128741216386785887485237903535011967264310628081772805865418646841771828566494422331918960006779437847077838135379228767020888503346221487625606989561464766098979814616703315507210693772610145034125230859983749588645164959428343521078143981388069312327590834893754353340597878304329885934981702091972861233235112086702432657823805894511173586979592324513656316233277138188041726134730813327170727412918157800292547137553599968223539455477182881451302345840397262591228570188639075069091419714621221515228131751997428030553750909780200688042650201197354138784925439398522202355311797657482825417733561554769273059694184443577163344456922849979375864932918277480426922461006181579979595175551025678507687767005418040140775325648296021675347480784155313735196658541390475382786035917854927202267064655575298997606246563002418114754520547111452824170514469941660203289142179436597789664589790956742643067949813175447395938470322035555857103944544604786924840964930980034735034061609181845620607831460713770653133250923402551987297201485964824552332216241734485472613124394181347966082706956069628488992962113147532752637402860615037429791398408738208746234630119201878416063821300338439170456002649326583102180899337539911353548588857412013494440755404189179361325795357642589917024984308672607471573171033808418901434055060528412896195791352491112725007239345428334121420462709200860489873502414825861884547983847540875126400700037685065643942537861409047019979211343772995432655685736347752224500995410494397219736129045486450195312500000000000000000000000000000000000000 (~7.9999050431205925, error of ~9.5e-5, time: ~49.2s, length: 54kbit)
The numbers getting that large is the implication of missing reduction of rationals as explained in the previous section, so yes, this definitely needs some improvements. The calculations then also get incredibly slow once multiple formulas are involved. This could potentially also be improved with better reducers (currently almost nothing is shared/parallelized).
Or, one could go a different route and just find a better formula
that converges faster. For example, for sqrt
and cbrt
:
# square root via x^0.5
pow-sqrt \pow (+0.5r) ⧗ Real → Real
# real square root
# Newton's/Heron's method, quadratic convergence
sqrt [[y [[[N.=?0 guess go]]] (1 0) 0]] ⧗ Real → Real
guess (+1.0q)
go [Q.div (Q.add 0 (Q.div 2 0)) (+2.0q)] (2 1 N.--0)
# real cube root
# Newton's/Heron's method, quadratic convergence
cbrt [[y [[[N.=?0 guess go]]] (1 0) 0]] ⧗ Real → Real
guess (+1.0q)
go [Q.div (Q.add (Q.mul (+3.0q) 0) (Q.div 2 (Q.mul 0 0))) (+4.0q)] (2 1 N.--0)
Comparison:
λ pow-sqrt (+2.0r) (+7)
89033367892848420325715129500331038667412886056368323573159431977310238895439218094980112725702840342993272134600096124077626686964733473513621053153947308252557325158480398438676319722582712899722911291367231026263483599926951994934434361052597331385548045264275022690705284987673823814184868974677551776642419208825170717874123852968703476348321320937471629780426469759321294332061735019180158114422592282352023648657975088457689473927068241499703325954672642420561944499393423107013113058488991800090967348390600054350006070192488093414571502040824332030927892137510876992387627104157234890750860521744152469018399002283691258127609527505782738117052145118472073066393963385416112881333398245475247711810954299587984415353368400124041917003009441829548789670468653120761682572674942660955340800000000000000000000000/62956098752232057387808350112841270755951480654727337171842489892946126359947885037380093400728241011317082442072604878261728100965154427911444757570897439621083387621109809746375368978651533367238019668081714567583330572204727730673367194112699820403536457063036501180427572233579383441923086911895434098875666345475230727812250574821804684019255992429473998197711015819847678940543489342616860950975760163509063079641664150888593726203230408507924318543814430554711334286943508216430250098869852834437753281169571221961391204635543888995125143434949553386638535183467482640471944118689238928081592507454317778459894864002701265511068264849228331714462408887817411373838522061660807386673012997199144341243681258637527781502269807089280818973625945249329852752982953473633532574022778880000000000000000000000000000000
(~1.4142135497189112, error: ~1.26e-8, time: ~3.6s, length: 14kbit)λ sqrt (+2.0r) (+7)
4946041176255201878775086487573351061418968498177/3497379255757941172020851852070562919437964212608 (~1.4142135623730951, error: ~2.89e-98, time: ~0.069s, length: 1071bit)
For and the previous method takes longer and is less precise! You can see how different formulas result in vastly different performance. This way, I believe, one could probably craft reasonably efficient implementations of many other mathematical functions.
Also fun:
# returns the derivative of a function
# (so much in that excellent formula!!1)
derive [[[[((3 (2 + 0)) - (3 2)) / 0] ((+1.0r) / conv) 0]]] ⧗ (Real → Real) → (Real → Real)
conv number→real (N.pow (+2) 0) # custom convergence ratio
# e.g. derive sqrt (+3.0r) n, n→∞
In REPL:
λ derive sqrt (+3.0r) (+5)
4548499328554333130253458031438461086963398946474394806600894005568776844296474084125321695485140926464/15756464483850753172492714095995675761947988863890026857797915483539179836951188627907703817126616760320 (~0.2886751233575412, error: ~1.12e-8, time: ~0.6s, length: 2201bit)
The results could be improved by using a better adapted convergence
ratio in derive
,
although this would also imply much longer results.
I’ve already shown you how to calculate using the Taylor series of . Other interesting constants:
π/2 [L.nth-iterate &[[[[[op]]]]] start 0 [[[[[3]]]]]] ⧗ Real
start [0 (+1) (+0.0q) (+1) (+1) (+1)]
op [0 N.++5 (Q.add 4 ((N.mul 3 2) : N.--1)) num-pow num-fac denom]
num-pow N.mul 3 (+2)
num-fac N.mul 2 (N.mul 5 5)
denom [N.mul 2 (N.mul 0 N.++0)] (N.mul (+2) 5)
# ratio of circle's circumference to its diameter
π π/2 ⋅ (+2.0r) ⧗ Real
# golden ratio from direct formula
φ ++(sqrt (+5.0r)) / (+2.0r) ⧗ Real
# conjugate golden ratio
ψ -(~φ) ⧗ Real
# golden ratio from fibonacci convergence
φ* ++[(L.nth-iterate &[[0 : (N.add 1 0)]] ((+0) : (+1)) 0) [[1 : N.--0]]] ⧗ Real
# real fibonacci becomes trivial
fib [((pow φ 0) - (pow ψ 0)) / (sqrt (+5.0r))] ⧗ Real → Real
Results:
λ π (+10)
58161346883519614707939305922472528176805953649166549778432000000000000000/18519311967783166520148611995625183290053887455979435458560000000000000000
(~3.140578169680337, error: ~0.001, time: ~0.141s, length: 1511bit)λ π (+20)
23770718233993033585696255811682815853064495117892501617611447251371134322557383888051869648991848387499657110718132777488467529319755142727319452657262872133986736812397012330767068148003469821569121725937860722210965191933852585434483132680685165347732376029802802822656597172962731748393261075522609795620777191496744960000000000000000000000000000000000000000000000000000000000000000000000000000000/7566456363918636918448056774706264206207563561731028263307703141167935514660420425628485969667487374826285184836250356538021409394395657113901240231996338685382563709519924768061981073010130735334419464639463771932114561342847927309380840121890195187487927033515312915490616963923270255231275266462271264256491368166195200000000000000000000000000000000000000000000000000000000000000000000000000000000
(~3.141591927675147, error: ~7.26e-7, time: ~8.3s, length: 8089bit)λ φ (+7)
138598980290857213752199295710131342640539745689910576085702017024/85658880625826545040363067218944712997310200964513375659560534016
(~1.618033988749895, error: 7.06e-54, time: ~0.17s, length: 1418bit)λ # another example of the convergence rate problem:
λ φ* (+500)
1206484255615496768210420703829205488386909032955899056732883572731058504300529011051/745648276861993189984204820755333242001144560869101973859680384188513887852280667476 (~1.618033988749895, error: 5.12e-85, time: ~0.17s, length: 1757bit)
Also see std/Math
for a list-based approximation of π
, it’s
actually more performant.
Some trigonometric functions (kinda arbitrarily chosen):
# hypotenuse
hypot [[sqrt ((0 ⋅ 0) + (1 ⋅ 1))]] ⧗ Real → Real → Real
# arctan by Taylor series
atan* [[[L.nth-iterate &[[[[op]]]] start 1] (1 0) [[[[3]]]]]] ⧗ Real → Real
start [0 1 [[0]] (Q.pow-n 1 (+3)) (+3.0q)]
op [0 ((3 Q.add Q.sub) 4 (Q.div 2 1)) \3 num denom]
num Q.mul 2 (Q.mul 5 5)
denom Q.add 1 (+2.0q)
# actual arctan for arbitrary x
atan [[[Q.lt? Q.|0 (+1.0) (atan* [1] 1) go] (1 0)]] ⧗ Real → Real
go Q.sub π/2-pm conj-atan
π/2-pm Q.mul (Q.>?0 (+1.0q) (-1.0q)) (π/2 1)
conj-atan atan* [Q.div (+1.0q) 1] 1
# 2-argument arctan
# this one's really weird
atan2 [[[[[Q.add a (Q.mul b (Q.add c d))]] (2 0) (1 0)]]] ⧗ Real → Real → Real
z (+0.0q)
a Q.=?0 z (atan [Q.div 2 1] 2)
b Q.sub (+1.0q) (Q.<?1 (+2.0q) z)
c Q.<?0 (π 2) z
d Q.=?0 (π/2 2) z
Formulas like sin
and cos
are similar
to atan
and are left as an exercise to the reader.
There are many ways to translate real approximations like Taylor
series to lambda calculus. Depending on your needs, my method using
L.nth-iterate
may not be the optimal approach.
If you use unary/Church numerals instead, you actually get repeated
iteration for free! A Church numeral is [[1 (1 (1 (... (1 0))))]]
,
so therefore we can just apply the iteration function as first argument
and the starting value as second. The previous exp
implementation could easily be translated:
# unary exp using Taylor series
unary-exp [[0 &[[[[[op]]]]] start [[[[1]]]]]] ⧗ Number → (Unary → Rational)
start [0 (+1) (+1) (+1.0q) (+1)]
op [[2 1 0 (Q.add 4 (1 : N.--0)) N.++3]] pow fac
pow N.mul 6 4
fac N.mul 3 1
Things like list indexing are also easier with Church numerals.
However, they also take more space, since they’re not only
,
but also because they literally duplicate the iteration function
times! They could still be a good option for rapidly converging series.
Also, if for whatever reason you’d want to golf an exp
function or
similar, you could use unary numerals at every step starting with
defining integers as a pair of two unary numerals. unary-exp
could
then get truly tiny.
I’ve demonstrated some conversions from Taylor series to recursive sequences. In theory, such recurrence relations could be resolved using a variadic fixed-point combinator.
In my experience, the use of this combinator doesn’t give performance
improvements over normal iteration using y
-recursion
(like I’ve used here) and can be quite confusing. Still, L.y*
can be fairly
small, so this could also be a good opportunity for golfing.
Complex numbers are, contrary to their name, quite simple. They are just a pair of real numbers!
Still, we need to lift the approximation arguments of the two real numbers to the outside, such that we only ever need to apply a single natural number:
:import std/Pair .
:import std/Combinator .
:import std/Math/Real R
:import std/Math/Rational Q
# new syntactic sugar with complex part as suffix!
:test ((+0.0+1.0i)) ([((+0.0r) 0) : ((+1.0r) 0)])
# imaginary unit
i (+0.0+1.0i)
# converts a balanced ternary number to a complex number
number→complex [[(R.number→real 1 0) : ((+0.0r) 0)]] ⧗ Number → Complex
:test (number→complex (+5)) ((+5.0+0.0i))
We can then extract the real/imaginary part using traditional pair selectors mixed with a lifted argument, therefore resembling a real number:
# returns real part of a complex number
real [[1 0 [[1]]]] ⧗ Complex → Real
:test (real (+5.0+2.0i)) ((+5.0r))
# returns imaginary part of a complex number
imag [[1 0 [[0]]]] ⧗ Complex → Real
:test (imag (+5.0+2.0i)) ((+2.0r))
Here we use φ
and b
again to lift
the approximator argument. After their application, the real numbers
inside the complex number will actually appear rational, such that we
can use rational operations directly.
# adds two complex numbers
add φ &[[&[[(Q.add 3 1) : (Q.add 2 0)]]]] ⧗ Complex → Complex → Complex
…+… add
# subtracts two complex numbers
sub φ &[[&[[(Q.sub 3 1) : (Q.sub 2 0)]]]] ⧗ Complex → Complex → Complex
…-… sub
# multiplies two complex numbers
mul φ &[[&[[p : q]]]] ⧗ Complex → Complex → Complex
p Q.sub (Q.mul 3 1) (Q.mul 2 0)
q Q.add (Q.mul 3 0) (Q.mul 2 1)
…⋅… mul
# divides two complex numbers
div φ &[[&[[p : q]]]] ⧗ Complex → Complex → Complex
p Q.div (Q.add (Q.mul 3 1) (Q.mul 2 0)) (Q.add (Q.mul 1 1) (Q.mul 0 0))
q Q.div (Q.sub (Q.mul 2 1) (Q.mul 3 0)) (Q.add (Q.mul 1 1) (Q.mul 0 0))
…/… div
# negates a complex number
negate b &[[(Q.negate 1) : (Q.negate 0)]] ⧗ Complex → Complex
-‣ negate
# inverts a complex number
invert b &[[p : q]] ⧗ Complex → Complex
p Q.div 1 (Q.add (Q.mul 1 1) (Q.mul 0 0))
q Q.div 0 (Q.add (Q.mul 1 1) (Q.mul 0 0))
~‣ invert
Several of the defined functions can now be combined to (attempt to) calculate :
# very similar to the real one!
exp [[L.nth-iterate &[[[[op]]]] start 0 [[[[1]]]] 0]] ⧗ Complex → Complex
start [0 (+1.0+0.0i) (+1.0+0.0i) (+1.0+0.0i) (+1)]
op [[[2 1 0 (4 + (1 / 0)) N.++3]] pow fac]
pow 6 ⋅ 4
fac 3 ⋅ (number→complex 1)
# complex natural logarithm
ln [[[p : q] (1 0)]] ⧗ Complex → Complex
p R.ln (&[[R.hypot [2] [1]]] 0) 1
q &[[\R.atan2 [2] [1]]] 0 1
# complex power function: complex^complex
pow [[exp (0 ⋅ (ln 1))]] ⧗ Complex → Complex → Complex
…**… pow
In REPL:
λ ln (+2.0+4.2i) (+3)
131333047804453819160134580326741441990719587929624714214939549619500790973864294006194984228668428249332105635341198746982341517036711701337460293632000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000/83484794802006666694721421339938768091461941633007787921780106775988824518415759807105367379415841350472645076302765692152680960498897197586158714880000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
: 110569850108556581649876480/108150992982294090821067600
(~1.57313734+1.02236556i, error: ~0.03+0.1i, time: ~1.5s, length: 6376bit)λ # ugh, numerator converges a lot faster than denominator, anyway...
λ (+2.0+0.0i) ** (+3.0+0.0i) (+2)
46279959818696076004546197821636492728892381317705253201507013980345318470977416794052934030379365364822245947734236637528515350739555014550295902286727588227228815359769649963518843577958400000000000000000000000000000000000000000000000000000000000000000000000000000000/8488066053046236477811349396282461894378738844356320153510331083656243451009448569006437421681264652074166946504077170970137877639991273678242184865776333754718651791926832629500105614950400000000000000000000000000000000000000000000000000000000000000000000000000000000
: 0/8488066053046236477811349396282461894378738844356320153510331083656243451009448569006437421681264652074166946504077170970137877639991273678242184865776333754718651791926832629500105614950400000000000000000000000000000000000000000000000000000000000000000000000000000000
(~5.45235623+0i, error: ~2.58+0i, time: ~2.9s, length: 10740bit)λ # yes, many abstractions make everything worse. Looking at you, Electron!
λ exp ([(R.π 0) : (+0.0)] ⋅ i) (+4)
-67479993202306472714667916921031425556099775398952539611521667570250447344805122008169714475335413570680821954970430751292829585907152642903675057639798448289757605848159425094599290482701423963369343202081203932066765475487409153721693127606181916420307244302193597548169505958661717082177471539283605044618375306864059617331387590135051834499754896430296524931141329329805419250287631708024288329733458404003715292655399963717259776754256684136966566721257811784729486080081920000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000/1361884868862060097549690181333056765046439560624409926616798792493328435787994562240567391560867539517753322628725069314972581967222209672892486307189350239216889094360499479590034032207920263636241298759966234243205570905766386207577079921475443490332191343090551009522870423225966591704205283036225548377573722514188921626990790748596113579717060293956259137133760787580672566454929838687772488012891423833655251357457695563011145177240090717870704091720599864440759816908963840000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
-2274464885548702030965181024866063200451887687181824158478906668527172710670660032680499836893937077839821263721442724551514516189004446201713686954381666886213283445931943125399617719268887457535960704836138403659962512688606713459006461086148384909887833506799178160812769954205264005149736354349093526426930883081681818503846693688041596303454328101341483063852229800372285465338831102296823361184831086127545753265883348598917220466614124079804220921276889958249694255644672000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000/1361884868862060097549690181333056765046439560624409926616798792493328435787994562240567391560867539517753322628725069314972581967222209672892486307189350239216889094360499479590034032207920263636241298759966234243205570905766386207577079921475443490332191343090551009522870423225966591704205283036225548377573722514188921626990790748596113579717060293956259137133760787580672566454929838687772488012891423833655251357457695563011145177240090717870704091720599864440759816908963840000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
(~-0.05-1.67i, error: 1.05+1.57i, time: ~19s, length: 27219bit) λ exp ([(+3.14159265358) : (+0.0)] ⋅ i) (+4)
5067240574827142493824884283487389404038327255808591266950606056410685157365136888074630405753850936889648437500000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000/40894549326657251992544008824259638146259021596051752567291259765625000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
: -82857269451682879813602954637549515087487631913625936626050655364750241460569668561220169067382812500000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000/40894549326657251992544008824259638146259021596051752567291259765625000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000
(~0.12-2.02i, error: 0.88+1.57i, time: ~80s, length: 46851bit)λ # oof
Sadly, these results are really disappointing. Unfortunately the complex functions converge rather slowly and became a victim of many layered abstractions and recursions.
Also, if you’re interested: The unevaluated REPL term with R.π
(after complete substitution of definitions) has a BLC length of
1,106,776bit – that’s around 138kB! Compressing the term with my BLoC tool results in a
size of 704B. Many duplicated terms could obviously be shared by
abstraction directly, the large size doesn’t seem to impact performance
significantly, though.
Otherwise the techniques shown for real numbers apply here, too, so I won’t demonstrate them again. However, feel free to apply your gained knowledge by PRing new functions to bruijn’s standard library!
In my experience many readers won’t immediately understand, but working with bruijn and pure lambda calculus in general is really satisfying and enjoyable. I have yet to find another language that allows the same freedom and elegancy. I can only recommend digging into all the definitions and trying to understand everything – it’s definitely possible and rewarding.
Especially in the end, we’ve seen that this technique comes with large performance disadvantages. If someone wanted to use this for actual programs, they would need to invest a huge amount of work to make all of these functions even remotely efficient and accurate. This article should mostly serve as a fun proof of concept and notebook. I’ll probably still work on some obvious performance and accuracy improvements – you will find them on GitHub.
Aside from showing bruijn’s weaknesses, this project also showed me how fast bruijn can actually be. If you think about it, there are literally millions of reductions happening for the last formulas and the result still appears within a few seconds! Still, having to optimize even the smallest operations is something I don’t experience in normal programming languages. If I were to implement something similar in a language like Python or Java, decomposing powers or sqrt into the iteration could even hurt performance, since the traditional functions are already hyper-optimized.
Here, we can’t trust any optimizations and need to consider the
usefulness of every additional abstraction and, especially, layered
y
-recursion –
it seems like I overestimated the performance of fixed-point recursion a
lot. While some of this could be improved with specially optimized
reducers, I think this goes against the spirit of the project and
language. For me, it’s more interesting trying to implement procedures
that run efficiently on off the shelf reducers. I’m obviously
not there yet, though.
Thanks for reading. Contact me via email. Support on Ko-fi. Subscribe on RSS. Follow on Mastodon. Discuss on hacker news. Program in bruijn.
You can find the source of every described function (and a lot more) here or on GitHub.
Imprint · AI Statement · RSS