7. Još programiranja

7.1. Mogućnosti ispisa rezultata:

Sage notebook sučelje može ispisati rezultat računa u različitim formatima:

  1. Defaultni ispis jedini omogućuje izravni copy/paste u input ćelije.
sage: var('x n theta')
(x, n, theta)
sage: s1 =sum(sin(theta+x)/log(n+x), x, 1, oo)
sage: s1
sum(sin(theta + x)/log(n + x), x, 1, +Infinity)
sage: s2 = matrix(ZZ, [[1, 2, 3], [4, 5, 6]])
sage: s2
[1 2 3]
[4 5 6]
  1. Ljepši ispis koji koristi interni LaTeX i tzv. jsMath dobiva se funkcijom show(), ili aktiviranjem “Typeset” gumba na vrhu notebook sučelja.
sage: show(s1)
sage: show(s2)
\[{\rm sum}\left(\frac{\sin\left(\theta + x\right)}{\log\left(n + x\right)}, x, 1, +\infty\right)\]
\[\begin{split}\left(\begin{array}{rrr} 1 & 2 & 3 \\ 4 & 5 & 6 \end{array}\right)\end{split}\]
  1. Ispis odgovarajućeg koda za copy/paste u LaTeX dokument.
sage: latex(s2)
\left(\begin{array}{rrr}
1 & 2 & 3 \\
4 & 5 & 6
\end{array}\right)
  1. Prikaz LaTeX->PDF inaćice u posebnom prozoru (Moguće je da ovo radi samo ukoliko su Sage klijent i server na istoj mašini ili uz pažljivo podešene X-windows autorizacije.)
sage: view(s2, viewer='pdf')

7.2. Još o funkcijama

Funkciju definiranu po dijelovima možemo dobiti na slijedeći način:

sage: def fpw(x):
....:     if x<-1:
....:         return -x
....:     elif x>1:
....:         return x
....:     else:
....:         return x^2
sage: plot(fpw, -2, 2)
../_images/cell_74_sage0.png

Python funkcije ne možemo općenito npr. simbolički integrirati. Posebno je opasno to što, ukoliko pokušamo, možemo dobiti pogrešan odgovor:

sage: integral(fpw(x), x, 1, 2)
7/3

No, uvijek možemo pribjeći numeričkoj integraciji koja u ovom slučaju daje točan rezultat 3/2:

sage: numerical_integral(fpw, 1, 2)[0]
1.5

Pri gornjem pokušaju simboličke integracije, prvo je izvrijednjen sam integrand fpw(x). Kako u tom trenutku x nema nikakvu vrijednost, on ne zadovoljava ni jedan od dva uvjeta (x>1 ili x<1) pa funkcija vraća simbolički izraz x^2 (is else bloka), što je pogrešno za integraciju u granicama od 1 do 2.

sage: fpw(x)
x^2

Iz istog razloga je prilikom crtanja gore bilo potrebno kao prvi argument staviti samu funkciju fpw, a ne izraz fpw(x) koji bi dao krivi crtež:

sage: plot(fpw(x), x, -2, 2)
../_images/cell_81_sage0.png

Naravno, ukoliko znamo da python funkcija uvijek vraća korektan simbolički izraz, smijemo je integrirati i diferencirati:

sage: def gun(f, x):
....:     return f(f(x))
sage: diff(gun(sin, x), x)
cos(x)*cos(sin(x))

7.3. Funkcionalno programiranje

Isti matematički problem je često moguće riješiti na različite načine, putem algoritama iza kojih stoje sasvim različiti načini razmišljanja. Razložimo to na primjeru funkcije faktorijel.

sage: factorial(7)
5040

Klasični način programiranja, upotrebljavan od dana prvih računala, je tzv. proceduralno (ili imperativno) programiranje kod kojeg na program gledamo kao na niz naredbi koje obično postupno mijenjaju vrijednosti nekih varijabli pohranjenih u memoriji računala:

sage: def factorialProc(n):
....:     fac = 1
....:     i = 1
....:     while i<n:
....:         i = i + 1
....:         fac = fac * i
....:     return fac
sage: factorialProc(7)
5040

Ovdje je očito riječ o standardnom algoritmu kakvog bismo isprogramirali u bilo kojem proceduralnom jeziku poput Fortrana ili C-a: Imamo petlju koja se prolazi n puta i svaki puta množimo rezultat prošlog prolaza sa za 1 većim brojem.

Međutim, postoje i drugačiji pristupi. Tako je kod funkcionalnog programiranja naglasak na izvrijednjavanju funkcija tj. primjeni funkcija na izraze. Faktorijel prirodno zamišljamo kao funkciju koja daje “umnožak svih brojeva od 1 do zadanog broja”.

Kao prvo, treba nam modul operator koji implementira funkcije koje odgovaraju standardnim matematičkim operacijama *, +, -,...

sage: import operator
sage: (operator.mul(2,3), operator.add(2,3))
(6, 5)

Zatim koristimo funkciju reduce(), koja kumulativno primjenjuje prvi argument (koji mora biti funkcija) na elemente drugog argumenta:

sage: f = function('f')
sage: reduce(f, range(1,5))
f(f(f(1, 2), 3), 4)
sage: reduce(operator.mul, range(1,8))
5040

Primijetite da nam ovdje nisu trebale pomoćne varijable (poput i i fac iz factorialProc), ali nam je trebalo više memorije jer smo trebali pohraniti čitavu listu [1, 2, ...] koju daje range().

Funkcionalno programiranje je stil programiranja koji stavlja naglasak na izvrijednjavanje izraza, a ne na sukcesivno izvršavanje komandi. Kod čistih funkcionalnih programa obično nema pomoćnih varijabli i nepotrebnih popratnih efekata izvrijednjavanja funkcija. U tom smislu funkcionalno programiranje je pomalo slično radu s tabličnim kalkulatorom (npr. MS Excel): redosljed izračunavanja ćelija je nebitan tj. očekujemo automatsku konzistenciju. Isto, vrijednosti ćelija su dane izrazima, a ne slijedovima komandi. (Misli se na Excel, a ne Sage ćelije.)

(Više informacija o funkcionalnom programiranju nudi Functional Programming FAQ ili WWW stranice Haskell funkcionalnog jezika.)

Takav stil programiranja često rabi nekoliko specijalnih funkcija (mogli bismo ih zvati “meta-funkcije”) čija je uloga upravljanje primjenivanjem drugih funkcija koje dolaze kao argumenti ovih meta-funkcija. Možda najvažnija takva funkcija je map() koja distribuira (mapira) funkciju po elementima neke liste:

sage: f = function('f')
sage: map(f, range(4))
[f(0), f(1), f(2), f(3)]

Ukoliko funkcija prima više argumenata, može se mapirati na više listi:

sage: map(f, range(4), range(3,7))
[f(0, 3), f(1, 4), f(2, 5), f(3, 6)]

Recimo da sad želimo ispisati tablicu s nizom prirodnih brojeva i njihovih faktorijela. Možemo prvo postupiti tako da prvo definiramo pomoćnu funkciju koja generira jedan red tablice i onda je pomoću map() primijenimo na niz brojeva:

sage: def auxfun(k):
....:     return (k, factorial(k))
sage: map(auxfun, range(7))
[(0, 1), (1, 1), (2, 2), (3, 6), (4, 24), (5, 120), (6, 720)]

Međutim, baš kao i pomoćne varijable, tako ni pomoćne funkcije nisu u duhu funkcionalnog programiranja. Zato postoje tzv. lambda-funkcije koje su bezimene i definiraju se i upotrebljavaju na slijedeći način:

sage: map(lambda k: (k, factorial(k)), range(7))
[(0, 1), (1, 1), (2, 2), (3, 6), (4, 24), (5, 120), (6, 720)]
sage: matrix(_)  # čitljiviji ispis bez uptrebe formatirajućih stringova
[  0   1]
[  1   1]
[  2   2]
[  3   6]
[  4  24]
[  5 120]
[  6 720]

Treba uočiti kako je često upotrebu funkcije map() i lambda-funkcija moguće izbjeći korištenjem obuhvaćanja liste. Npr, gornji primjer se elegantnije realizira ovako:

sage: [(k, factorial(k)) for k in range(7)]
[(0, 1), (1, 1), (2, 2), (3, 6), (4, 24), (5, 120), (6, 720)]

Kako su Sage/Python funkcije punopravni objekti, možemo i sami definirati funkcije koje primaju druge funkcije kao argumente. Evo funkcije koja implementira kompoziciju funkcija:

sage: def compose(f, n, x):
....:     """Uzastopno komponiranje funkcije sa samom sobom n puta."""
....:     if n <= 0:
....:         return x
....:     x = f(x)
....:     for i in range(n-1):
....:         x = f(x)
....:     return x

Npr. tzv. zlatni omjer \((\sqrt{5}+1)/2\) se može izraziti kao beskonačni razlomak

\[\frac{\sqrt{5}+1}{2} = 1+\frac{1}{1+\frac{1}{1+\frac{1}{1+\frac{1}{\cdots}}}} = 1.618034\ldots\]

Takav razlomak možemo izvrijedniti na slijedeći način:

sage: gr = compose(lambda x: 1+1/x, 5, x); gr
1/(1/(1/(1/(1/x + 1) + 1) + 1) + 1) + 1
sage: gr.subs(x=1).n()
1.62500000000000

Ili, uz povećanje točnosti:

sage: numerical_approx((1+sqrt(5))/2)
1.61803398874989
sage: compose(lambda x: 1+1/x, 32, x).subs(x=1).n()
1.61803398874986

7.3.1. Primjer razvoja funkcionalnog programa

Promotrimo razvoj programa koji ispisuje tablicu frekvencija pojavljivanja elemenata u listi. Načinimo prvo jednostavnu testnu listu

sage: lst = ['a', 'c','b', 'a', 'c', 'a', 'd', 'b', 'c','d', 'c']

Metoda liste koja daje broj pojavljivanja nekog elementa je count()

sage: lst.count('c')
4

Da bismo ispisivali tablicu trebamo listu oblika [(element1, frekvencija elementa1), (element2, frekvencija elementa2),...]. Element te liste (red tablice) se može dobiti ovako:

sage: def el(x):
....:     return (x, lst.count(x))
sage: el('c')
('c', 4)

Tablicu frekvencija za različite elemente dobijemo korištenjem lambda-funkcije analogne el() i njenom primjenom pomoću funkcije map() na popis elemenata:

sage: map(el, ['a', 'b', 'c', 'd'])
[('a', 3), ('b', 2), ('c', 4), ('d', 2)]
sage: map(lambda x: (x, lst.count(x)), ['a', 'b', 'c', 'd'])
[('a', 3), ('b', 2), ('c', 4), ('d', 2)]

To je sad to, jedino još želimo izbjeći da sami moramo ručno identificirati koji se sve elementi pojavljuju u listi. To nam može raditi funkcija uniq() koja izbacuje duplikate iz liste:

sage: res = map(lambda x: (x, lst.count(x)), uniq(lst)); res
[('a', 3), ('b', 2), ('c', 4), ('d', 2)]

Kao zadnju stvar, poželjno je listu sortirati po frekvencijama. Funkcija sorted() je već definirana tako da zna sortirati liste različitih objekata, no ovdje bi po defaultu sortirala parove po prvom elementu i to po abecednom redu, ...

sage: sorted(res)
[('a', 3), ('b', 2), ('c', 4), ('d', 2)]

... dok nama treba sortiranje po drugom elementu i to po veličini. Stoga moramo putem opcionalnog argumenta key funkciji sorted() proslijediti funkciju koja će vraćati onaj dio elementa liste po kojem želimo sortiranje:

sage: def keyfun(item):
....:     """Return last element of item."""
....:     return item[-1]
sage: sorted(res, key=keyfun, reverse=True)
[('c', 4), ('a', 3), ('b', 2), ('d', 2)]

(Opcionalni argument reverse daje silazno sortiranje umjesto defaultnog uzlaznog.)

Naravno, i ovu pomoćnu funkciju keyfun() možemo zamijeniti lambda funkcijom:

sage: sorted(res, key=lambda it: it[-1], reverse=True)
[('c', 4), ('a', 3), ('b', 2), ('d', 2)]

I to je to. Sad definiramo kompletnu funkciju:

sage: def frequencies(lst):
....:     """Elements of list lst with their frequencies."""
....:     return sorted(map(lambda x: (x, lst.count(x)), set(lst)),
....:            key=lambda it: it[-1], reverse=True)
sage: for it in frequencies(lst):
....:     print "%s: %d" % it
c: 4
a: 3
b: 2
d: 2

Primijenimo to na frekvenciju pojavljivanja slova u Shakespeareovom Mletačkom trgovcu (zapravo koristimo engleski original The Merchant of Venice, sa WWW stranica projekta Gutenberg)

sage: import urllib
sage: venice = urllib.urlopen(
....:   'http://www.gutenberg.org/ebooks/2243.txt.utf-8').read()[16400:]
sage: len(venice)
120784
sage: print venice[:370]
The Merchant of Venice

Actus primus.

Enter Anthonio, Salarino, and Salanio.

  Anthonio. In sooth I know not why I am so sad,
It wearies me: you say it wearies you;
But how I caught it, found it, or came by it,
What stuffe 'tis made of, whereof it is borne,
I am to learne: and such a Want-wit sadnesse makes of
mee,
That I haue much ado to know my selfe
sage: for it in frequencies(venice)[:10]:
....:     print "%s:  %d" % it
 :  20927
e:  12250
o:  7260
t:  7227
a:  6233
h:  5560
n:  5518
s:  5317
r:  5232
i:  5206

Vidimo da je “e” daleko najčešće slovo (poslije razmaka), činjenica vrlo važna pri dešifriranju engleskih tekstova.

Kao netrivijalniji primjer upotrebe gore definirane funkcije compose() možemo nacrtati tzv. Mandelbrotov skup, definiran kao skup svih točaka \(c\) kompleksne ravnine za koje iteracija

\[z_0 = 0\;; \qquad z_{n+1} = z_{n}^2 + c\]

ne divergira.

sage: complex_plot(compose(lambda z: z^2+x, 8, x), (-2, 1), (-1, 1))
../_images/cell_33_sage0.png

(Analizirajte sami kako radi ovaj “program”.)

Zadatak 1

Definirajte funkciju allequal() koja testira da li su svi elementi neke liste jednaki i iskoristite je da pronađete neprekidni niz od 6 istih znamenki u prvih 1000 decimala broja \(\pi\). (Za pretvaranje broja u listu znamenaka možete iskoristiti funkciju str().) Koristite ili obuhvaćanje liste ili map() tj. nije dopušteno korištenje petlji osim eventualno unutar funkcije allequal(), gdje to isto nije nužno (Naputak: all()).

Zadatak 2

Logističko preslikavanje je dano iteracijskom formulom \(x_{n+1}= \lambda x_n(1-x_n)\). Za male vrijednosti parametra \(\lambda\), to preslikavanje za veliki \(n\) konvergira k jednoj vrijednosti. Kad \(\lambda\) raste, negdje blizu \(\lambda=3\), dolazi do bifurkacije i iteracije više ne konvergiraju već skaču između dvije vrijednosti. S daljnjim rastom \(\lambda\) dolazi do nove bifurkacije i sustav se za veliki \(n\) ponavlja s periodom četiri, itd. Rezultat se može prikazati kao tzv. Feigenbaumov bifurkacijski dijagram (vidi sliku). Nacrtajte ga tako da prvo definirate funkciju composelist() koja je analogna gornjoj funkciji compose(), ali ne vraća samo kranji rezultat kompozicije funkcija već i listu sa svim međukoracima. Nakon toga iskoristite tu funkciju za iteriranje logističkog preslikavanja. Eliminirajte prvi dio liste (dok se iteracije ne stabiliziraju) i nacrtajte drugi dio pomoću list_plot(...., pointsize=1). Uočite da svaka vrijednost apscise lambda traži posebno iteriranje.

Zadatak 3

Odredite najmanji prirodni broj koji se ne može dobiti iz brojeva 2, 3, 4 i 5 korištenjem operacija zbrajanja, oduzimanja i množenja, a gdje se svaki broj smije koristiti najviše jednom. (Npr. 1=3-2, 2=3-5+4, ..., 6=2*3, ..., 14=4*5-3*2, ..., 30=(2+4)*5, .... Naputak: Od koristi se možda može pokazati već ranije korišteni modul operator, te funkcije permutations() i combinations.

Zadatak 4

Prim-brojevi blizanci su parovi prim-brojeva koji se razlikuju za dva, poput (3,5) ili (17,19). V. Brun je 1919. dokazao da suma recipročnih vrijednosti prim-brojeva blizanaca konvergira

\[B=\left(\frac{1}{3}+\frac{1}{5}\right) + \left(\frac{1}{5}+\frac{1}{7}\right) + \left(\frac{1}{11}+\frac{1}{13}\right) + \left(\frac{1}{17}+\frac{1}{19}\right) + \cdots = 1.902160583104\]

Ovako preciznu vrijednost je vrlo teško dobiti, no napišite funkciju brun(n) koja izračunava \(B\) koristeći listu od prvih n prim-brojeva. Pokušajte je isprogramirati i unutar paradigme funkcionalnog i unutar proceduralnog programiranja. Inače, zanimljivo je da je upravo računalno određivanje ove konstante ukazalo na bug u prvoj generaciji Pentium procesora.

Literatura