1/1
Technické výpočty rychle a jednoduše
Pro (elektro)technické a vědecké výpočty a vizualizaci dat se velice často používá prostředí MATLAB. Jde o poměrně drahý program a proto zde budu využívat program GNU Octave a programovací jazyk Python. Jak s pomocí jednoho, tak s pomocí druhého lze dosáhnout téhož. Pokusím se zde nastínit výhody a nevýhody.
MATLAB platí v této oblasti za nepsaný standard, proto se Octave drží jeho syntaxe a kód vytvořený pro MATLAB lze spustit v Octave a naopak.
Pro první pokusy nabízím následující zdrojový kód. Octave stačí
stáhnout a nainstalovat. Kód můžete interaktivně spouštět nebo použít
rovnou celý skript. Potom ho umístěte do pracovního adresáře a jako příkaz
uveďte jméno skriptu bez přípony. Určitě se vám budou hodit příkazy
pwd a cd. Octave má také český web octave.cz.
- prvni.m
1 #!/usr/bin/octave -q 2 3 % vytvořím hodnoty (vektor) x 4 x = [1, 2, 3, 4, 5, 6, 7 ] 5 6 % vytvořím hodnoty (vektor) y = x*x 7 y = x .^ 2 8 % nebo 9 y = x .* x 10 11 figure(1) 12 plot(x,y) 13 stem(x,y) 14 plot(x,y,'*') 15 grid('on') 16 title('y=f(x)') 17 xlabel('x') 18 ylabel('y=x^2') 19 20 % to, co končí středníkem se nevypisuje 21 zacatek = 0; 22 konec= 2*pi; 23 pocet = 100; 24 x = linspace(zacatek, konec, pocet); 25 y = sin(x); 26 27 figure(2) 28 plot(x,y) 29 grid('on') 30 title('y=f(x)') 31 xlabel('x') 32 ylabel('y=sin(x)') 33 34 pause()- `--> stáhnout
Další možností je použít programovací jazyk Python. Python je kompletní vysokoúrovňový programovací jazyk se vším všudy. (Více v sekci Python.) Čím hlouběji do něj pronikám, tím jsem si více a více jistý, že Python je ten pravý jazyk, jak pro technické výpočty a vizualizaci dat, tak pro výuku programování na střední škole. Python je maximálně jednoduchý a existuje pro něj řada kvalitních knihoven, umožňujících jednoduše provádět numerické a technické výpočty a vykreslování. Vše opět ala Matlab.
Aby bylo možné sílu Pythonu využít je nutné nainstalovat nejen Python jako takový, ale i výše zmíněné knihovny. Na mnou doporučovaném operačním systému Debian (v Ubuntu) Linux zadáte příkaz
# aptitude install python-numpy python-scipy python-matplotlib
spyder
případně
# aptitude install python-doc python-matplotlib-doc
Na Windows lze Python, všechny potřebné knihovny a související programy poměrně jednoduše nainstalovat jako distribuci Python(x,y). Stačí opět stáhnout a instalovat. Při instalaci je možné vybrat konkrétní balíčky, které se budou instalovat. Já jsem volil autory doporučovaný seznam balíčků (tedy výchozí nastavení). Pro samotnou práci doporučuji prostředí spider.
- prvni.py
1 #!/usr/bin/python 2 # -*- coding: utf8 -*- 3 # Soubor: prvni.py 4 # Datum: 19.04.2011 21:34 5 # Autor: Marek Nožka, nozka <@t> spseol <d.t> cz 6 # Licence: GNU/GPL 7 # Úloha: numerické výpočty v pythonu -- první seznámení a pokusy 8 # 9 10 # import knihoven 11 import pylab as lab 12 from pylab import pi 13 14 ############# první graf ################### 15 16 # vytvořím hodnoty (vektor) x 17 x = lab.array([1, 2, 3, 4, 5, 6, 7 ]) 18 # vytvořím hodnoty (vektor) y = x*x 19 y = x ** 2 # operátor ** je mocnina 20 #nebo 21 y = x * x 22 23 lab.figure(1) 24 lab.plot(x,y) 25 lab.stem(x,y) 26 lab.plot(x,y,'*') 27 lab.grid(True) 28 lab.title('y=f(x)') 29 lab.xlabel('x') 30 lab.ylabel('y=x^2') 31 32 ############# druhý graf ################### 33 zacatek = 0; 34 konec= 2*pi; 35 pocet = 100; 36 x = lab.linspace(zacatek, konec, pocet); 37 #x = lab.arange(zacatek, konec, krok); 38 39 y = lab.sin(x); 40 41 lab.figure(2) 42 lab.plot(x,y) 43 lab.grid(True) 44 lab.title('y=f(x)') 45 lab.xlabel('x') 46 lab.ylabel('y=sin(x)') 47 48 lab.show()- `--> stáhnout
Časový průběh napětí
Další jednoduchý příklad: Chci zobrazit časový průběh napětí o frekvenci 50Hz, efektivní hodnotě 230V v čase od 18ms do 64ms.
- sitove-napeti.py
1 #!/usr/bin/env python 2 # -*- coding: utf8 -*- 3 # Soubor: sitove-napeti.py 4 # Datum: 20.04.2011 20:05 5 # Autor: Marek Nožka, nozka <@t> spseol <d.t> cz 6 # Licence: GNU/GPL 7 # Úloha: časový průběh síťového napětí 8 # 9 10 import pylab as l 11 12 f = 50 13 Um = 230 * l.sqrt(2) # sqrt = odmocnina 14 15 t = l.linspace(18, 64, 200) 16 17 u = Um * l.sin(2*l.pi*f*t*1e-3) 18 19 l.plot(t,u) 20 l.grid(True) 21 l.xlabel('t[ms]') # všimněte si, že čas je v ms 22 l.ylabel('u(t)[V]') 23 l.text(20, 330, '$u(t)=sin(2\pi ft)$', {'fontsize' : 20}) 24 25 l.show() 26- `--> stáhnout
Výkon střídavého proudu
Na tomto místě chci ukázat, jak lze s pomocí jednoho nástroje -- Pythonu, vytvořit klasickou aplikaci i interaktivní webovou aplikaci.
Online
Offline
- vykony.py
1 #!/usr/bin/python 2 # -*- coding: utf8 -*- 3 # Soubor: vykony.py 4 # Datum: 20.04.2011 19:55 5 # Autor: Marek Nožka, nozka <@t> spseol <d.t> cz 6 # Licence: GNU/GPL 7 # Úloha: okmžitá hodnota výkonu je dána součinem okamžitých hodnot 8 # napětí a proudu 9 10 11 # import modulů 12 import pylab as l 13 from pylab import pi # nechce se mi psát l.pi 14 15 # fázový posuv mezi napětím a proudem 16 angle = 50 17 18 # čas -- arange (od, do, krok ) 19 t = l.arange(-0.1,1.1,0.001) 20 21 # napětí 22 u = l.sin(2*pi*2*t) 23 # proud 24 i = l.sin(2*pi*2*t+(angle*pi/180)) 25 26 l.figure() # nový obrázek 27 l.grid(True) # mřížka zapnuta 28 # vynesení grafů 29 l.plot(t, u,label="u") 30 l.plot(t, i,"--", label="i") 31 l.plot(t, u*i, label="p", linewidth=2) 32 l.axis([-0.1,1.1,-1.2,1.2]) # meze os 33 # popisky 34 l.title("Okamzite hodnoty napeti, proudu a vykonu") 35 l.xlabel("t[s]") 36 l.ylabel("u[V], i[A], p[W]") 37 l.legend() # legendu zobrazit 38 l.show() # vykreslit graf 39 40 41 # vim:nospell:- `--> stáhnout
Sériový rezonanční obvod
- sro.py
1 #!/usr/bin/env python 2 # -*- coding: utf8 -*- 3 # Soubor: sro.py 4 # Datum: 21.04.2011 10:15 5 # Autor: Marek Nožka, nozka <@t> spseol <d.t> cz 6 # Licence: GNU/GPL 7 # Úloha: frekvenční charakteristika sériového rezonančního obvodu 8 # 9 10 import pylab as l 11 import numpy as n 12 from pylab import pi 13 14 from matplotlib import rcParams #, rc 15 #rc('font',**{'family':'sans-serif','sans-serif':['Helvetica']}) 16 ## for Palatino and other serif fonts use: 17 #rc('font',**{'family':'serif','serif':['Palatino']}) 18 19 rcParams['text.usetex']=True 20 rcParams['text.latex.unicode']=True 21 L = .8e-3 22 C = 42e-12 23 24 #### první graf + příprava ###### 25 R = 10 26 27 # rezonanční kmitočet 28 f0 = 1.0/( 2*pi*l.sqrt(L*C) ) 29 # jakost 30 Q = 2*pi*f0*L/R 31 # šířka pásma 32 B = f0/Q 33 34 print "f0 = {0:.5} kHz".format(f0/1000.) 35 print "Q = {0:.5}".format(Q) 36 print "B = {0:.5} kHz".format(B/1000) 37 38 39 f_min = f0 - 2*B 40 f_max = f0 + 2*B 41 f= l.linspace(f_min, f_max, 10000) 42 w=2*pi*f 43 44 45 Z = R + 1j*w*L - 1j/(w*C) 46 l.plot(f/1000,n.abs(Z),label="Q = {0:8.5}".format(Q)) 47 48 ###### druhý graf ############ 49 R=20 50 Q = 2*pi*f0*L/R 51 Z = R + 1j*w*L - 1j/(w*C) 52 l.plot(f/1000,n.abs(Z),label="Q = {0:8.5}".format(Q)) 53 54 ###### třetí graf ############ 55 R=50 56 Q = 2*pi*f0*L/R 57 Z = R + 1j*w*L - 1j/(w*C) 58 l.plot(f/1000,n.abs(Z),label="Q = {0:8.5}".format(Q)) 59 60 # nadpis 61 l.title(u'Frekvenčni závislost modulu impedance SRO') 62 # meze os 63 l.axis([f_min/1000,f_max/1000,0,max(n.abs(Z))]) 64 # popisky os 65 l.xlabel(r'$f[kHz]$',{'fontsize' : 12}) 66 l.ylabel(r'$|Z|[\Omega]$',{'fontsize' : 12}) 67 # legenda a mřížka 68 l.legend() 69 l.grid(True) 70 # ukaž se! 71 l.show()- `--> stáhnout
Paralelní rezonanční obvod
Tento program zde vystavuji pouze jako polotovar. Pohrajte si s hodnotami R, L, C tak, aby výstup vypadal přibližně takto.
- pro.py
1 #!/usr/bin/env python 2 # -*- coding: utf8 -*- 3 # Soubor: pro.py 4 # Datum: 21.04.2011 10:15 5 # Autor: Marek Nožka, nozka <@t> spseol <d.t> cz 6 # Licence: GNU/GPL 7 # Úloha: frekvenční charakteristika paralelního rezonančního obvodu 8 # 9 10 import pylab as l 11 import numpy as n 12 from pylab import pi 13 14 L = 80e-3 15 C = 42e-12 16 R = 0.1 17 18 f0 = (1.0/2.0/pi) * l.sqrt( 1.0/L/C - (R/L)**2 ) 19 # jakost 20 Q = 2*pi*f0*L/R 21 # šířka pásma 22 B = f0/Q 23 24 print "f0 = {0:.5} kHz".format(f0/1000.) 25 print "Q = {0:.5}".format(Q) 26 print "B = {0:.5} kHz".format(B/1000) 27 28 29 f_min = 86e3 30 f_max = 87e3 31 f= l.linspace(f_min, f_max, 1000) 32 w=2*pi*f 33 34 ########## první graf ################ 35 Y = ( 1j*w*C ) + ( 1./ ( R + 1j*w*L ) ) 36 Z = 1./Y 37 l.plot(f/1000,n.abs(Z)/1e6,label="Q = {0:8.5}".format(Q)) 38 grafmax = (max(n.abs(Z))) 39 40 ########## druhý graf ################ 41 R=1 42 f0 = (1.0/2.0/pi) * l.sqrt( 1.0/L/C - (R/L)**2 ) 43 Q = 2*pi*f0*L/R 44 Y = ( 1j*w*C ) + ( 1./ ( R + 1j*w*L ) ) 45 Z = 1./Y 46 l.plot(f/1000,n.abs(Z)/1e6,label="Q = {0:8.5}".format(Q)) 47 48 ########## první graf ################ 49 R=10 50 f0 = (1.0/2.0/pi) * l.sqrt( 1.0/L/C - (R/L)**2 ) 51 Q = 2*pi*f0*L/R 52 Y = ( 1j*w*C ) + ( 1./ ( R + 1j*w*L ) ) 53 Z = 1./Y 54 l.plot(f/1000,n.abs(Z)/1e6,label="Q = {0:8.5}".format(Q)) 55 56 ######### vzhled ################## 57 l.title('Frekvencni zavislost modulu impedance PRO') 58 # meze os 59 l.axis([f_min/1000,f_max/1000, 0, 1.1*grafmax/1e6]) 60 # popisky os 61 l.xlabel(r'$f[kHz]$',{'fontsize' : 20}) 62 l.ylabel(r'$|Z|[M\Omega]$',{'fontsize' : 20}) 63 # legenda a mřížka 64 l.legend() 65 l.grid(True) 66 # ukaž se! 67 l.show()- `--> stáhnout
