[Edellinen] [Seuraava] [Sisällys]

Ohjelmointi

Lyhyen oppaan ohjelmointi-osa (Hiukan sekavassa tilassa toistaiseksi)
Laajemman oppaan loogiset ja vertailuoperaattorit

Tästä lyhyen oppaan kohdasta voit tutustua/palauttaa mieleesi m-tiedostojen perusasiat. Olemme käyttäneet tässäkin oppaassa m-tiedostoon kirjoitettuja komentojonoja jo edellä. Kun siirrymme ohjelmointiin, kirjoitamme funktiota (ohjelmia) funktio-m-tiedostoihin, kuten yllä viitatussa kohdassa esitetään.

Käsittellä "ohjelmointi" tarkoitetaan yleensä kahta asiaa:

  1. Komentojonon suoritusta, jota tarvittaessa hallitaan sopivilla ehdoilla säädeltävillä ohjausrakenteilla.
  2. Aliohjelmien (funktioiden) kirjoittamista. Olennaista on parametrien (argumenttien) välittäminen ylemmän tason aliohjelmasta, "pääohjelmasta", jollaisena toimii Matlab-tulkin komentotila tai m-skripti.
Ohjausrakenteita Matlabissa on Matlabin tapauksessa tulee välttää turhaa ohjausrakenteiden (erityisesti pitkien for-luuppien) käyttöä, usein ne voidaan korvata "vektoriajattelulla", joka on sekä ohjelman suorituksen että ohjelmoinnin kannalta usein huomattavasti tehokkaampaa.
Toki on paljon tilanteita, joissa ohjausrakenteita ei voida välttää, eikä "puhdasoppisuudessa" pidä pyrkiä äärimmäisyyksiin (vaikka joskus mieli tekee).

Esimerkkiohjelma 1

Tämä on itse asiassa "vain" skripti, mitään parametrejä ei välitetä. Ts. kyseessä on komentojono, jollaisia olemme jo nähneet ja tehneet.
Kuitenkin tietoa ohjelmalle voidaan välittää vuorovaikutteisesti kyselemällä. Tällainen tapa on enemmän leikkimistä varten, mutta voi toki joskus olla hyötykäytössäkin.
Lisäksi näemme tässä merkkijonoilla operointia. 
% Tiedostossa tervehdys.m oleva koodi kysyy nimesi
% ja vastaa tervehdyksellä ja päivämäärällä.
% Samalla esitellään vähän merkkijono-operointia.
disp('Hei, kuka olet?')
nimi=input(['Kirjoita nimesi hipsukoissa: ',' Esim: ''Hemuli''  '] );
pvm=date;
vastaus=['Hei ',nimi,'. Tänään on ',pvm,' .   '];
disp(vastaus)


 Merkkijonoista lyhyesti  


 Merkkijo on vektori,
joka koostuu ASCII-merkeistä. Esim: 
  >> mjono='Olen merkkijono'.


Merkkijonoja liitetään peräkkäin, kuten vektoreita. Esim:

>> mjono='Olen merkkijono'
mjono =
Olen merkkijono
>> mjono2=[mjono,'. ', mjono,'.']
mjono2 =
Olen merkkijono. Olen merkkijono.


Jos merkkijonoon halutaan hipsukka, se kahdennetaan.  


 Esimerkkejä omista funktioista 


Lyhyessä oppaassa on joukko esimerkkejä alkaen keskiarvo-funktiosta. Lisätään 
yksinkertaisten esimerkkien valikoimaa aloittamalla Pythagoraan lauseesta:


 Tehtävä: 
Kirjoita funktio-m-tiedosto hypo.m, jonka otsikkorivi voisi olla:

 function c=hypo(a,b)


(Huom: Ennen editorin avaamista on syytä valita Matlab-istuntoikkunasta 
polku omaan työhakemistoon.) 


Argumenttina annetaan samankokoiset matriisit a ja b, tuloksena palautetaan
edelleen samankokoinen matriisi c, jossa on vastaavien kateettiarvojen hypotenuusa-arvot. 

Kirjoita alkukommentit ja sisällytä niihin muutama esimerkki, kuten

 % hypo(3,4)
 % C=hypo((1:10),(10:-1:1))
 % hypo(hilb(3),inv(hilb(3)))


Tässäkin tapauksessa suurin osa funktiosta on kommentteja, varsinaista 
koodiahan on tasan yksi lyhyt rivi.


Kun ryhdyt aikaansaannostasi testaamaan,
 suorita ensin  >> help hypo , testaile sitten ainakin help-esimerkit läpi.

Huom: Anna komento  >> which hypo  ja myös which hypot . (Tässä syy, miksi en ehdottanut jälkimmäistä nimeä.)






Komentorivifunktiot, funktio-osoittimet (vanhanaikaiset inlinet)


Lyhyen oppaan inline-kohdassa
käsitellään näitä. 

Kyseessä on lyhyet, komentorivillä määriteltävät funktiot. Uusissa versioissa
ei kannata enää aktiivisesti opetella muuta kuin funktio-osoitin ("function handle")-tapaa.


Jos haluaisimme määritellä vaikkapa funktion f(x) = esin(x), 
voisimme m-tiedoston sijasta kirjoittaa komentoriville:

>> f=@(x) exp(sin(x))
>> ezplot(f,[-pi,pi])


Huomaa: Funktion ezplot kutsussa f ei ole merkkijono, kuten ennen,
vaan "osoitin-tyyppinen olio". (Vanha merkkijonotyylikin toimii, mutta silloin
funktion määrittely on tehtävä myös "vanhanaikaisesti", inline:lla.)


Jos funktiomääritystä ei muualla tarvita, se voidaan antaa myös suoraan
nimettömänä määrityksenä kutsuvalle funktiolle, tässä ezplot:lle

>> ezplot(@(x) exp(sin(x)),[-pi,pi])


Määrittely  @(x) f(x)  vastaa matemaattista merkintää:
 x->f(x) , joka itse asiassa on Maple-ohjelmassa otettu vastaavan
määrittelyn syntaksiksi. 




 Lue: AT(x) arvo(x)  




Argumentteja voi tietysti olla n kappaletta:

>> g=@(x,y,z) x.^2+y.^2+z.^2
g = 
    @(x,y,z)x.^2+y.^2+z.^2





 Yhden rivin iteraatiot, nuolinäppäinalgoritmit 


 Motto: Pieni on kaunista  ( Vaikka nykyisin ei aina siltä
näytä)


Jatketaan funktio-osoittimien parissa. 

Varsin näppärää, nokkelaa, opettavaista,havainnollistavaa,
vuorovaikutteista (ja mitä vielä) on useissa tilanteissa kirjoittaa yhdelle riville iteraatioaskel ja siihen 
liittyvät muuttujien päivitykset. Askelta voidaan toistaa nuoli-ylös/ENTER-yhdistelmällä. Alkuarvon hylkäys/hyväksymispäätös ja lopetusehto määräytyvät 
vuorovaikutteisesti.



 
Huom! Uusissa Matlab:n versioissa on välineitä, joiden ansiosta on 
mukavaa ja suositeltavaa siirtyä entistä
enemmän työskentelemään komentoikkunan sijasta komentoja sisältävään m-tiedostoon, skriptiin. Nuolinäppäinalgoritmeille saadaan silloin "jalostuneempi" muoto.
Lyhyesti: Kahden kaksoisprosenttimerkin (%%) väliin jäävä m-tiedoston kappale 
suorittuu, kun osoitin viedään johonkin kohtaan kappaletta ja painetaan:
CTR-ENTER. Toki joissakin tilanteissa vanhan hyvän ajan nuolinäppäintyyli
voi olla nopein ja näppärinkin. Sitä nyt tässä esittelemme.




Tyypillinen algoritmi on Newtonin menetelmä epälineaarisen yhtälön 
 ratkaisemiseksi. Menetelmän ominaisuuksiin kuuluu nopea (kvadraattinen) 
suppeneminen silloin, kun se toimii. Ongelmatapauksissa menetelmä taas 
huitelee äärettömyyden rajamaille tai jää ikuisesti kiertelemään kehää.

Kaikessa lyhykäisyydessään menetelmän impelementointi käy tällä funktiolla:


  >> newt1=@(f,df,x) x-f(x)/df(x)


Kyseessä on siis yhtälön f(x)=0 juuren haku. Argumentteina on funktio f, f:n derivaatta df ja alkuarvaus x.

Testataanpa helpolla esimerkillä. Haetaan sin:n nollakohtaa:

>> f=@sin,df=@cos
f = 
    @sin
df = 
    @cos              % sin:n derivaatta
>> x0=1.2             % Alkuarvaus
x0 =
    1.2000
>> x0=newt1(f,df,x0)  % 1. iteraatio
x0 =
   -1.3722
>> x0=newt1(f,df,x0)  % Nuoli ylös ja ENTER: 2. iteraatio
x0 =
    3.5956
>> x0=newt1(f,df,x0)  % ...
x0 =
    3.1076
>> x0=newt1(f,df,x0)
x0 =
    3.1416
>> format long
>> x0=newt1(f,df,x0)
x0 =
   3.141592653589793
>> x0=newt1(f,df,x0)
x0 =
   3.141592653589793  % Ei muutu Matlabin laskentatarkkuudella,'
                      %    turha jatkaa.


 Huom! Voitaisiin tietysti kirjoittaa suoraan: 
 >> x0=newt1(@sin,@cos,x0)  


Jos haluaisimme tutkia lähemmin iteraatiojonon ominaisuuksia, tarvittaisiin 
vektori, johon jonon alkioita kumuloitaisiin. Sekin hoituu nuolinäppäintekniikalla:

>> x0=1.3;X=[];   % Kumuloidaan vektoriin X, joka alustetaan tyhjäksi.
>> X=[X x0],x0=newt1(@sin,@cos,x0)
X =
    1.3000
x0 =
   -2.3021
>> X=[X x0],x0=newt1(@sin,@cos,x0)
X =
    1.3000   -2.3021
x0 =
   -3.4166
>> X=[X x0],x0=newt1(@sin,@cos,x0)
X =
    1.3000   -2.3021   -3.4166
x0 =
   -3.1344
>> X=[X x0],x0=newt1(@sin,@cos,x0)
X =
    1.3000   -2.3021   -3.4166   -3.1344
x0 =
   -3.1416
>> X=[X x0],x0=newt1(@sin,@cos,x0)
X =
    1.3000   -2.3021   -3.4166   -3.1344   -3.1416
x0 =
   -3.1416  


Huomattiin, että pieni alkuarvon muutos vei tuloksen pi:n sijasta arvoon -pi.
Jos jatketaan vielä pari kierrosta ja halutaan tutkia iteraatiojon 
peräkkäisten termien erotusta, niin mikäs siinä, funktio diff on 
työkalupakissamme:

>> format long
>> diff(X)'
ans =
  -3.602102447967979
  -1.114488696328861
   0.282146985976172
  -0.007148617036663
   0.000000121767538
                   0


Kun suppeneminen pääsee vauhtiin, oikeiden numeroiden lukumäärä likimain kaksinkertaistuu joka askeleella, sitä se kvadraattisuus on.


Lisää nuolinäppäinalgoritmeja nähnemme jatkossa. Samoin eri versioita teemasta 
Newtonin menetelmä.



Ohjausrakenteet 



 if,else,elseif,end 



>> help if
 IF Conditionally execute statements.
    The general form of the IF statement is
 
       IF expression
         statements
       ELSEIF expression
         statements
       ELSE
         statements
       END
 
    The statements are executed if the real part of the expression 
    has all non-zero elements. The ELSE and ELSEIF parts are optional.
    Zero or more ELSEIF parts can be used as well as nested IF's.
    The expression is usually of the form expr rop expr where 
    rop is ==, <, >, <=, >=, or ~=.
 
    Example
       if I == J
         A(I,J) = 2;
       elseif abs(I-J) == 1
         A(I,J) = -1;
       else
         A(I,J) = 0;
       end



Koodia voit testata kirjoittamalla se m-tiedostoon, antamalla i:lle ja j:lle 
arvoja Matlab-istunnossa ja ajamalla m-tiedosto. (Toki koodin voi kirjoittaa
 suoraan istuntoonkin, mutta näinkin monen rivin pläjäys on tuskaista 
korjailtavaa ja uudelleen ajettavaa.)


Huom! Yllä oleva tapa sisennyksineen on yleensä suositeltava.
Lyhyissä if-lauseissa on joskus mukava kirjoittaa samalle riville, se
on sallittua, mutta silloin on käytettävä erotinta, joko , tai ;
(Tällaisia on helppo kirjoittaa suoraan istuntoonkin.)


Esim:

>> x=3;y=2;if x < y, pieni=x; else pieni=y, end
pieni =
        2


Seuraavan voit leikata/liimata suoraan istuntoon, koska siinä ei ole muuteltavia parametreja. Samalla tulee pieni merkkijonharjoitelma, funktio 
num2str on tällaisissa tarpeen. (Nämä ovat usein hyödyllisiä niksejä
kuvien otsikko-ym. teksteissä.) 

Tässä selvitetään laskemalla, kumpi on suurempi: e2 vai 2e 
ja kirjoitetaan opettavaisesti näytölle. Teepä, jos huvittaa! 

e=exp(1);
if 2^e > e^2
  disp(['2^e=',num2str(2^e),' siis suurempi,koska ','e^2=',num2str(e^2)])
else
  disp(['e^2=',num2str(e^2),' siis suurempi,koska ','2^e=',num2str(2^e)])
end





 for i=1:n, komennot, end 




Esimerkki, for ja if:


function [y] = askel(t)
%askel -  Heavisiden yksikköaskelfunktio
%   t on vektori
%   Funktion arvo =1, jos t(k) >=0, muuten 0.
% Esim: 
%      askel(-10:10)
n=length(t);
y=zeros(1,n);
for k=1:n
    if t(k)>=0
        y(k)=1;
    end    
end



Kuten edellä todettiin, usein voidaan ohjausrakenteita välttää 
vektoriajattelulla. Tässä on dramaattinen esimerkki: 
Yllä oleva funktio voidaan kaikessa lyhykäisyydessään yhtä hyvin
määritellä  näin:


>> vaude=@(t) t>=0 ;



Testaa ja selvitä, miksi tämä toimii!



 while-rakenne  



doc while antaa mm. tyyppiesimerkin koodista, jolla lasketaan
"kone-epsilon", eli liukulukulaskennan suhteellinen tarkkuusraja. Kyseessä on
pienin positiivinen liukuluku eps, jolla 1+eps > 1.
Jotta laskennassa ei suotta muutettaisi Matlabin omaa eps-muuttujaa (oikeammin 
funktiota, joka käyttäjälle näkyy muuttujana), lisään nimen perään suomalaisittain i:n.

epsi = 1;
while (1+epsi) > 1
    epsi = epsi/2;
end
epsi = epsi*2



Toinen esimerkki while-rakenteesta:

Lyhyessä oppaassa on tyyppiesimerkki  iteratiivisesta algoritmista, neliöjuuren laskennasta Newtonin menetelmällä.

 *** tee a name tälle kohdalle lyhyeen *** 


Todettakoon, että ääretön ohjelmasilmukka on helppo toteuttaa:

       while 1,...,end


Ehtona on aina 1, eli totuusarvo 'true'. Poistumisehto voidaan rakentaa 
vaikka if-lauseen ja break-käskyn avulla. (Tai GUI-ohjelmoinnissa
napin painalluksella.)




Newtonin menetelmä kokonaisena ohjelmana:


Ohjelman ydin voisi mennä näin:


   while abs(x-xedell) > eps*abs(x)   % Testataan suhteellista virhettä.
      xedell=x;
      x=x-f(x)/df(x);
   end



 Tehtävä:  Kirjoita tämä täydelliseksi ohjelmaksi kommentteineen
m-tiedostoon. Palauta tuloksena koko iteraatiojono X (Vrt. yllä olevaan 
rivi-iteraatioon.) Ota valinnaisiksi argumenteiksi toleranssi tol, jonka 
oletusarvona on eps, sekä maxiter, jonka oletus voisi olla vaikka 25.


Voit ottaa mallia lyhyen oppaan neliöjuuri 2-ohjelmasta.
X:n kumuloinnin voit tehdä kahdella eri tavalla, joko nuolinäppäin-esimerkissä 
käytetyllä: X=[X x] tai X(i+1)=... ,
kuten esim. sqrt(2)-ohjelmassa. Tämä koodinpätkä tukee paremmin edellistä,
mutta helppohan sitä on modifioida. 
Todettakoon, että muistin hallinnan kannalta jälkimmäinen tapa, jossa allokoidaan sopivan pituinen vektori etukäteen: (X=zeros(1,n)), on tehokkaampi. Pienehköissä tehtävissä tällä ei ole merkitystä.



"Projektitehtäviä"




 Funktioargumentit, alifunktiot, sisäkkäiset funktiot




 Alifunktiot 



Funktio-m-tiedosto alkaa pääfunktiolla. Sen jälkeen voi tulla yksi tai useampia
funktiomäärityksiä, alifunktioita. Ne näkyvät vain pääfunktiolle, eivät sen 
ulkopuolelle.


Esim

function maxvirhe=poly1virhe(n)
% polyvirhe   Virheapproksimaatio lineaarisessa interpolaatiossa
% poly1virhe(n) laskee lineaarisessa interpolaatiossa välillä [0,1]
% n:ssä tasavälisesti valitussa otospisteessä lasketun 
% maksimivirheen.
% Interpoloitava funktio annetaan alifunktiona f

% Muokattu Kirjan [Hig-Hig] poly1err-esimerkistä s. 148
%     19.8.2010  HA

f0=f(0);f1=f(1);
x=linspace(0,1,n);
p=f0*(1-x)+f1*x;
plot(x,f(x),x,p,x,p,'o')
maxvirhe=max(abs(f(x)-p));

% Alifunktio:
function y=f(x)
y=x.*sin(x);


Koodi on valmiina myös tässä. Ota käyttöön, kokeile, muokkaa.


Toinen, meille jo tuttu tapa olisi ottaa funktio f argumenttilistaan, jolloin
otsikkorivi olisi:

   function poly1virhe2(f,n)


ja kutsu:  >> poly1virhe2(@(x) x.*sin(x),n) 


Tässä tapauksessa jälkimmäinen tapa näyttäisi yksinkertaisemmalta, mutta 
alifunktiossa voi olla useampia kuin 1 rivi ja niitä voi olla useita.
Myös pääfunkion argumenttilista tulee näin lyhyemmäksi.


Sanotaanpa vain, että riippuu tilanteesta, molempia tapoja tarvitaan, jälkimmäistä erityisesti kutsuttaessa valmiita funktioita, joiden argumenttilistaan 
sisältyy funktioita. Valmiita kirjastofunktioitahan käyttäjä ei pääse muokkaamaan, joten alifunkioita ei niihin voi kirjoittaa.



 Sisäkkäiset ("nested") funktiot 



Kts. [Hig-Hig] 10.7 ss. 151-152


Sisäkkäiset funktiot asustavat toistensa sisällä, nimensä mukaisesti.
Niitä käytettäessä kaikki funktiomääritykset on lopetettava end-sanaan. 

Perusominaisuus:




 Sisäkkäisellä funktiolla on pääsy kaikkien sitä ympäröivien funktioiden
työtilaan, ts. jos f2 on f1:n sisällä, niin f2 voi käyttää f1:n muuttujia.








 Esimerkki sisäkkäisten funktioiden käytöstä 


Määritellään ensin derivaatan differenssiapproksimaation toteuttava
erotusosamääräfunktio  eros . (Vrt. [Hig-Hig].)



function y=eros(f,x,h)
% Erotusosamääräfunktio (fd_deriv Hig-Hig s. 145)

if nargin < 3, h=sqrt(eps); end 
      % Numeerinen derivointi on häiriöaltis tehtävä, askeleena
      % ei voi käyttää pienempää kuin sqrt(eps). (eps on aivan
      % liian pieni)
 y=(f(x+h)-f(x))/h; 



Koodi on valmiina myös tässä. 


Testataan aluksi laskemalla sinifunktion derivaattoja, joita verrataan cos-arvoihin:

>> [eros(@sin,1:10)' cos(1:10)']
ans =
    0.5403    0.5403
   -0.4161   -0.4161
   -0.9900   -0.9900
   -0.6536   -0.6536
    0.2837    0.2837
    0.9602    0.9602
    0.7539    0.7539
   -0.1455   -0.1455
   -0.9111   -0.9111
   -0.8391   -0.8391



Lyhyellä näyttötarkkuudella (4 numeroa) täsmälleen samat. Tähän asti kaikki hyvin. 

 No nyt se sisäkkäisen funktion käyttö: 


Otetaan tehtäväksi laskea yleisen 1. asteen rationaalifunktion derivaattojen
differenssiapproksimaatioita, kun kertoimet a,b,c,d ovat muuteltavia 
parametreja.


function dy=ratder(x)
% 1. asteen rationaalifunktion derivaatan differenssiapproksimaatio.
% Sisäkkäisfunktioiden valaisemiseksi.
a=1;b=1;c=1;d=-1;  % Muuteltavat parametrit
dy=eros(@rat1,x);

    function r=rat1(x)
    % Määritellään ao. rationaalifunktio
    %   (Tämä siis pääsee käyttämään yllä asetettuja arvoja a=1,...,d=-1.)
    r=(a+b*x)./(c+d*x);
    end
end




>> ratder(-2)
ans =
    0.2222
>> format rat
>> ratder(-2)
ans =
       2/9     
>> format short
>> ratder(2:10)'
ans =
    2.0000
    0.5000
    0.2222
    0.1250
    0.0800
    0.0556
    0.0408
    0.0313
    0.0247


(Tarkistettu  "Symbolic toolboxilla".)


Työtapa on siis sellainen, että avataan editoriin tiedosto ratder.m,
editoidaan riviä "Muuteltavat parametrit" ja kutsutaan (komentoikkunassa tai 
erillisessä ajotiedostossa).


Tehtäviä:
Tähän tehtäviä, joissa 

    

  1.  ajetaan ratder-esimerkkiä eri parametreillä,

  2.  otetaan joitakin muita parametreistä riippuvia funktiojoukkoja,
kuten korkeamman asteen rationaalifkt., trig.-, exp-, ym. funktioiden 
lineaarikombinaatioita, ym. ym.





 Globaalit muuttujat 


Näiden käyttöä ei suositella, mutta joskus on kätevää (kunhan käytät salaa ja 
vähän huonolla omallatunnolla).


Funktion sisäiset muuttujat ovat aina paikallisia, ellei global-komennolla toisin määrätä.


Valaistaan globaaleiden käyttöä parilla esimerkillä.


Esim 1:

function y=omafun(x)
global A,B
y=1./(A+(x-B).^2);



Esim 2:


Komentorivillä voidaan kirjoittaa:

 >> global A,B
 >> A=0.01;B=0.5;  % Muutellaan tarpeen mukaan.
 >> fplot(omafun,[0 1])



Eräs hyötykäyttö voisi olla seuraavanlainen: Ajatellaan, että kutsumme
kirjastofunktiota ratkaisija, joka puolestaan kutsuu ratkaisussa tarvittavaa
funktiota f. Esim. quad, joka integroi numeerisesti f:ää.
Haluaisimme selvittää, montako kertaa ratkaisija kutsuu f:ää.
Emme pääse "ratkaisijan" koodiin kirjoittamaan. 

Voisimme kirjoittaa oman integrandifunktiomme f näin:

  function y=f(x)
  global Laskuri
  Laskuri=Laskuri+1;
  y=... ;


Istunnossa näin:

>> global Laskuri
>> Laskuri=0;
>> q=ratkaisija(@f,a,b,...)
>> Laskuri    % Katsotaan, mikä on arvo nyt.



Globaalien muuttujien käyttö voidaan yleensä korvata anonyymeillä tai sisäkkäisillä funktioilla ja kirjastofunktioden omilla monipuolisilla mekanismeilla.


 Ohjelmointioppia esimerkkien voimalla 


Hyvä tapa oppia, on tutkia ammattilaisten kirjoittamaa koodia, Mathworks'n
omat koodit lienevät esimerkillisimpiä. Koodia voit tutkia komennolla 
type funnimi, mutta "rauhallisempaa" on 
edit funnimi.
Polkua ei tarvitse tietää, ellei halua hakea funktiota omaan mielieditoriin.
Annetaan varmuudeksi alla myös polku. Selityksenä mainitaan tähdennettävä piirre. (Lista on [Hig-Hig]-kirjasta s. 156)


















[Edellinen]
[Seuraava]
[Alkusivu]