Seminarie MATLAB D. Deses 2012-2013
Contents 1 Voorwoord
2
2 Het 2.1 2.2 2.3
command venster Getallen, bewerkingen en variabelen . . Scripts . . . . . . . . . . . . . . . . . . Vectoren en matrices . . . . . . . . . . 2.3.1 Bewerkingen . . . . . . . . . . . 2.3.2 Stelsels oplossen . . . . . . . . . 2.3.3 Eigenwaarden en eigenvectoren 2.4 Grafieken . . . . . . . . . . . . . . . . 2.4.1 2D-Grafieken . . . . . . . . . . 2.4.2 3D-Grafieken . . . . . . . . . . 2.4.3 Animaties . . . . . . . . . . . . 2.5 Veeltermen . . . . . . . . . . . . . . . 2.6 Symbolische wiskunde . . . . . . . . .
. . . . . . . . . . . .
. . . . . . . . . . . .
. . . . . . . . . . . .
. . . . . . . . . . . .
. . . . . . . . . . . .
. . . . . . . . . . . .
. . . . . . . . . . . .
. . . . . . . . . . . .
. . . . . . . . . . . .
. . . . . . . . . . . .
. . . . . . . . . . . .
. . . . . . . . . . . .
. . . . . . . . . . . .
2 2 4 5 6 7 8 9 9 11 14 15 16
3 Programmeren 18 3.1 Functies . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19 3.2 Verdere oefeningen . . . . . . . . . . . . . . . . . . . . . . . . 21
1
1
Voorwoord
Deze cursus bestaat uit dertien uren seminarie. Tijdens deze uren is het de bedoeling dat je deze nota’s doorneemt en de voorbeelden en oefeningen onder begeleiding op computer uitvoert. Gaandeweg leer je met het wiskundepakket MATLAB werken. Na deze dertien uren volgt een evaluatie onder de vorm van een taak. Voor deze taak zal je twee verschillende stellingen uit verschillende cursussen kiezen en deze aan de hand van voorbeelden in MATLAB illustreren. Er is geen examen MATLAB . Zorg er dus voor dat je taak in orde is. Deze tekst is geschreven voor MATLAB R2010b. Na het opstarten krijg je verschillende vensters. Het centrale venster is het command venster. Hierin kan je opdrachten ingeven die MATLAB zal uitvoeren. Verder draait deze tekst eerst en vooral om wiskunde. De verschillende voorbeelden en oefeningen zullen het gebruik vergen van de leerstof uit andere cursussen. Denk daarom eerst na en begin dan pas MATLAB te gebruiken om het rekenwerk te verlichten, MATLAB zal niet voor jou denken!
2
Het command venster
2.1
Getallen, bewerkingen en variabelen
Je kan gewoon getallen en bewerkingen gebruiken. >> 2+3 ans = 5 De uitkomsten worden in de variabele ans gestoken, je kan aldus telkens het laatste antwoord verder gebruiken. >> sqrt(-36) ans = 0 + 6.0000i >> ans^2 ans = 2
-36 Opmerking 1. MATLAB gebruikt standaard complexe getallen, de uitdrukking 2i^2 is dus equivalent aan (2i)2 en wordt dus -4. Dit in tegenstelling tot 2*i^2 dat gelijk is aan 2 · (i2 ) = −2. Verder kent MATLAB ook het getal π,voor e wordt dan weer de functie exp gebruikt. >> exp(i*pi) ans = -1.0000 + 0.0000i Naast de getallen kent MATLAB ook alle nodige wiskundige functies: sin, cos, tan, sqrt, exp, sinh, cosh, tanh, ... >> sqrt(4*(cos(pi/4)+i*sin(pi/4))) ans = 1.8478 + 0.7654i Een van de andere meest gebruikte functies van MATLAB is de helpfunctie. Je kan altijd de helpfunctie gebruiken om meer informatie over een commando te krijgen. Gebruik je help commando, dan krijg je een korte beschrijving. Voor een meer uigebreide uitleg kan je doc commando gebruiken. Je kan ook kortweg doc ingeven om de help bladzijden te openen. >> help sqrt SQRT Square root. SQRT(X) is the square root of the elements of X. Complex results are produced if X is not positive. See also sqrtm, realsqrt, hypot. Overloaded methods: sym/sqrt Reference page in Help browser doc sqrt 3
MATLAB kent ook variabelen, om een waarde toe te kennen gebruik je gewoon =. In het workspace venster wordt er bijgehouden welke variabelen je al hebt gebruikt en welke waarde ze hebben. >>a=1;b=-1;c=-1; >>D=b^2-4*a*c D = 5 >> x1=(-b+sqrt(D))/(2*a),x2=(-b-sqrt(D))/(2*a) x1 = 1.6180
x2 = -0.6180 Opmerking 2. Met een komma kan je verschillende commando’s op eenzelfde lijn ingeven. Met een puntkomma kan je voorkomen dat een resultaat wordt weergegeven. Oefening 1. De NMBS gebruikt stukken spoor van een lengte van 50m, die in de zomer 1cm uitzetten. Verklaar waarom de treinen in de zomer vertraging oplopen. (hint: Een spoor dat uitzet buigt en komt van de grond, berken de hoogte door te benaderen met driehoeken)
2.2
Scripts
Je hebt het reeds gemerkt, vaak moeten een aantal commando’s achtereenvolgens worden uitgevoerd. Men kan zo’n reeks commando’s saven als een script. Je maakt een nieuwe script aan via het file-menu. We geven hier het voorbeeld van het oplossen van een vierkantsvergelijking ax2 +bx+c = 0. a=1; b=-1; c=-1; D=b^2-4*a*c; 4
x1=(-b+sqrt(D))/(2*a) x2=(-b-sqrt(D))/(2*a) Opmerking 3. Let op het gebruik van ; om het wel of niet tonen van de opdrachten uit een script te regelen.
2.3
Vectoren en matrices
Alles is een vector in MATLAB . Ook functies zijn eigenlijk vectoren van functiewaarden. >> x=-2*pi:.1:2*pi; >> y=sin(x) y = Columns 1 through 10 0.0000
0.0998
0.1987
0.2955
0.3894
0.4794
...
0.9320
0.9636
0.9854
0.9975
...
Columns 11 through 20 0.8415
0.8912
... Opmerking 4. Wanneer je werkt met grote vectoren is het gebruik van ; zeer nuttig. Ook getallen zijn eigenlijk 1 × 1-vectoren. >> size(x) ans = 1
126
>> size(pi) ans = 1
1 5
Naast vectoren kent MATLAB ook matrices. Die voer je als volgt in. >> A=[1,2,3;4,5,6;7,8,9] A = 1 4 7
2 5 8
3 6 9
Je kan rijen en kolommen uit een matrix halen. >> A(2,1:3) ans = 4
5
6
>> A(1:3,2) ans = 2 5 8 Opmerking 5. Let op! Het gebruik van , en ; is hier anders. We gebruiken , om elementen van een matrix te scheiden, terwijl ; de rijen van elkaar scheidt. 2.3.1
Bewerkingen
>> A=[3,-1;5,-2] A = 3 -1 5 -2 >> B=[2,4;3,-5] B = 2 4 3 -5 Je kan de gewone bewerkingen gebruiken voor matrices. 6
>> A*B,A^-1,B^2 ans = 3 17 4 30 ans = 2.0000 -1.0000 5.0000 -3.0000 ans = 16 -12 -9 37 Opmerking 6. Een bewerking die we niet gewoon zijn in de wiskunde is het delen van matrices. In MATLAB gaat dit wel. A/B betekent hier A ∗ B −1 . >> B/A ans = 24.0000 -19.0000 >> B*A^-1 ans = 24.0000 -19.0000
-14.0000 12.0000
-14.0000 12.0000
Opmerking 7. De bewerkingen *, ^ bestaan ook in een ”gepunte” versie .*, .^ Deze versies werken element per element. >> A.*B,A.^-1,B.^2 ans = 6 -4 15 10 ans = 0.3333 -1.0000 0.2000 -0.5000 ans = 4 16 9 25 2.3.2
Stelsels oplossen
Eenmaal je matrices kan invoeren, kunnen natuurlijk ook stelsels lineaire vergelijkingen worden opgelost. Daartoe gebruik je het commando linsolve.
7
Let wel, gezien het feit dat MATLAB numerieke methoden gebruikt, zullen enkel stelsels met een reguliere matrix tot een oplossing leiden. Het stelsel ( 2x + 4y = 3 3x + 8y = 7 los je als volgt op. >> C=[2,4;3,8]; D=[3;7]; >> linsolve(C,D) ans = -1.0000 1.2500 We hadden dit ook korter gekund door C^-1*D te gebruiken. Oefening 2. Gebruik een Vandermonde matrix (vander) om een veelterm van graad n − 1 te bepalen door n punten. 2.3.3
Eigenwaarden en eigenvectoren
We hernemen de matrix A van voorheen. en gaan opzoek naar eigenwaarden en eigenvectoren. Het commando [V,D]=eig(A) geeft de diagonaalmatrix D, die de eigenwaarden bevat, terug, samen met de overgangsmatrix V , die de eigenvectoren als kolommen heeft. We weten dat A = V DV −1 . >> A=[3,-1;5,-2] A = 3 -1 5 -2 >> [V,D]=eig(A) V = 0.5862 0.2664 0.8101 0.9639 D = 1.6180 0 0 -0.6180 >> V*D*V^-1 ans = 3.0000 -1.0000 5.0000 -2.0000 8
De eigenvectoren uit V zijn genormeerd. We kunnen nagaan dat de eerste eigenvector, wel degelijk bij de eerste eigenwaarde hoort. >> E=V(1:2,1) E = 0.5862 0.8101 >> norm(E) ans = 1.0000 >> A*E ans = 0.9485 1.3108 >> D(1,1)*E ans = 0.9485 1.3108 Oefening 3. Bepaal de eigenwaarden en eigenvectoren van een draaing, een homothetie, een spiegeling en een puntspiegeling in het vlak.
2.4 2.4.1
Grafieken 2D-Grafieken
De grafiek van een functie maken, herleidt zich tot het berekenen van twee vectoren. We geven eerst het voorbeeld van de sech functie. >> x=-6:.01:6; >> y=sech(x); >> plot(x,y) Voor de functie f : R → R : x 7→ 21/x wordt dit het volgende. Let op het gebruik van .^ en ./! >> x=-6:.01:6; >> y=2.^(1./x); >> plot(x,y) Je ziet onmiddellijk dat de grafiek niet bruikbaar is, hiervoor moeten we eerst de assen aanpassen. >> axis([-6,6,-2,5]) 9
Je kan dus de eigenschappen van de grafiek aanpassen nadat je de grafiek hebt gemaakt. Om meerdere grafieken op een enkele tekening te krijgen gebruik je hold all. Hiermee worden alle eigenschappen van een grafiek behouden en zullen de verdere grafieken bijgevoegd worden in een andere stijl. Een andere mogelijkheid is om bij plot meerdere x,y vectoren mee te geven. >> >> >> >> >> >> >> >>
x=-6:.01:6; y1=x.^2-6*x+8; plot(x,y1,’LineWidth’,2) axis([-6,6,-2,10]) hold all y2=abs(x.^2-6*x+8); y3=abs(x).^2-6*abs(x)+8; plot(x,y2,x,y3)
Oefening 4. Elke functie f is de som fe + fo van een even en een oneven functie: f (x) + f (−x) fe (x) = 2 f (x) − f (−x) fo (x) = 2 Gebruik dit om een aantal functies op te splitsen in een even en een oneven functie. Maak tevens de grafieken. Oefening 5. Gebruik MATLAB om de hypergoniometrische functies te onderzoeken. Hierbij is sinh het oneven deel van ex en cosh het even deel. Verder definieert men de hyperbolische tangens, cotangens, secans en cosecans analoog als in de goniometrie. Oefening 6. Als je een ketting tussen twee punten (−1, 1) en (1, 1) laat hangen zou je denken dat de ketting volgens een parabool hangt. Niets is minder waar, de werkelijke kromme is een catenaire of kettinglijn. De vergelijking is een verschuiving van een hyperbolische cosinus: y = cosh(x) + c Bepaal c zodat de kromme door de gewenste punten gaat. Maak de grafiek en bepaal het minimum. Vergelijk met een parabool die door de eindpunten en door het minimum gaat. Maak hiervoor een grafiek van beide krommen. Oefening 7. Beschouw de functies f (x), f (|x|) en |f (x)|. Wat is het verband tussen hun grafieken? 10
Je kan in MATLAB ook grafieken maken in poolco¨ordinaten. >> >> >> >> >>
t=0:.01:2*pi; r=t/2; [x,y]=pol2cart(t,r); plot(x,y) axis equal
Oefening 8. Onderzoek de volgende polaire grafieken. r = sin(aθ), θ ∈ [0, 2π[ 1 r = 1 + ,θ > 0 θ Ook parameterkrommen kunnen getekend worden. >> >> >> >>
t=0:.01:2*pi; x=2*cos(t); y=3*sin(t); plot(x,y)
Oefening 9. Onderzoek de Lissajous-krommen. ( x = sin(at) , t ∈ [0, 2π[ y = sin(bt) 2.4.2
3D-Grafieken
Om 3D-grafieken te maken vertrek je van een meshgrid. Dit geeft twee matrices die x en y co¨ordinaten bevatten voor een rechthoek van punten in het xy-vlak. Je kan een of twee parameters meegeven aan meshgrid (voor een vierkantig, resp. rechthoekig domein). Daarna gebruik je mesh, surf of contour voor de passende 3D-grafiek. >> >> >> >> >>
[x,y]=meshgrid(-3:.1:3); z=x.*exp(-x.^2-y.^2); mesh(x,y,z) surf(x,y,z) contour(x,y,z)
Soms is het handiger om bij mesh een zwart-wit grafiek te tekenen. Met dezelfde optie kan je ook bij surf de wireframe weglaten.
11
>> mesh(x,y,z,’EdgeColor’,’black’) >> surf(x,y,z,’EdgeColor’,’none’) Met colormap kan je kiezen het kleurepalet aan te passen. >> >> >> >>
surf(x,y,z) colormap hot colormap cool colormap gray
Je kan tenslotte ook kiezen om de kleur te laten afhangen van een functie. >> c=sqrt(x.^2+y.^2); >> surf(x,y,z,c) Oefening 10. Maak de grafiek van de hyperboloide gegeven door f (x, y) = x2 − y 2 . Bepaal het raakvlak in (2, 3, f (2, 3)) en teken die op dezelfde grafiek. Gebruik eerst een surf-plot zonder wireframe, daarna een mesh-plot in zwartwit. Gebruik ten laatste hidden off. Wat doet dit? √ Oefening 11. Maak de grafiek van f (x, y) = 3 xy. Wat is hier het probleem? Hoe komt dit? Waarom biedt z=sign(x.*y).*(abs(x.*y)).^(1/3) een oplossing? p Oefening 12. De grafiek van f (x, y) = 4 − x2 − y 2 is een halve bol. Maak deze grafiek. Los de problemen die je tegenkomt zelf op. Je kan ook in 3D andere co¨ordinatensystemen gebruiken. Voor sferische coordinaten doe je dit als volgt. >> >> >> >>
[theta,phi]=meshgrid(0:.1:2*pi,-pi/2:.1:pi/2); r=1; [x,y,z]=sph2cart(theta,phi,r); surf(x,y,z)
Ook cylindrische co¨ordinaten kan je gebruiken. >> >> >> >>
[theta,z]=meshgrid(0:.1:2*pi,-6:.5:6); r=1; [x,y,z]=pol2cart(theta,r,z); surf(x,y,z)
Oefening 13. Waar ken je de vorm r = θ in sferische co¨ordinaten van?
12
Oefening 14. Om het omwentelingslichaam gegeven door een 2D-kromme y = f (x) te bekomen kan je gebruik maken van cylindrische co¨ordinaten. Je gebruikt hiervoor het oppervlak gegeven door r = f (z). Doe dit voor volgende omwentelingslichamen: een kegel, een bol en ellipsoide. Oefening 15. De Belgische natuurkundige Joseph Plateau (1801-1883) heeft ge¨experimenteerd met zeepfilms. Hij zocht welke oppervlakte met minimale spanning men bekwam door zeepfilms tussen twee cirkels te maken. Het antwoord hierop is een catenoide, het omwentelingslichaam van een kettinglijn. Gebruik cylindrische co¨ordinaten om een catenoide te maken. We kunnen tenslotte MATLAB ook gebruiken om oppervlakken in parameterco¨ordinaten te tekenen. >> >> >> >> >> >>
[u,v]=meshgrid(0:.1:2*pi); R=3;r=1; x=(R+r*cos(u)).*cos(v); y=(R+r*cos(u)).*sin(v); z=R*sin(u); mesh(x,y,z)
Oefening 16. Een regeloppervlak wordt geparametriseerd door x = (1 − u)px (t) + uqx (t) y = (1 − u)py (t) + uqy (t) , u ∈ [0, 1] z = (1 − u)pz (t) + uqz (t) waarbij p(t) = (px (t), py (t), pz (t)) en q(t) = (qx (t), qy (t), qz (t)) twee krommen in R3 zijn. Enkele voorbeelden: • Cylinder: twee evenwijdige cirkels. p(t) = (cos(t), sin(t), −1), q(t) = (cos(t), sin(t), 1), t ∈ [0, 2π] • Hyperboloide: twee evenwijdige cirkels doorlopen met een faseverschil. p(t) = (cos(t), sin(t), −1), q(t) = (cos(t+π/2), sin(t+π/2), 1), t ∈ [0, 2π] • M¨obiusband: p(t) = ((2 + cos(t)) cos(2t), (2 + cos(t)) sin(2t), sin(t)) q(t) = ((2 − cos(t)) cos(2t), (2 − cos(t)) sin(2t), − sin(t)) t ∈ [0, 2π] 13
• Helicoide: p(t) = (cos(2t), sin(2t), t), q(t) = (0, 0, t), t ∈ [0, 2π] Je kan tenslotte ook krommen in de ruimte tekenen. Dit doe je met plot3. >> >> >> >> >>
t=0:.01:1; x=cos(2*pi*t); y=sin(2*pi*t); z=t; plot3(x,y,z)
Oefening 17. Maak de grafieken uit de vorige oefening opnieuw, maar nu met de randen bijgetekend. Oefening 18. Onderzoek een 3D versie van de Lissajous-krommen. 2.4.3
Animaties
Het is ook mogelijk om met MATLAB animaties te maken. Dit wordt gedaan door een aantal grafieken uit te rekenen en die achtereenvolgens te tekenen en op te slaan met getframe. Het afspelen zelf gebeurd met movie. We gebruiken voor dit alles een for-lus, die uitgevoerd zal worden als je ze sluit met end. >> for k=1:50 t=0:.05:2*pi; e=-2+k/50*4; r=1./(1+e.*cos(t)); [x,y]=pol2cart(t,r); plot(x,y) axis([-3,3,-3,3]); M(k)=getframe; end >> movie(M) Je mag niet vergeten om na elke grafiek de assen juist te zetten. Voordat je movie aanroept is het best om elk bestaand grafisch venster te sluiten. Met movie(M,30) herhaal je het filmpje dertig keer. Oefening 19. Maak een animatie die het oppervlak z = f (x, y) continu omvormt tot het oppervlak z = g(x, y). Werk een aantal voorbeelden uit. Oefening 20. Maak een script die voor de regeloppervlakken uit vorige oefeningen een animatie maakt waarbij het regeloppervlak rechte per rechte verschijnt. 14
2.5
Veeltermen
Naast numerieke computerprogramma’s bestaan er ook symbolische programma’s. MATLAB werkt numeriek. Om te laten zien wat symbolische berekeningen zijn zullen we in eerste instantie gaan kijken naar veeltermen. Voor MATLAB is een veelterm een vector van co¨effici¨enten. >> v=[1,0,-2]; >> polyval(v,[-1,0,1,sqrt(2),2]) ans = -1.0000 -2.0000 -1.0000 >> roots(v) ans = 1.4142 -1.4142 >> poly([sqrt(2),-sqrt(2)]) ans = 1.0000 0 -2.0000
0.0000
2.0000
Met poly kan ook je de characteristieke veelterm van een matrix vinden. En daarmee dus ook de eigenwaarden. >> A=[3,-1;5,-2]; >> p=poly(A) p = 1.0000 -1.0000 >> roots(p) ans = 1.6180 -0.6180
-1.0000
Met conv kan je het product van veeltermen uitrekenen. >> poly([1,2,3]) ans = 1 -6 11 -6 >> conv(conv([1,-1],[1,-2]),[1,-3]) ans = 1 -6 11 -6 >> roots(ans) ans = 3.0000 2.0000 1.0000 15
MATLAB kan tenslotte ook veeltermen integreren en afleiden. >> p=poly([1,2,3]) p = 1 -6 11 >> q=polyint(p) q = 0.2500 -2.0000 >> polyder(q) ans = 1 -6 11
-6
5.5000
-6.0000
0
-6
Oefening 21. Gebruik polyval om de grafiek van een veelterm te maken. Oefening 22. In de oude Islamkultuur worstelde men met volgend probleem: vind getallen m, n zodat m + n = S en mn = P indien S en P gegeven zijn. Men wist dat m en n de oplossingen zijn van de vierkantsvergelijking x2 − Sx + P = 0. Gebruik deze methode om m en n te bepalen voor S = 10, P = 16 en S = 2000, P = 802864. Waarom had men problemen voor S = 10, P = 90? Controleer je antwoorden. Oefening 23. Maak met MATLAB een derdegraadsveelterm V (x) die in x = 1, 2 en 3 de waarde 6 heeft. Maak een grafiek van V (x) en de rechte y = 6. (hint: beschouw V(x)-6) Oefening 24. Bepaal de oppervlakte tussen de parabolen y = x2 + 3x en y = −x2 + x + 4. Oefening 25. Gebruik polyfit om een veelterm van graad n − 1 te bepalen die door n punten gaat. Vergelijk met de oefening op Vandermonde matrices. Wat gebeurt er als je de graad van de veelterm doet dalen? Maak ook telkens de nodige grafieken. Oefening 26. Implementeer het algoritme van Newton-Raphson om nulwaarden te vinden voor veeltermen. Gebruik dit om de nulwaarden te benaderen van 2x3 − x2 − 3x − 1.
2.6
Symbolische wiskunde
De functionaliteit van MATLAB kan uitgebreid worden dmv toolboxes. Een ervan is de Symbolic Math Toolbox dat MATLAB in staat stelt om echte symbolische berekeningen te doen. Deze toolbox is gebaseerd op Mupad dat een programma is om symbolische wiskunde te doen en dat gratis te downloaden is op internet. 16
Een variabele in MATLAB moet eerst symbolisch gemaakt worden alvorens ze te kunnen gebruiken in symbolische berekeningen. >> syms x y z >> (x+y)*x-z+(x-z)*(1-y) ans = x*(x + y) - z - (x - z)*(y - 1) >> simplify(ans) ans = x^2 + x - 2*z + y*z Soms is het ook nodig om getallen in symbolische vorm om te zetten. >> sym(’1/sqrt(2)’) ans = 2^(1/2)/2 >> pretty(ans) 1/2 2 ---2 Verder kan MATLAB hiermee ook berekeningen uitvoeren. >> syms a b c x >> p=a*x^2+b*x+c p = a*x^2 + b*x + c >> solve(p,x) ans = -(b + (b^2 - 4*a*c)^(1/2))/(2*a) -(b - (b^2 - 4*a*c)^(1/2))/(2*a) >> diff(p,x) ans = b + 2*a*x >> diff(p,x,2) ans = 2*a >> int(p,x) ans = (a*x^3)/3 + (b*x^2)/2 + c*x Ook matrices kunnen in symbolische vorm voorkomen. Hiermee kan je de exacte waarden van de eigenwaarden vinden ipv een benadering! 17
>> A=sym(’[3,-1;5,-2]’) A = [ 3, -1] [ 5, -2] >> [V,E]=eig(A) V = [ 1/2 - 5^(1/2)/10, 5^(1/2)/10 + 1/2] [ 1, 1] E = [ 1/2 - 5^(1/2)/2, 0] [ 0, 5^(1/2)/2 + 1/2] De symbolische toolbox kent ook de functies ezplot, ezplot3, ezsurf, ezcontour en ezmesh om grafieken van symbolische uitdrukkingen te maken. >> syms x y >> f=exp(-x^2-y^2) f = 1/exp(x^2 + y^2) >> int(int(f,x,-inf,inf),y,-inf,inf) ans = pi >> ezmesh(f,[-3,3,-3,3]) Oefening 27. Gebruik factor en expand om een aantal veeltermen te ontbinden in factoren en uit te werken. Oefening 28. Gebruik de symbolische toolbox om de functie f : R → R : x 7→ 21/x zo volledig mogelijk te onderzoeken. Oefening 29. Bepaal de oppervlakte tussen de grafieken van de hyperbolische sinus, hyperbolische cosinus en de y-as. Bepaal ook de oppervlakte onder de hyperbolische secans. Oefening 30. Gebruik taylor om achtereenvolgende benaderingen te maken van de sinus functie. Maak de nodige grafieken. Oefening 31. Zoek het verband tussen de volumes van een cilinder, een kegel en een paraboloide met eenzelfde grondvlak en gelijke hoogten.
3
Programmeren
Er bestaan verscheidene soorten programmeertalen: 18
• 1960: spaghetti programeren (BASIC) • 1970: procedureel programmeren (PASCAL) • 1990: object georienteerd programmeren (JAVA) MATLAB beschikt over een volwaardige procedurele programmeertaal. Iedere programmeertaal beschikt over volgende commando structuren: 1. variabelen (a,b,c, ...) en elementaire bewerkingen (+,*,sin, ¡) 2. voorwaardelijke structuren: als a < b dan (doe dit) anders (doe dat) 3. lussen: voor k = 1 to 10 (doe dit) 4. commando’s inherent aan de taal zelf: Java: internet MATLAB : wiskunde Variabelen en bewerkingen heb je ondertussen al onder de knie, en de for-lus ben je al tegengekomen. De if-structuur ziet er als volgt uit. >> x=35; >> relpriem=[]; >> for k=1:x if gcd(k,x)==1 relpriem=[relpriem, k]; end end >> relpriem relpriem = Columns 1 through 17 1 2 3 4 Columns 18 through 24 26 27 29 31
3.1
6
8
9
32
33
34
11
...
Functies
Soms is het handig om zelf een functie te kunnen schrijven die een uitkomst teruggeeft. Dit doe je ook via het file-menu.
19
function [ c ] = sommod( a,b,n ) % Som modulo n % sommod(a,b,n) geeft a+b mod n weer. c=mod(a+b,n); end De twee lijnen met commentaar zijn optioneel, maar worden wel gebruikt door het helpsysteem van MATLAB . Probeer maar eens help sommod of doc sommod. Het aanroepen gebeurt zoals je zou verwachten. >> sommod(3,4,5) ans = 2 Oefening 32. Schrijf een functie die n! uitrekent. Doe dit zowel iteratief als recursief. Oefening 33. Schrijf een functie die het nde Fibonacci getal uitrekent. Doe dit zowel iteratief als recursief. Ga na dat de verhouding van twee opeenvolgende Fibonacci getallen naar de gulden snede convergeert. Oefening 34. Schrijf een functie die nagaat of een getal in een bepaalde verzameling zit. Een functie doorgeven aan een andere functie, gebeurt via een function handle die je aanmaakt met @. Onderstaande functie gaat na of de bewerking b op de verzameling S commutatief is.
20
function [ res ] = iscomm( S, b ) % iscomm(S,@b) gaat na of b een comm bewerking is res=true; for k=1:length(S) for l=1:length(S) if not(b(S(k),S(l))==b(S(l),S(k))) res=false; break; end end if not(res) break; end end end Nu maken we een bewerking aan. function [ c ] = bew( a,b ) c=sommod(a,b,5); end We kunnen nu de functie iscomm oproepen. >> S=0:4 S = 0 1 2 >> iscomm(S,@bew) ans = 1
3
4
Oefening 35. Maak een functie die voor een bepaalde bewerking op een verzameling een Cayley-tabel maakt. Oefening 36. Schrijf een verzameling functies die nagaan of een gegeven bewerking op een gegeven verzameling een commutatieve groep is.
3.2
Verdere oefeningen
Oefening 37. Maak een script die de driehoek van Pascal genereert (< 50 lijnen, wat gebeurt er voor meer lijnen?). Zet dan alle p-vouden op 1 en alle andere op 1 voor een bepaald priemgetal p (begin met 2). Maak dan de grafiek, welk beeld krijg je? 21
Oefening 38. Onderstaand script maakt de Mandelbrot fractaal. Ga na wat elke regel doet. Leg het algoritme uit. clear for k=1:500 x=-2.5+k/500*4; for l=1:500 y=-1.5+l/500*3; z=x+y*i; c=z; for it=1:50 z=z^2+c; if abs(z)>2 M(k,l)=it; break end end end end image(M) axis equal Oefening 39. Maak een animatie van p een kromme in poolcoordinaten die afhangt van een parameter, vb: r = a2 − cos(2θ). Oefening 40. Maak een aantal scripts om π te benaderen. (hint: arctangens (taylorreeks), oppervlakte cirkel, Monte-Carlo methode, integratie ... bedenk maar iets) Oefening 41. Schrijf een procedure die het beeld van een kromme in C maakt voor de Moebius-transformatie z 7→
4z − 6 2z − 4
Bepaal hiermee de beelden van de krommen a(t) = 2 + 2eIπt , t = 0..1 en b(t) = 4 − (2 + 4I)t, t = 0..1 en c(t) = (2 − 4I)t , t = 0..1. Maak een tekening van de krommen en van hun beelden. Oefening 42. Maak een functie die, gegeven een dataset, de statistieken ervan teruggeeft. Doe dit eventueel ook voor twee dimensionale statistieken. Oefening 43. Onderstaand script maakt de Koch kromme. Het eerste deel maakt een ”programma” aan het tweede deel gebruikt dit om een turtlegrafiek te maken. Ga na wat elke regel doet. Leg het algoritme uit. 22
% Koch kromme L=[1]; S=[1,60,1,-120,1,60,1]; h=1/3; n=5; % substituties s=1; for k=1:n NL=[]; s=s*h; for l=1:length(L) if L(l)==1 NL=[NL,S]; else NL=[NL,L(l)]; end end L=NL; end % turtle grafiek x=0;y=0; t=0; for l=1:length(L) x1=x;y1=y; if L(l)==1 x=x+s*cos(t*pi/180);y=y+s*sin(t*pi/180); line([x1,x],[y1,y],’Color’,’black’) else t=t+L(l); end end axis([0,1,-.1,1])
23
