Tillämpad matematik TMA 990 KB2 - 2002/2003

Kölista

Datorstudio 1: Fourierserier.

Syftet med denna övning är att skriva ett Matlabprogram som beräknar approximationer till Fourierserier av en godtycklig periodisk funktion. Genom att plotta olika delsummor av serien kan man också få en visuell förståelse för konvergensen hos Fourierserier. Titta i kapitel 3.1-3 i Jan Peterssons Fourieranalysbok för mer detaljerade definitioner och mer teori.

Problemformulering.

Givet en funktion f(t) med perioden T ska ni skriva en Matlabfunktion som beräknar delsumman

S_N(t) = 1/2*a_0 + sum_n( a_n*cos(n*W*t)+b_n*sin(n*W*t) ), 1 <= n <= N

där W=2*pi/T betecknar grundvinkelfrekvensen. S_N(t) ger då en approximation till f(t) som konvergerar punktvis (i alla punkter där f är dervierbar) då N går mot oändligheten. Vi har att

a_n = 2/T int_0^T f(t)*cos(n*W*t)dt
b_n = 2/T int_0^T f(t)*sin(n*W*t)dt

f(t) kan vara besvärlig att skriva med ett enkelt uttryck och ska därför implementeras i en egen Matlabfil som ges ett passande namn; i instruktionerna har den fått namnet f.m.
Er Matlabfunktion ska som indata ta namnet på en funktion 'f', funktionens period T, en vektor t med t-värden samt antalet termer i delsumman N. Utdata ska vara en vektor S_N med delsummans värden för varje t-värde.
Koefficienterna a_n och b_n beräknas genom att loopa från 0 till N och räkna ut integralerna mha av enkel kvadratur. Termerna i summan adderas lämpligen ihop allt eftersom koefficienterna beräknas.

Programstruktur.

Nedan ges en struktur att följa när ni skriver ert program:

Indata: 'f', t, T, N
Utdata: S_N

Struktur:
Initialisera
S_N ska ha samma storlek som t
noQpnts antalet kvadraturpunkter/delintervall
W grundvinkelfrekvensen

För varje n = 0:N
beräkna a_n och b_n mha kvadratur
addera motsvarande termer till S_N

Att tänka på:
f.m måste skrivas så att den kan ta en vektor t som indata och returnera funktionsvärden för varje komponent i t
Antalet kvadraturpunkter vi använder är känsligt, varför? Hur ska vi välja detta?

Uppgifter.

Skriv programmet. Testa först att det ger rätt svar om ni ger en sinus- eller cosinusfunktion som indata (serien ska ju då endast innehålla denna komponent och sammanfaller då med funktionen).

Prova sedan med en fyrkantsvåg enligt Exempel 1 på sidan 3:07 i boken. Fyrkantsvågen kan implementeras som 


function y = fyrkant(x)
xx = mod(x,2*pi);
y = (abs(xx-pi/2) < pi/2)-(abs(xx-3*pi/2) < pi/2);

Plotta 'fyrkant' och S_N i samma figur. Kan ni återskapa bliderna i Figur 7 på sidan 3:09?

Prova med andra funktioner. Fundera över följande frågor: Vad är S_0? Vad händer när ni ökar N? Vad händer om antalet kvadraturpunkter är för få?

Utöka programmet så att du plottar S_N varje gång du adderat nya termer i loopen. Börja med att plotta 'f' inne i programmet, plotta sedan S_N varje gång i loopen (ta bort den gamla S_N innan du plottar den nya, se tipsen nedan).

Räkna ut och plotta ett amplitudspektrum.

Om man använder ett stort N och många kvadraturpunkter så blir programmet väldigt långsamt. Försök skriva om integralberäkningen så att ni inte behöver använda en loop. Finns det fler sätt att få programmet att gå snabbare?

Tips!

1)
Kommandot xx = mod(x,T) gör att xx blir en periodisk vektor med värden i intervallet [0,T), dvs x = xx + m*T för lämpliga heltal m. Är användbart när man ska skapa periodiska funktioner.
Exempel:


>> x=0:.1:1 

x =

         0    0.1000    0.2000    0.3000    0.4000    0.5000    0.6000    0.7000    0.8000    0.9000    1.0000

>> xx=mod(x,.3)

xx =

         0    0.1000    0.2000         0    0.1000    0.2000         0    0.1000    0.2000         0    0.1000

2)
d = plot(x,y) ger ett handtag d till plotten så att du kan ta bort just denna graf senare mha delete(d).
Exempel:


x=0:.1:2;
plot(x,x.^2,'r')
hold
d1=plot(x,sin(x))
d2=plot(x,sin(2*x),'g')
delete(d1)

3)
stegfunktioner kan skrivas mha boolska operationer (större än, mindre än etc) Exempel:


x=-2:.1:2;
H=(x>=0);
plot(x,H,'r')

(Se även funktionen fyrkant ovan.)

4)
subplot kan användas för att plotta flera figurer i samma fönster; läs mer genom att skriva help subplot

Lycka till!!!

P.S.
Så här skulle en lösning kunna se ut. Ladda hem och jämför med hur det ser ut när ni kör ert program. Här finns förmodligen en del "kodningstrick" som ni inte är vana vid så försök inte kopiera min kod.


Tillbaka till kurshemsidan.


Editor: Rickard Bergström
Last modified: 2002-10-23