Matematički softver

Sadržaj kolegija

  • Python
  • Numpy / Scipy
  • Pandas
  • Sympy
  • Sage
  • Markdown
  • $\LaTeX$

Kolegij će se bazirati na korištenju web platforme Sagemath cloud. Materijali za kolegij se nalaze na web stranici kolegija, na Sagemath oblaku te u Git repozitoriju kolegija.

Plan kolegija

  1. Uvod u Python; upoznavanje s web platformom
  2. Markdown, $\LaTeX$ notacija za matematičke formule
  3. Python kao alat za znanstvenike (Matplotlib, Pandas, ...)
  4. Sage
  5. $\LaTeX$

Polaganje kolegija

Ocjena se formira na osnovu dva kolokvija te domaćih zadaća.

Kolokviji open-book tipa.

Kolokviji nose najviše 80 bodova.

Domaće zadaće u obliku eseja, sa zadanim elementima.

Predaja ispravno napravljenih svih domaćih zadaća nužan je uvjet za polaganje kolegija. Domaća zadaća je uspješno predana samo ako ima sve tražene elemente.

Domaće zadaće nose do 20 bodova, ali na sljedeći način: studenti čije domaće zadaće odskaču kvalitetom dobit će do 20 dodatnih bodova.

Popravni ispit nije predviđen.

Literatura

GIYF

Kolegij se dosta oslanja na Lectures on scientific computing with Python.

Python

Prvi primjeri, korištenje Jupyter (IPython) notebook-a

Jupyter notebook ima puno mogućnosti, koje ćemo upoznavati tijekom kolegija. Jedna od njih je i jednostavan pristup datotečnom sustavu.


In [4]:
ls


astro.png         mojmodul.py   terminal.term  Uvod.ipynb
custom.css        mojmodul.pyc  thumb2.png     Uvod.slides.html
integral_demo.py  __pycache__/  thumb.png      verzije.py

In [5]:
pwd


Out[5]:
'/projects/95cb524c-f1ac-4a19-908e-5f47c10f47c1/01. Uvod; Python'

In [6]:
%matplotlib inline

In [ ]:
%load integral_demo.py

In [8]:
integral_plot(func)



In [9]:
from IPython.display import Image
Image(url='http://python.org/images/python-logo.gif')


Out[9]:

In [10]:
from IPython.display import HTML
HTML('<iframe width="560" height="315" src="//www.youtube.com/embed/-SqzarfTqRY" frameborder="0" allowfullscreen></iframe>')


Out[10]:

In [11]:
from IPython.display import YouTubeVideo
from IPython.display import display
display(YouTubeVideo("GMKZD1Ohlzk",width="560",height="315", start=3))


Jedan primjer korištenja Jupytera. Verzija sa svim pomoćnim datotekama na Sagemath oblaku.

Primjer korištenja Pythona

Krenimo s nečim korisnim. Učitajmo sliku, prikažimo je, prikažimo dio slike te smanjimo šum na tom dijelu slike.


In [12]:
#Scikit-image je paket za manipulaciju slika
from skimage import data
import matplotlib.pyplot as plt
coins = data.coins()
plt.imshow(coins,cmap='gray');


Što se dešava u gornjem kodu?

Prvi redak koda služi da se iz biblioteke scikit-image učita pod-bibiloteka data.

U drugom retku iz biblioteke matplotlib učitavamo pod-biblioteku pyplot. Kako ćemo pyplot često koristiti dajemo joj kraće ime plt.

U trećem retku iz biblioteke data učitavamo objekt coins. To je u ovom slučaju slika novčića.

Na koncu prikazujemo sliku.


In [13]:
print(coins.shape)


(303, 384)

Varijabla coins.shape sadrži par brojeva.


In [14]:
a,b=coins.shape
print(a)
print(b)


303
384

In [15]:
from skimage import restoration
coins_zoom = coins[10:80, 300:370]
tv_coins = restoration.denoise_tv_chambolle(coins_zoom, weight=0.1)
plt.figure()
plt.subplot(1,2,1)
plt.imshow(coins_zoom,cmap='gray'); 
plt.subplot(1,2,2)
plt.imshow(tv_coins, cmap='gray');


Što točno radi funkcija denoise_tv_chambolle? To možemo saznati na sljedeći način.


In [16]:
help(restoration.denoise_tv_chambolle)


Help on function denoise_tv_chambolle in module skimage.restoration._denoise:

denoise_tv_chambolle(im, weight=0.1, eps=0.0002, n_iter_max=200, multichannel=False)
    Perform total-variation denoising on n-dimensional images.
    
    Parameters
    ----------
    im : ndarray of ints, uints or floats
        Input data to be denoised. `im` can be of any numeric type,
        but it is cast into an ndarray of floats for the computation
        of the denoised image.
    weight : float, optional
        Denoising weight. The greater `weight`, the more denoising (at
        the expense of fidelity to `input`).
    eps : float, optional
        Relative difference of the value of the cost function that
        determines the stop criterion. The algorithm stops when:
    
            (E_(n-1) - E_n) < eps * E_0
    
    n_iter_max : int, optional
        Maximal number of iterations used for the optimization.
    multichannel : bool, optional
        Apply total-variation denoising separately for each channel. This
        option should be true for color images, otherwise the denoising is
        also applied in the channels dimension.
    
    Returns
    -------
    out : ndarray
        Denoised image.
    
    Notes
    -----
    Make sure to set the multichannel parameter appropriately for color images.
    
    The principle of total variation denoising is explained in
    http://en.wikipedia.org/wiki/Total_variation_denoising
    
    The principle of total variation denoising is to minimize the
    total variation of the image, which can be roughly described as
    the integral of the norm of the image gradient. Total variation
    denoising tends to produce "cartoon-like" images, that is,
    piecewise-constant images.
    
    This code is an implementation of the algorithm of Rudin, Fatemi and Osher
    that was proposed by Chambolle in [1]_.
    
    References
    ----------
    .. [1] A. Chambolle, An algorithm for total variation minimization and
           applications, Journal of Mathematical Imaging and Vision,
           Springer, 2004, 20, 89-97.
    
    Examples
    --------
    2D example on astronaut image:
    
    >>> from skimage import color, data
    >>> img = color.rgb2gray(data.astronaut())[:50, :50]
    >>> img += 0.5 * img.std() * np.random.randn(*img.shape)
    >>> denoised_img = denoise_tv_chambolle(img, weight=60)
    
    3D example on synthetic data:
    
    >>> x, y, z = np.ogrid[0:20, 0:20, 0:20]
    >>> mask = (x - 22)**2 + (y - 20)**2 + (z - 17)**2 < 8**2
    >>> mask = mask.astype(np.float)
    >>> mask += 0.2*np.random.randn(*mask.shape)
    >>> res = denoise_tv_chambolle(mask, weight=100)

Analogno smo mogli koristiti i restoration.denoise_tv_chambolle?

Još jedan primjer korištenja biblioteka.


In [17]:
import math
print ('pi je (otprilike): {}'.format(math.pi))
print ('pi je (otprilike): {:f}'.format(math.pi))
print ('pi je (otprilike): {:1.8f}'.format(math.pi))


pi je (otprilike): 3.141592653589793
pi je (otprilike): 3.141593
pi je (otprilike): 3.14159265

In [18]:
import math as m
print ('drugi korijen iz 7 je: {}'.format(m.sqrt(7))) # umjesto math.sqrt


drugi korijen iz 7 je: 2.6457513110645907

In [19]:
from math import sqrt
print ('drugi korijen iz 9 je: {}'.format(sqrt(9))) # bez korištenja točke


drugi korijen iz 9 je: 3.0

Ukoliko želimo učitati sve iz biblioteke koristimo *.


In [20]:
from math import *
x = cos(2 * pi)
print(x)


1.0

Varijable se definiraju sa znakom za jednakost =.


In [21]:
masa = 90.5
godine = 50
ime = 'Osoba #1' # može i ime = "Osoba #1"

Brojevi


In [22]:
2/3 # U pythonu 2 rezultat bi bio 0


Out[22]:
0.6666666666666666

In [23]:
2//3


Out[23]:
0

In [24]:
float(2)/3


Out[24]:
0.6666666666666666

In [25]:
8%4


Out[25]:
0

In [26]:
5**2


Out[26]:
25

In [27]:
# kompleksni brojevi
x = 1.0 - 1.0j
print(x.real, x.imag)


1.0 -1.0

Razlomci


In [28]:
from fractions import Fraction
print(Fraction(6,7))
print(Fraction(2))
print(Fraction('-3/7'))
print(Fraction(2.25))


6/7
2
-3/7
9/4

Stringovi, liste, rječnici


In [29]:
s = "Hello world"
print(len(s))
s2 = s.replace("world", "test")
print(s2)
print(s[:5])
print(s[6:])
print(s[::2])
s3 = "Hello" + "world"


11
Hello test
Hello
world
Hlowrd

In [30]:
l = [1,2,3,4]
print(l[1:3])
l2 = [1, 'a', 1.0, 1-1j]
l3 = [1, [2, [3, [4, [5]]]]]
s2=list(s2)
s2


[2, 3]
Out[30]:
['H', 'e', 'l', 'l', 'o', ' ', 't', 'e', 's', 't']

In [31]:
l = []
l.append("A")
l.append("d")
l.append("d")
l


Out[31]:
['A', 'd', 'd']

In [32]:
l[1]='A'
l


Out[32]:
['A', 'A', 'd']

In [33]:
l.remove("A")
l


Out[33]:
['A', 'd']

In [34]:
# nizovi (tuples) su isto što i liste, ali ne mogu se mijenjati
point = (10, 20)
point[0] = 20


---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
<ipython-input-34-678f722d1c34> in <module>()
      1 # nizovi (tuples) su isto što i liste, ali ne mogu se mijenjati
      2 point = (10, 20)
----> 3 point[0] = 20

TypeError: 'tuple' object does not support item assignment

In [36]:
# rječnici
params = {"parametar1" : 1.0,
          "parametar2" : 2.0,
          "parametar3" : 3.0,}
print("parametar2 = " + str(params["parametar2"]))


parametar2 = 2.0

Skupovi


In [37]:
korpa = ['jabuka', 'malina', 'jabuka', 'kupina', 'dunja', 'banana']
voce = set(korpa)
voce


Out[37]:
{'banana', 'dunja', 'jabuka', 'kupina', 'malina'}

In [38]:
'dunja' in voce, 'ananas' in voce


Out[38]:
(True, False)

In [39]:
skup = {'a','b','c'}
skup


Out[39]:
{'a', 'b', 'c'}

In [40]:
{x for x in 'abracadabra' if x not in 'abc'}


Out[40]:
{'d', 'r'}

In [41]:
a = set('abracadabra')
b = set('alacazam')
print(a - b) # razlika
print(a | b) # unija
print(a & b) # presjek
print(a ^ b) # simetrična razlika


{'r', 'd', 'b'}
{'r', 'd', 'a', 'l', 'b', 'c', 'z', 'm'}
{'c', 'a'}
{'r', 'd', 'l', 'b', 'z', 'm'}

Funkcije


In [42]:
def zero():
    return 0
zero()


Out[42]:
0

In [43]:
def powers(x):
    """
    Potencije od x.
    """
    return x ** 2, x ** 3, x ** 4
x2, x3, x4 = powers(3)
print(x3)


27

Blokovi se u Pythonu označavaju uvlačenjem koda. Praznine imaju sintaktičku ulogu.

Kompliciraniji primjer.


In [44]:
from skimage import io, transform
def thumbnail(slika, sirina=100, ime='thumb.png'):
    """Ova funkcija kao ulazne parametre prima sliku, širinu (s defaultnom vrijednošću 100) 
    te ime za thumbnail (s defaultnom vrijednošću thumb.png)"""
    visina = int(slika.shape[1] * float(sirina) / slika.shape[0])
    io.imsave(ime, transform.resize(slika,(sirina, visina)))

In [45]:
astro = data.astronaut()
astro.view()


Out[45]:
array([[[154, 147, 151],
        [109, 103, 124],
        [ 63,  58, 102],
        ..., 
        [127, 120, 115],
        [120, 117, 106],
        [125, 119, 110]],

       [[177, 171, 171],
        [144, 141, 143],
        [113, 114, 124],
        ..., 
        [127, 118, 112],
        [124, 115, 108],
        [121, 116, 105]],

       [[201, 194, 193],
        [182, 178, 175],
        [168, 165, 164],
        ..., 
        [128, 120, 117],
        [126, 116, 112],
        [124, 114, 109]],

       ..., 
       [[186, 170, 176],
        [186, 170, 177],
        [183, 168, 170],
        ..., 
        [  0,   0,   0],
        [  0,   0,   1],
        [  0,   0,   0]],

       [[183, 169, 170],
        [182, 167, 171],
        [185, 164, 176],
        ..., 
        [  0,   0,   1],
        [  1,   1,   1],
        [  0,   0,   0]],

       [[184, 167, 172],
        [183, 165, 169],
        [180, 162, 171],
        ..., 
        [  0,   0,   0],
        [  1,   1,   1],
        [  0,   0,   0]]], dtype=uint8)

In [46]:
plt.imshow(astro);



In [47]:
io.imsave('astro.png', astro)

In [48]:
astro_s_diska = io.imread('astro.png')
thumbnail(astro_s_diska)


/projects/anaconda3/lib/python3.5/site-packages/skimage/util/dtype.py:110: UserWarning: Possible precision loss when converting from float64 to uint8
  "%s to %s" % (dtypeobj_in, dtypeobj))

In [49]:
Image('thumb.png')


Out[49]:

In [50]:
thumbnail(astro_s_diska, sirina=200, ime='thumb2.png')
Image('thumb2.png')


/projects/anaconda3/lib/python3.5/site-packages/skimage/util/dtype.py:110: UserWarning: Possible precision loss when converting from float64 to uint8
  "%s to %s" % (dtypeobj_in, dtypeobj))
Out[50]:

In [51]:
# Anonimne funkcije
f1 = lambda x: x**2

# Gdje anonimne funkcije mogu biti korisne?
map(lambda x: x**2, range(-3,4))
list(map(lambda x: x**2, range(-3,4)))


Out[51]:
[9, 4, 1, 0, 1, 4, 9]

Funkcije su objekti kao i svi drugi.


In [52]:
def linearna (a,b):
    def rezultat (x):
        return a*x+b
    return rezultat
f=linearna(0.5,2)
f(1e3)


Out[52]:
502.0

In [53]:
def kompozicija(f,g):
    return lambda x: f(g(x))
def g(x):
    return x**2
fg = kompozicija(f,g)
fg(1)


Out[53]:
2.5

In [54]:
from operator import pow
pow(5,2)


Out[54]:
25

In [55]:
def parcijalno(metoda, parametar):
    return lambda x: metoda(x,parametar)
h = parcijalno(pow,2)
h(20)


Out[55]:
400

In [56]:
#kompozija proizvoljnog broja funkcija
from functools import reduce
def compose(*funcs): 
    return lambda x: reduce(lambda v, f: f(v), reversed(funcs), x)
ffgh = compose(f,f,g,h)
ffgh(1)


Out[56]:
3.25

Kontrola toka


In [57]:
izjava1 = False
izjava2 = False

if izjava1:
    print("izjava1 je True")
    
elif izjava2:
    print("izjava2 je True")
    
else:
    print("izjava1 i izjava2 su False")


izjava1 i izjava2 su False

In [58]:
izjava1 = izjava2 = True

if izjava1:
    if izjava2:
        print("i izjava1 i izjava2 su True")


i izjava1 i izjava2 su True

In [59]:
for x in range(4): #range počinje od 0
    print(x)


0
1
2
3

In [60]:
for word in ["Znanstvenici", "vole", "koristiti", "python"]:
    print(word)


Znanstvenici
vole
koristiti
python

In [61]:
for key, value in params.items():
#   print(key + " = " + str(value)) isto ispisuje ali je teže čitati kod
    print("{} = {}".format(key,str(value)))


parametar3 = 3.0
parametar2 = 2.0
parametar1 = 1.0

In [62]:
# kreiranje liste na elegantniji način pomoću for petlje
l1 = [x**2 for x in range(0,5)]

print(l1)


[0, 1, 4, 9, 16]

In [63]:
i = 0

while i < 5:
    print(i)
    
    i = i + 1
    
print("gotovo")


0
1
2
3
4
gotovo

Klase

Smisao klasa (tj. objektnog programiranja) je da omogući grupiranje varijabli i pripadnih funkcija koje rade s njima.


In [64]:
class Point:
    """
    Jednostavna klasa koja služi za rad s točkama u Euklidskoj ravnini.
    """    
    def __init__(self, x, y):
        """
        Kreiranje točke s koordinatama x i y.
        """
        self.x = x
        self.y = y        
    def translate(self, dx, dy):
        """
        Translacija točke za dx u smjeru x-osi i dy u smjeru y-osi.
        """
        self.x += dx
        self.y += dy        
    def __str__(self):
        """
        Prikaz točke.
        """
        return("Point at [{:f}, {:f}]".format(self.x, self.y))

Funkcija oblika __imeFunkcija__ obično implementira neku standardnu metodu za odgovarajući tip podataka.

  • __init__služi za kreiranje objekta dane klase
  • __str__služi za ispisivanje objekta
  • s matematičkog aspekta je važno da se mogu definirati i matematičke operacije nad objektima, kao npr. zbrajanje, množenje,...

In [65]:
# ako smo zaboravili definirati metodu pri kreiranju klase
def add(self, other):
    return Point(self.x + other.x,self.y + other.y)
Point.__add__ = add

In [66]:
p1 = Point(1,1)
p1.translate(2,3.5)
p2 = Point (2,3.5)

In [67]:
print(p1)
print(p2)
print(p1 + p2)


Point at [3.000000, 4.500000]
Point at [2.000000, 3.500000]
Point at [5.000000, 8.000000]

Moduli

Moduli (ono što učitavamo s import) strukturom su ili datoteke (s nastavkom .py) koje na početku imaju doc string (komentar) koji opisuje modul.

Primjer modula uz korištenje Jupyter magije.


In [68]:
%%file mojmodul.py

# -*- coding: utf-8 -*-
"""
Primjer modula. Sadrži varijablu my_variable,
funkciju my_function te klasu MyClass.
"""
my_variable = 0
def my_function():
    """
    Primjer funkcije
    """
    return my_variable
class MyClass:
    """
    Primjer klase
    """
    def __init__(self):
        self.variable = my_variable        
    def set_variable(self, new_value):
        """
        Daje novu vrijednost varijabli self.variable
        """
        self.variable = new_value        
    def get_variable(self):
        return self.variable


Overwriting mojmodul.py

In [69]:
import mojmodul
help(mojmodul)


Help on module mojmodul:

NAME
    mojmodul

DESCRIPTION
    Primjer modula. Sadrži varijablu my_variable,
    funkciju my_function te klasu MyClass.

CLASSES
    builtins.object
        MyClass
    
    class MyClass(builtins.object)
     |  Primjer klase
     |  
     |  Methods defined here:
     |  
     |  __init__(self)
     |      Initialize self.  See help(type(self)) for accurate signature.
     |  
     |  get_variable(self)
     |  
     |  set_variable(self, new_value)
     |      Daje novu vrijednost varijabli self.variable
     |  
     |  ----------------------------------------------------------------------
     |  Data descriptors defined here:
     |  
     |  __dict__
     |      dictionary for instance variables (if defined)
     |  
     |  __weakref__
     |      list of weak references to the object (if defined)

FUNCTIONS
    my_function()
        Primjer funkcije

DATA
    my_variable = 0

FILE
    /projects/95cb524c-f1ac-4a19-908e-5f47c10f47c1/01. Uvod; Python/mojmodul.py


Paketi

To su kompleksniji moduli, koji su strukturirani kao dorektoriji, s datotekom __init__.py.

Lovljenje grešaka


In [70]:
raise Exception("opis greške")


---------------------------------------------------------------------------
Exception                                 Traceback (most recent call last)
<ipython-input-70-fc41adc288f0> in <module>()
----> 1 raise Exception("opis greške")

Exception: opis greške

In [71]:
try:
    # varijabla varijabla nije definirana
    print(varijabla)
except:
    print("Opa!")


Opa!

Jupyter interactive

Jupyter notebook podržava i interaktivan rad. Dokumentaciju možete ovdje pogledati. Za sada ćemo samo napraviti jedan primjer.


In [80]:
from ipywidgets import interact
import matplotlib.pyplot as plt
import networkx as nx

In [81]:
def random_lobster(n, m, k, p):
    return nx.random_lobster(n, p, p / m)

def powerlaw_cluster(n, m, k, p):
    return nx.powerlaw_cluster_graph(n, m, p)

def erdos_renyi(n, m, k, p):
    return nx.erdos_renyi_graph(n, p)

def newman_watts_strogatz(n, m, k, p):
    return nx.newman_watts_strogatz_graph(n, k, p)

def plot_random_graph(n, m, k, p, generator):
    g = generator(n,m, k, p)
    nx.draw(g)
    plt.show()

In [79]:
interact(plot_random_graph, n=(2,30), m=(1,10), k=(1,10), p=(0.0, 1.0, 0.001),
        generator={'lobster': random_lobster,
                   'power law': powerlaw_cluster,
                   'Newman-Watts-Strogatz': newman_watts_strogatz,
                   u'Erdős-Rényi': erdos_renyi,   });


<matplotlib.figure.Figure at 0x7f99d1120ef0>

In [75]:
from verzije import *
from IPython.display import HTML
HTML(print_sysinfo()+info_packages('matplotlib,IPython,scikit-image,networkx'))


Out[75]:
Python verzija3.5.3
kompajlerGCC 4.8.2 20140120 (Red Hat 4.8.2-15)
sustavLinux
broj CPU-a8
interpreter64bit
matplotlib verzija1.5.3
IPython verzija5.3.0
scikit-image verzija0.12.3
networkx verzija1.11

Zadaci za vježbu

  1. Napišite funkciju koja prima n, a vraća listu neparnih brojeva od 1 do n.
  2. Napišite funkciju koja rješava kvadratnu jednadžbu.
  3. Napišite funkciju koja numerički računa integral funkcije koristeći trapeznu formulu $$ \int_a^b f(x)\approx \frac{h}{2}\sum_{i=1}^n (f(x_{i-1})+f(x_i)). $$ Funkcija treba ovako izgledati: trapezint(f,n,a,b)
  4. Napišite funkciju za numeričko deriviranje oblika diff(f,x,h=1e-6)