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

| navigace |

Č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

| navigace |

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

vykony.cgi

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

| navigace |

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

obrázek

| navigace |

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

| navigace |

Licence Creative Commons Valid XHTML 1.0 Strict Valid CSS! Antispam.er.cz Blog: Tlapicka.net