[haskell] optimalisatie fibonacci-getallen

[Knutselsmurf]

Naar aanleiding van een topic van een tijdje terug ben ik begonnen mijzelf te verdiepen in de wondere wereld van het functioneel programmeren. Als oefening daarin had ik voor mijzelf een programmaatje geschreven voor het uitrekenen van de fibonacci-getallen. Omdat de definitie daarvan recursief is, heb ik ook een recursieve definitie aangehouden:

code:

fib 0 = 1
fib 1 = 1
fib (n+2) = fib (n) + fib (n+1)

De code werkt, maar is gruwelijk traag voor grotere getallen, omdat er 100.000 keer hetzelfde uitgerekend moet worden. Nu is dit te optimaliseren, met behulp van bijvoorbeeld een lijst, maar in mijn ogen stap je dan af van het principe van FP, namelijk dat je beschrijft wat je wilt bereiken en niet direct hoe dat moet gebeuren. Wie heeft er ideeën over hoe dit te optimaliseren, zonder van zo'n mooie korte definitie af te stappen?

[Twilight Burn]

Bij Inleiding Functioneel Programmeren op de TU/e, hebben we dit bestandje (als voorbeeld) gekregen met allerlei manieren om Fibonacci getallen uit te rekenen:

code:

module Fib (fib, fibt, fibg, fibm, fibw) where



-- Fibonacci met inductietechniek

fib :: Integer -> Integer
fib 0 = 0
fib 1 = 1
fib (n+2) = fib n + fib (n+1)


-- Fibonacci met tuplingtechniek

fibt :: Integer -> Integer
fibt = fst . tfib 

tfib :: Integer -> (Integer, Integer)           
tfib 0 = (0,1)
tfib (n+1) = (b, a+b) where (a,b) = tfib n


-- Fibonacci met generalisatietechniek

fibg :: Integer -> Integer
fibg n = gfib 1 0 n

gfib :: Integer -> Integer -> Integer -> Integer
gfib a b 0 = b
gfib a b (n+1) = gfib b (a+b) n


-- Fibonacci via matrices

fibm :: Integer -> Integer
fibm n = mfib n 0 1 0 1

mfib :: Integer -> Integer -> Integer -> Integer -> Integer -> Integer
mfib n a b x y
  | n == 0         = x
  | odd(n)         = mfib (n-1) a b x' y'
  | otherwise      = mfib (n `div` 2) a' b' x y
   where
      a' = a^2 + b^2
      b' = 2 * a * b + b^2
      x' = a * x + b * y
      y' = b * x + (a + b) * y


-- Fibonacci via wortelformule

fibw :: Floating a => Integer -> a 
fibw n = c * ((power p n) - (power q n))
   where
      c = 1/wortel5
      p = (1 + wortel5)/2
      q = (1 - wortel5)/2
      wortel5 = sqrt 5


power :: Num a => a -> Integer -> a
power x n = gpower 1 x n

gpower :: Num a => a -> a -> Integer -> a
gpower m x k
  | k == 0       = m
  | odd(k)       = gpower (m*x) x (k-1)
  | otherwise    = gpower m (x*x) (k `div` 2)

[DRvDijk]

Originele code (en dus traag bij getallen boven de 20):

code:

fib :: Integer -> Integer
fib n
  | n == 0         = 1
  | n == 1         = 1
  | otherwise      = fib (n - 1) + fib (n - 2)

Waarom deze code zo traag is is omdat er veel te vaak hetzelfde getal uitgerekend wordt. Wat je eigenlijk wil is dat een bepaalde waarde van fib maar 1 keer uitgerekend wordt en daarna gewoon 'hergebruikt'.
Hoe los je dit op? Je kan een hulpfunctie maken die als resultaat een tupel van 2 integers oplevert:

code:

hulpFib :: Integer -> (Integer, Integer)
hulpFib 0 = (1, 1)
hulpFib n = (b, a + b)
  where (a,b) = hulpFib (n-1)

Deze functie is niet dubbelrecursief, maar enkelrecursief. Om het echt goeie antwoord op te leveren kun maak je dan nóg een functietje:

code:

snelFib :: Integer -> Integer
snelFib n = a
  where (a, b) = hulpFib n

[Verwijderd]

Laatst hoorde ik dat Haskell automatisch functies aanroepen cached. Dus hij rekent eenmalig fib 4 uit, en gebruikt bij verdere aanroepen steeds het antwoord wat ie al heeft uitgerekend. Aangezien er geen side-effects zijn, zou dit ook wel eens kunnen.
Nu weet ik het niet zeker of dit echt gebeurd, maar als dat wel het geval is, veroorzaakt iets anders de traagheid.

[MaxxRide]

Inderdaad, een functionele programmeertaal rekent maar een keer een bepaalde functiewaarde uit. In principe zijn de 'variabelen' in een programmeertaal niet aanpasbaar, hebben ze een keer een waarde toegekend gekregen dan veranderd deze nooit weer. Clean en Haskell 'cachen' inderdaad de resultaten.

Tevens kennen ze lazy evaluation, wat inhoudt dat ze iets pas berekenen als ze het nodig hebben.

Als ze een bepaald resultaat niet meer nodig hebben wordt dit trouwens wel 'weggegooid'.

Althans zo heb ik het allemaal geleerd

[DRvDijk]

Haskell kent inderdaad lazy evaluation. Zo kent het dus ook een oneindige lijst van alle integers [1 ..], die je gewóón kan gebruiken zónder dat er een oneindige loop ontstaat.

Dat fib 4 echter maar 1 keer uitgerekent wordt is niet waar. Waar zou de traagheid ánders vandaan moeten komen? Het is een redelijk straight-forward fuctietje niet? fib 0 = 1, fib 1 = 0, fib x = fib (x-1) + fib (x-2). De traagheid móét in de laatste regel zitten, en die rekent dus meermalig hetzelfde uit.

[Glimi]

OMG, nooit gedacht dat ik dit nog nodig zou hebben

De algemene vorm van recursie is vaak (en ook zo bij fibonacci) a_n = Aa_n-1 + Ba_n-2 ...
Dit is om te rekenen naar een polynoom met graad x waarbij x gelijk is aan de diepte van de recursieve betrekking, in dit geval 2. Oftwel w² -Aw -B = 0

Nadat je hier de nulpunten van berekend (kan ook complex zijn) vul je deze in in de algemene formule, waarbij w1 en w2 de gevonden nulpunten zijn:
a_n = Cw₁² + Dw₁²
Dit kun je invullen voor a₀ en a₁ (de startwaardes). Je krijgt dan een specefieke oplossing voor de recursie.

Bij fibonacci met startwaardes 1 en 1 wordt dat:
a_n = (1/ ( sqrt(5) * 2 ⁿ⁺¹ ) ) * ( ( 1+ sqrt(5) )ⁿ⁺¹ - ( 1- sqrt(5) )ⁿ⁺¹ )

Dus als je echt alleen fibonacci wilt berekenen met startwaardes 1 en 1, ben je met deze berekening echt een stuk sneller klaar

[Zoijar]

Lineaire twee orde homogene recurrente betrekking met constante coefficienten

Glimli's uitleg is perfect, maar nog ter opheldering: (ik zeg eingelijk hetzelfde maar ik denk iets duidelijker)

We proberen een oplossing van de vorm k*r^n, dan krijgen we:
a_n-a_n-1-a_n-2 = k*r^n - k*r^(n-1) - k*r^(n-2)

= k*r^(n-2)*(r^2-r-1)

= 0

Dus geldt r^2-r-1=0, nulpunten zijn r₁=(1+sqrt(5))/2 en r₂=(1-sqrt(5))/2

We hebben nu twee onafhankelijke algemene oplossingen die een lineaire basis vormen voor de specifieke oplossingen; dus:

a_n=k₁*r₁ + k₂*r₂

Voor a₀ wordt dit:

a₀ = k₁+k₂
Dus: k₂=1-k₁

Nu voor a₁:
a₁ = k₁*sqrt(5)+(1/2)-(1/2)*sqrt(5) = 1

Dus:

k₁ = (1/2)+sqrt(5)/10
k₂ = (1/2)-sqrt(5)/10

En daaruit volgt de uiteindelijke oplossing:

a_n = ((1/2)+sqrt(5)/10)*((1+sqrt(5))/2)^n + ((1/2)-sqrt(5)/10)((1-sqrt(5))/2)^n

[Glimi]

En of u het iets duidelijker zegt Ik had alleen ff geen tijd om een hele berekening (die overigens niet zo moeilijk is, het probleem begint pas echt als je een 3 of hogere graads polynoom hebt ) op te schrijven. Ik zie dat jij ook geen zin meer had om sup en sub te gebruiken op het eind Maar ik wil u toch bedanken voor de verduidelijking

[HCMarX]

[ontopic]
Zoijar
ISO/IEC 14882:1998
gaf inderdaad het goede antwoord, wij berekenen nu ook polynomen op school uit rijen. Is wel grappig dat je het weet, het leren is wat minder.....

Gang-ster

[ Voor 0% gewijzigd door D2k op 13-10-2002 20:07 . Reden: D2k: opschonen ]

[D2k]

heren: ontopic blijven
vragen van n00bs over Java horen niet in een haskell topics

[Okyah]

mkae

[D2k]

Okyah schreef op 13 oktober 2002 @ 20:13:
[...]



ik snap de onvriendelijkheid niett
ik vroeg juist wat ik eventueel zou kunnen gebruiken.
waarom je zo kwaad wordt vat ik echt niet...
tis niet alsof ik t hele topic verziek

en hiervoor is dus icq/mail
niet binnen een topic
en ja het werd vervelend gevonden
en ik was het er mee eens

---

nu weer ontopic

[CyberSnooP]

Knutselsmurf schreef op 12 oktober 2002 @ 22:53:
Wie heeft er ideeën over hoe dit te optimaliseren, zonder van zo'n mooie korte definitie af te stappen?

Zoals je hierboven duidelijk geworden zal zijn is het heel goed mogelijk een (/jouw) oplossing te optimaliseren met een boel extra regels code. Het is helaas niet (direct) mogelijk een dergelijke optimalisatie te bereiken met je korte definitie zoals je misschien zou willen.

Dat is ook niet helemaal het idee achter functioneel programmeren. Je begint met je korte definitie zoals jij die geeft en vormt zo een overduidelijk correct programma. Door nu wat (correcte) transformaties to te passen zoals genoemd in de 2e post (tuplingtechniek etc..) krijg je een sneller programma wat (omdat de transformatie eerder bewezen is) nog steeds hetzelfde bereikt. Die transformatie zorgt er echter wel voor dat je definitie groter wordt.. jammer maar helaas

Ik meen dat sommige transformatie echter (relatief) eenvoudig te omschrijven zijn waardoor ze ook door software uitgevoerd kunnen worden. Er zal dus vast ergens iets een soort van Haskell-parser te verkrijgen zijn die op je functie poogt van die optimalisatie technieken toe te passen zonder dat jij dat weet. Op die manier bereik je weer je ideale situatie waarin je met een eenvoudige definitie ook nog een snel programma bouwt.

[Verwijderd]

Ik denk dat het uitschrijven c.q. oplossen van de polynomen wel een erg extreme optimalisatie is voor iemand die net met functioneel programmeren begonnen is! De oplossing met tuplen lijkt me meer in de "geest" van het leren functioneel programmeren passen.

CyberSnoop's hint lijkt me ook erg nuttig om eens te bekijken.

[Glimi]

Verwijderd schreef op 13 oktober 2002 @ 21:04:
Ik denk dat het uitschrijven c.q. oplossen van de polynomen wel een erg extreme optimalisatie is voor iemand die net met functioneel programmeren begonnen is! De oplossing met tuplen lijkt me meer in de "geest" van het leren functioneel programmeren passen.

CyberSnoop's hint lijkt me ook erg nuttig om eens te bekijken.

Het zoeken van een specefieke oplossing met Polynomen is ook niet altijd mogelijk. Kijk bij polynomen van graad 2 is het nog makkelijk te doen. Voor graad 3 is er ook een formule en voor graad 4 ook. Wat echter als er een graad 5, 6, 7 op komt draven? Dan moet je aan de slag met gokwerk en hopen dat iets een oplossing is. Dat is niet echt makkelijk te doen dmv een functioneel programma

Echter de vraag was het 'versnellen' van fibonacci. De reeks van fibonacci is een vaststaande reeks met vaste begingetallen! Het het omschrijven naar een polynoom en oplossen is helemaal geen moeilijke constructie die ook maar iets met functioneel programmeren te maken heeft, maar eerder met wiskunde! Daarom is het voor een optimalisatie binnen een functionele programmeertaal ook niet 'afstappen van het functionele princiepe'.

[Knutselsmurf]

Om eerlijk te zijn was ik niet eens op de hoogte van het feit dat je een reeks van gehele getallen om kunt schrijven naar een polynoom. Wat mij daaraan nog verwonderd is het feit dat er niet-rationele coefficienten gebruikt worden en dat het blijkbaar toch netjes uitkomt. De oplossing met tupels komt idd nog het dichtst bij de oorspronkelijke definitie, waarbij geen dubbele recursie meer voorkomt. Het jammere aan andere oplossingen is, dat het verband met de oorspronkelijke, lekker kort omschreven definitie niet meer goed zichtbaar is.

[CyberSnooP]

Knutselsmurf schreef op 13 oktober 2002 @ 21:35:
Het jammere aan andere oplossingen is, dat het verband met de oorspronkelijke, lekker kort omschreven definitie niet meer goed zichtbaar is.

Nee, maar dat moet je dus volgens mij niet echt meer willen. Zodra je aan het transformeren slaat doe je dat om snelheid te winnen.. je code is in eerste instantie altijd het begrijpelijkst, dat neemt dus alleen maar af.

Het verband is er echter wel, namelijk de gebruikte transformatie. Wil je dus meer te weten komen moet je die transformaties leren begrijpen.

BTW.. Het omschrijven tot een polynoom valt eigenlijk niet meer binnen het functioneel programmeren, maar meer binnen de wiskunde toch?

[Zoijar]

offtopic:
Je kan elke reeks schrijven als polynoom. En dat van die niet-rationele coefficienten en toch netjes uitkomen is een algebraisch iets...maar gaat allemaal wat verder en dieper dan hier van toepassing is.

[Glimi]

Dat van de niet- reële wortels die je uit de nulpunten van de polynoom krijgt wordt geëlimineerd in de C en de D van mijn bovenstaand voorbeeld. i² is namelijk gewoon -1 + i0

[Soultaker]

Zelf vind ik het in dit geval wel aardig, om gewoon de lijst van alle Fibonnacci getallen te generen en daar de gewenste elementen uit te pakken:

Clean:

fibs a b = [ a : fibs b (a+b) ]

'fibs 1 1' is dan bijvoorbeeld de Fibonacci rij die begint met twee enen (zoals 'ie meestal bekend staat) maar niets weerhoudt je er van om 'm met andere getallen te beginnen. Uit deze definitie is ook duidelijk dat elke functie-'aanroep' een element aan de lijst toevoegt. De complexiteit van dit algoritme is dan ook O(N).

Mijn functie voor het berekenen voor het zoveelste Fibonacci getal, zou dus zoiets worden:

Clean:

fib n = (fibs 1 1) !! n

Effectief komt dit overeen met Twilight Burn's "Fibonacci met generalisatietechniek" maar ik vind dit nog wat mooier omdat het nog declaratiever is. Dezelfde methode is geschikt om allerlei functies op andere rijen te definiëren (zoals priemgetallen e.d.).

_{Omdat Parse 't er schijnbaar niet meer bij wil zetten: mijn voorbeeldcode is in Clean; '!!' is volgens mij specifiek van Clean en geeft het zoveelste (hier n-de) element uit een lijst.}

[Knutselsmurf]

Soultaker schreef op 13 oktober 2002 @ 23:14:
Zelf vind ik het in dit geval wel aardig, om gewoon de lijst van alle Fibonnacci getallen te generen en daar de gewenste elementen uit te pakken:

Clean:

fibs a b = [ a : fibs b (a+b) ]

'fibs 1 1' is dan bijvoorbeeld de Fibonacci rij die begint met twee enen (zoals 'ie meestal bekend staat) maar niets weerhoudt je er van om 'm met andere getallen te beginnen. Uit deze definitie is ook duidelijk dat elke functie-'aanroep' een element aan de lijst toevoegt. De complexiteit van dit algoritme is dan ook O(N).

Mijn functie voor het berekenen voor het zoveelste Fibonacci getal, zou dus zoiets worden:

Clean:

fib n = (fibs 1 1) !! n

Effectief komt dit overeen met Twilight Burn's "Fibonacci met generalisatietechniek" maar ik vind dit nog wat mooier omdat het nog declaratiever is. Dezelfde methode is geschikt om allerlei functies op andere rijen te definiëren (zoals priemgetallen e.d.).

_{Omdat Parse 't er schijnbaar niet meer bij wil zetten: mijn voorbeeldcode is in Clean; '!!' is volgens mij specifiek van Clean en geeft het zoveelste (hier n-de) element uit een lijst.}

Deze methode met een lijst ken ik en heb ik zelf ook toegepast in bijvoorbeeld Pascal. Ik had echter geen idee hoe deze methode in Haskell te gebruiken.....

[Soultaker]

Pinda schreef op 13 oktober 2002 @ 11:48:Clean en Haskell 'cachen' inderdaad de resultaten.

Over Haskell kan ik weinig zeggen (vooral omdat er zoveel Haskell-achtige systemen zijn) maar Clean cached in ieder geval niets. Miranda doet dit voor zover ik weet wel; Amanda misschien ook wel.

Het is een leuke feature, maar je kunt beter goede code schrijven, denk ik. Je zal in het algemeen immers niet precies weten hoe groot je cache is en hoe 'duur' je aanroepen zijn.

[MaxxRide]

Oh Soultaker, bij ons werd op school vertelt dat Clean dat wel deed.


Maar idd je kunt beter goed programmeren