1  Introducció

1.1 Introducció al python

El python és un llenguatge interpretat que permet orientació a objectes. És un dels llenguatges més usats actualment per la seva flexibilitat, potència i quantitat de paquets desenvolupats que permeten efectuar multitud d’operacions, des de càlculs numèrics o simbòlics, fins a aplicacions web o pràcticament qualsevol aplicació al mercat.

El python consta d’una consola que invocarem cridant python a la línia de comandes. De totes formes nosaltres usarem una versió un xic millorada que s’anomena ipython . Aquesta consola ens permetrà treballar de forma molt semblant amb què fem amb les comandes de Matlab, per exemple. A la consola podeu escriure:

El python és un llenguatge interpretat que permet orientació a objectes. És un dels llenguatges més usats actualment per la seva flexibilitat, potència i quantitat de paquets desenvolupats que permeten efectuar multitud d’operacions, des de càlculs numèrics o simbòlics, fins a aplicacions web o pràcticament qualsevol aplicació al mercat.

El python consta d’una consola que invocarem cridant python a la línia de comandes. De totes formes nosaltres usarem una versió un xic millorada que s’anomena ipython . Aquesta consola ens permetrà treballar de forma molt semblant amb què fem amb les comandes de Matlab, per exemple. A la consola podeu escriure:

#| eval:false

2+3 # evident...
exit() # per sortir. També CRTL-D
2*3 
2**3 # Atenció en com escrivim les potències
2/3 # Compte que quan són dos enters...
2/-3 # ...pren la part sencera!
2.0/3.0 # Operació amb float.
1234567890 # podem representar fins 2**31
type(1); type(1.0) # ens dóna el tipus d'una variable.
2+4.0j # podem representar complexos

És important poder assignar valors a variables:

x=0 # El signe = és assingació
0==1 # El signe == és comprovació
x=y=z=0 # Abusem de la notació
a,b,c=1,2,3 # Múltiples assignacions
a=3.0+3.0j 
a.real # El punt indica que és un atribut de a
a.imag
abs(a)

Hi ha més tipus de dades com llistes, tuples o diccionaris que veurem més endavant. Per ara veiem només les cadenes:

cadena='El meu nom es Quim'
cadena2=' Puig i Sadurni'
cadena+cadena2 #L'operació + a cadenes ho ajunta.
print(cadena+cadena2)
len(cadena) # Això ens dóna la longitud
nom='Quim'
edat= 32
'En %s te  %s anys' % (nom,edat)

Per mitjà de llistes podem definir vectors, matrius i objectes més complicats:

vector1=[1,2,3,4,5,6] 
vector1=range(6) # versio (quasi-)alternativa
vector1[0];  # accedim com en C
matriu= [[1,2],[3,4]]  
matriu[0][1]
mixt=['pi','es',3,.141592]
mixt[3]

1.2 El paquet scipy

El paquet scipy conté nombroses rutines numèriques de gran qualitat i comparables a les del Matlab. Està molt relacionat amb el paquet numpy.

from numpy import matrix # Del numpy nomes importem la classe matrix
from scipy.linalg import inv, det, eig  # del scipy.linalg tambe importem aquestes

A=matrix([[1,1,1],[4,4,3],[7,8,5]]) 
b = matrix([1,2,1]).transpose()    

print det(A)     
print inv(A)*b   
print eig(A)
type(A)

També ho podríem haver fet com

from numpy import * # Del numpy ho importem tot
from scipy import *  # idem

A=[[1,1,1],[4,4,3],[7,8,5]] # Usem els tipus natius del python

print det(A)     
print inv(A)   
print eig(A)
type(A) # Ara és una llista

De cara a l’eficiència numèrica és millor treballar amb els tipus matrix del scipy. Anem ara a integrar equacions diferencials de primer ordre amb el scipy. Tot i que veurem formes més sofisticades (a través de funcions on definirem el nostre camp), el següent mètode serà suficient:

from scipy.integrate import odeint
odeint(lambda y, t: y, 1, [0, 2])
help odeint # Aixi veiem totes les opcions.

El que estem fent és integrar l’equació diferencial escalar \[y'(t)= y(t), \quad y(0)=1\] Per tant el resultat serà \(y(2)=e^2\). El vector [0,2] ens indica que integrem des de 0 (a on tenim la condició inicial) fins a \(t=2\), que serà la segona component del vector. És interessant com definim funcions de forma anònima usant lambda:

producte=lambda numero1, numero2: numero1*numero2
producte(2,3)

Per dibuixar les solucions també podem demanar que ens registri

import numpy as N # aixo és pq cridarem les rutines com N.elquesigui
x=N.linspace(0,1,10) # vector de 0 a 1 amb increment (1-0)/9. 
f=odeint(lambda y, t: y, 1, x) # ens retorna el vector de y en els punts de x.

1.3 Gràfiques amb pylab

El paquet pylab és un paquet de molta qualitat que permet fer gràfiques a l’estil del Matlab, però sense haver de gastar-nos milers d’euros. Veiem un exemple fàcil:

from pylab import *
import numpy  # manera alternativa d'importació
t = arange(0.0, 10, 0.1) # manera alternativa amb N.linspace
ft = t*cos(t) # el vector ft ara conté les imatges
plot(t,ft) 
show() # ara el mostrem per pantalla

Ens apareixerà una finestra amb el dibuix i tota una sèrie d’icones on podem modificar l’àrea de visualització i/o. Per modificar els atributs i/o afegir llegendes ho fem així:

from pylab import *
import numpy  # manera alternativa d'importació
t = arange(0.0, 10, 0.1)
ft = t*cos(t) 
gt = t*sin(t) 
plot(t,ft,'-', t,gt,'ro') # Això ens ho dibuixa amb diferents estils la gt i la ft
xlabel('temps')
ylabel('poblacio') # com posem accents i caràcters especials?
title('Evolucio de $\dot x=x$, $x(t)=e^t$') # Accepta LaTeX!
show()

Podeu veure tots els estils fent

plot?

1.4 Posant-ho tot junt

Anem a representar les solucions de l’exponencial:

from numpy import *
from pylab import *
from scipy.integrate import odeint
x=arange(0,10,0.1) #valors dels temps
fx=odeint(lambda y, t: y, 1, x) # ens retorna el vector de y en els punts de x.
plot(x,fx) #dibuixem la solucio
xlabel('temps')
ylabel('poblacio') 
title('Evolucio de $\dot x=x$, $x(0)=1$') # Accepta LaTeX!
show()

Per veure millor l’evolució també podem representar per a un rang de valors inicials:

from numpy import *
from pylab import *
from scipy.integrate import odeint
x=arange(0,10,0.1) #valors dels temps
ci=[0,0.5,1] #valors de les condicions inicials
fx=odeint(lambda y, t: y, ci, x) # ens retorna el vector de y en els punts de x.
plot(x,fx) #dibuixem la solucio
xlabel('temps')
ylabel('poblacio') 
title('Evolucio de $\dot x=x$, $x(0)=1$') # Accepta LaTeX!
show()

on ara ci conté les condicions inicials i fx serà una matriu amb les ordenades corresponents a aquestes condicions inicials. Ho podem fer de forma més automatitzada com segueix:

from numpy import *
from pylab import *
from scipy.integrate import odeint
x=arange(0,10,0.1) #valors dels temps
ci=arange(0,2,0.1) # valors de les condicions inicials
fx=odeint(lambda y, t: y, ci, x) # ens retorna el vector de y en els punts de x.
plot(x,fx) #dibuixem la solucio
xlabel('temps')
ylabel('poblacio') 
title('Evolucio de $\dot x=x$, $x(0)=1$') # Accepta LaTeX!
show()

1.5 Exercicis

Feu el mateix per a l’equació logística i veieu el tipus de comportament de \[x'=x\cdot(1-x)\] depenent del valor de la condició inicial.

1.6 Referències